Conference turris::decc

2203.0. "OS96 converter kinks" by KERNEL::SMITH () Thu May 29 1997 06:34

A cust has noted the following differences between dec and sun. The program
    runs ok on the sun prducing smooth curves, while on the dec alpha, the
    curves are more angular than curved.
    
    Wondering if anyone else has come across this before ?
    Is there a way to smooth out his angles ?
    
    ewan - uk csc
    
    
================================================================================

compilers: 
both standard dec cc (whatever is installed) and gcc2.7.2.1 (dec)
Sun Sparcworks cc and gcc2.7.2.1 (sun)
machines: 
OSF1 imaw08 V3.2 148 alpha
OSF1 imas02 V3.2 148 alpha
SunOS boudin 5.5 Generic_103093-02 sun4m sparc SUNW,SPARCstation-20
Compiler flags: -g -Wall -O -pedantic -pg (gcc)
Compiler flags: -g ( dec & sun cc)

Symptoms: the OS96 backwards converter converts data from OS96 to OS93 citf 
formats. During this process, parallel railway tracks are generated from a 
center line. The sun executable produces the correct parallels, however the 
dec has kinks which are visible. Examples of the data produced are available 
but will be relatively meaningless to outsiders. Relevant portions of the 
source code are attached below.
I have managed to minimise the problem through defensive coding but I would 
still appreciate any information on known problems with floating point 
calculation on DEC unix systems.

LIST GenerateRailTracks(FEATURE *iGenFeature, DMU *iDmu)
{
    LIST oFeatureList = NULL;
    char **Attr;
    int AttrCount;
    int i;
    int OutFsn1, OutFsn2;
    FEATURE *Track1, *Track2;
    char *ClipAttr = NULL;
    GEOM1 *Start, *End;
    POLAR2D Polar;                                                              
    CART2D Offset;
    long dx, dy;

    GetAttrValues(&iGenFeature->attr, "O3", &Attr, &AttrCount);
    if (AttrCount == 2) {
OutFsn1 = atol(Attr[0]);
OutFsn2 = atol(Attr[1]);
    } else if (AttrCount == 1) {
OutFsn1 = iGenFeature->ftype.point.serial;
OutFsn2 = atol(Attr[0]);
    } else {
dbg(WARNING, "O3 count %d for fsn %ld", AttrCount,
iGenFeature->ftype.point.serial);
return (NULL);
    }

    Track1 = Malloc(sizeof(FEATURE));
    Track1->Action = iGenFeature->Action;
    Track1->CitfType = CITF96;
    Track1->type = LINE_REC;
    Track1->OutFc = RAILTRACK_93FC;
    Track1->attr.CitfType = CITF96;
/* shallow copy */
    Track1->attr = iGenFeature->attr;
    Track1->ftype.line.recdesc = LINE_REC;
    Track1->ftype.line.serial = OutFsn1;
    Track1->OutFsn = OutFsn1;
    strcpy(Track1->ftype.line.udate, iGenFeature->ftype.point.udate);
    Track1->ftype.line.ustatus = iGenFeature->ftype.point.ustatus;
    Track1->ftype.line.divider = iGenFeature->ftype.point.divider;
    Track1->ftype.line.divider2 = iGenFeature->ftype.point.divider2;

    Track1->geom = Malloc(sizeof(GEOMREC));
    Track1->geom->recdesc = iGenFeature->geom->recdesc;
    Track1->geom->gtype = 2;
    Track1->geom->CitfType = CITF96;
    Track1->geom->numcoord = iGenFeature->geom->numcoord;

    Track2 = Malloc(sizeof(FEATURE));
    Track2->Action = iGenFeature->Action;
    Track2->CitfType = CITF96;
    Track2->type = LINE_REC;
    Track2->OutFc = RAILTRACK_93FC;
/* shallow copy */
    Track2->attr = iGenFeature->attr;
    Track2->attr.CitfType = CITF96;
    Track2->ftype.line.recdesc = LINE_REC;
    Track2->ftype.line.serial = OutFsn2;
    Track2->OutFsn = OutFsn2;
    strcpy(Track2->ftype.line.udate, iGenFeature->ftype.point.udate);
    Track2->ftype.line.ustatus = iGenFeature->ftype.point.ustatus;
    Track2->ftype.line.divider = iGenFeature->ftype.point.divider;
    Track2->ftype.line.divider2 = iGenFeature->ftype.point.divider2;

    Track2->geom = Malloc(sizeof(GEOMREC));
    Track2->geom->recdesc = iGenFeature->geom->recdesc;
    Track2->geom->gtype = 2;
    Track2->geom->CitfType = CITF96;
    Track2->geom->numcoord = iGenFeature->geom->numcoord;

/* generate new geometry */
    Polar.dist = RAILTRACK_WIDTH / 2.0;
    if (Track1->geom->recdesc == GEOMETRY1_REC) {
Track1->geom->g1crd = Malloc(sizeof(GEOM1) * iGenFeature->geom->numcoord);
Track2->geom->g1crd = Malloc(sizeof(GEOM1) * iGenFeature->geom->numcoord);
    } else {
Track1->geom->g2crd = Malloc(sizeof(GEOM2) * iGenFeature->geom->numcoord);
Track2->geom->g2crd = Malloc(sizeof(GEOM2) * iGenFeature->geom->numcoord);
    }

    for (i = 0; i < iGenFeature->geom->numcoord-1; i++) {
Start = (iGenFeature->geom->g1crd) ? &iGenFeature->geom->g1crd[i] : (GEOM1 *)
&iGenFeature->geom->g2crd[i];
End = (iGenFeature->geom->g1crd) ? &iGenFeature->geom->g1crd[i + 1] : (GEOM1 *)
&iGenFeature->geom->g2crd[i + 1];
dx = End->x - Start->x;
dy = End->y - Start->y;
Polar.theta = (dx) ? atan(dy / dx) : M_PI_2;
Polar.theta = M_PI_2 - Polar.theta;
ComputePolarGeomOffset(&Polar, &Offset);

if (Track1->geom->recdesc == GEOMETRY1_REC) {
    Track1->geom->g1crd[i].x = Start->x - Offset.x;
    Track1->geom->g1crd[i].y = Start->y + Offset.y;
} else {
    Track1->geom->g2crd[i].x = Start->x - Offset.x;
    Track1->geom->g2crd[i].y = Start->y + Offset.y;
    Track1->geom->g2crd[i].z = iGenFeature->geom->g2crd[i].z;
}

if (Track2->geom->recdesc == GEOMETRY1_REC) {
    Track2->geom->g1crd[i].x = Start->x + Offset.x;
    Track2->geom->g1crd[i].y = Start->y - Offset.y;
} else {
    Track2->geom->g2crd[i].x = Start->x - Offset.x;
    Track2->geom->g2crd[i].y = Start->y + Offset.y;
    Track2->geom->g2crd[i].z = iGenFeature->geom->g2crd[i].z;
}
if (Track1->geom->recdesc == GEOMETRY1_REC) {
    Track1->geom->g1crd[i + 1].x = End->x - Offset.x;
    Track1->geom->g1crd[i + 1].y = End->y + Offset.y;
} else {
    Track1->geom->g2crd[i + 1].x = End->x - Offset.x;
    Track1->geom->g2crd[i + 1].y = End->y + Offset.y;
    Track1->geom->g2crd[i + 1].z = iGenFeature->geom->g2crd[i + 1].z;
}

if (Track2->geom->recdesc == GEOMETRY1_REC) {
    Track2->geom->g1crd[i + 1].x = End->x + Offset.x;
    Track2->geom->g1crd[i + 1].y = End->y - Offset.y;
} else {
    Track2->geom->g2crd[i + 1].x = End->x - Offset.x;
    Track2->geom->g2crd[i + 1].y = End->y + Offset.y;
    Track2->geom->g2crd[i + 1].z = iGenFeature->geom->g2crd[i + 1].z;
}
    }

    ClipAttr = GetAttrValue(&iGenFeature->attr, "PF");
    if (ClipAttr && atoi(ClipAttr)) {
ClipFeatureIntoDmu(Track1, iDmu);
ClipFeatureIntoDmu(Track2, iDmu);
    }
    oFeatureList = Consl(oFeatureList, Track2);
    oFeatureList = Consl(oFeatureList, Track1);

    return (oFeatureList);
}


static const double TDEG_IN_RAD = M_PI / 1800.0;

void
ComputePolarGeomOffset(POLAR2D *iPolar, CART2D *oPoint)
{
    oPoint->x = (long) rint(sin(iPolar->theta) * (double)iPolar->dist);
    oPoint->y = (long) rint(cos(iPolar->theta) * (double)iPolar->dist);
}


double
CitfOrientToLibm(double theta)
{
    return ((M_PI_2) - (theta * TDEG_IN_RAD));

}

It's pretty hard to tell from this what the problem is.  Ideally, your customer
could reduce this to a single statement where the inputs are the same on Sun
and Digital systems, but the outputs are different.

Note that "standard dec cc" on "OSF1 V3.2" is the MIPS-derived compiler,
not DECC.  To try with DEC C, use "cc -migrate".

As a wild guess, perhaps the program depends on support for denormal numbers.
Try compiling with one of the "-ieee" switches -- preferably
"-ieee_with_no_inexact" (if that's the right spelling, it's not in the V4
man page I have access to).  I don't know whether you can do that with gcc.

It's also possible that there's a problem with your math library calls.
For example:

    long dx, dy;
Polar.theta = (dx) ? atan(dy / dx) : M_PI_2;

depending on how "atan" is declared, this might pass a "long int" instead
of a "double".  The easiest way to be sure it gets the right argument type
would be to compile with "-std" (ANSI mode) and be sure to "#include <math.h>"

Looking again at this call -- are you sure you want an integer divide here?
It seems much more likely the program wants to say
Polar.theta = (dx) ? atan((double)dy / (double)dx) : M_PI_2;