IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 5090 for trunk/psLib/src/astro


Ignore:
Timestamp:
Sep 21, 2005, 4:47:16 PM (21 years ago)
Author:
gusciora
Message:

The nX terms for the polynomials now refer to the order, not the number
of terms.

File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/psLib/src/astro/psCoord.c

    r4991 r5090  
    1010*  @author GLG, MHPCC
    1111*
    12 *  @version $Revision: 1.86 $ $Name: not supported by cvs2svn $
    13 *  @date $Date: 2005-09-11 22:18:40 $
     12*  @version $Revision: 1.87 $ $Name: not supported by cvs2svn $
     13*  @date $Date: 2005-09-22 02:47:16 $
    1414*
    1515*  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
    … …  
    9090
    9191    A = transform->x->coeff[0][0];
    92     if (transform->x->nX >= 2) {
     92    if (transform->x->COOL_2D_nX >= 1) {
    9393        B = transform->x->coeff[1][0];
    9494    }
    95     if (transform->x->nY >= 2) {
     95    if (transform->x->COOL_2D_nY >= 1) {
    9696        C = transform->x->coeff[0][1];
    9797    }
    9898    D = transform->y->coeff[0][0];
    99     if (transform->y->nX >= 2) {
     99    if (transform->y->COOL_2D_nX >= 1) {
    100100        E = transform->y->coeff[1][0];
    101101    }
    102     if (transform->y->nY >= 2) {
     102    if (transform->y->COOL_2D_nY >= 1) {
    103103        F = transform->y->coeff[0][1];
    104104    }
    … …  
    154154    PS_ASSERT_PTR_NON_NULL(transform->y, 0);
    155155
    156     for (psS32 i=0;i<(transform->x->nX);i++) {
    157         for (psS32 j=0;j<(transform->x->nY);j++) {
     156    for (psS32 i=0;i<(1 + transform->x->COOL_2D_nX);i++) {
     157        for (psS32 j=0;j<(1 + transform->x->COOL_2D_nY);j++) {
    158158            if (transform->x->coeff[i][j] != 0.0) {
    159159                if (!(((i == 0) && (j == 0)) ||
    … …  
    166166    }
    167167
    168     for (psS32 i=0;i<(transform->y->nX);i++) {
    169         for (psS32 j=0;j<(transform->y->nY);j++) {
     168    for (psS32 i=0;i<(1 + transform->y->COOL_2D_nX);i++) {
     169        for (psS32 j=0;j<(1 + transform->y->COOL_2D_nY);j++) {
    170170            if (transform->y->coeff[i][j] != 0.0) {
    171171                if (!(((i == 0) && (j == 0)) ||
    … …  
    208208}
    209209
     210// XXX: Verify the order/nterms poly changes
    210211psPlaneTransform* psPlaneTransformAlloc(int n1, int n2)
    211212{
    … …  
    214215
    215216    psPlaneTransform *pt = psAlloc(sizeof(psPlaneTransform));
    216     pt->x = psPolynomial2DAlloc(n1, n2, PS_POLYNOMIAL_ORD);
    217     pt->y = psPolynomial2DAlloc(n1, n2, PS_POLYNOMIAL_ORD);
     217    pt->x = psPolynomial2DAlloc(n1-1, n2-1, PS_POLYNOMIAL_ORD);
     218    pt->y = psPolynomial2DAlloc(n1-1, n2-1, PS_POLYNOMIAL_ORD);
    218219
    219220    psMemSetDeallocator(pt, (psFreeFunc) planeTransformFree);
    … …  
    276277}
    277278
     279// XXX: Verify the order/nterms poly changes
    278280psPlaneDistort* psPlaneDistortAlloc(int n1, int n2, int n3, int n4)
    279281{
    … …  
    284286
    285287    psPlaneDistort *pt = psAlloc(sizeof(psPlaneDistort));
    286     pt->x = psPolynomial4DAlloc(n1, n2, n3, n4, PS_POLYNOMIAL_ORD);
    287     pt->y = psPolynomial4DAlloc(n1, n2, n3, n4, PS_POLYNOMIAL_ORD);
     288    pt->x = psPolynomial4DAlloc(n1-1, n2-1, n3-1, n4-1, PS_POLYNOMIAL_ORD);
     289    pt->y = psPolynomial4DAlloc(n1-1, n2-1, n3-1, n4-1, PS_POLYNOMIAL_ORD);
    288290
    289291    psMemSetDeallocator(pt, (psFreeFunc) planeDistortFree);
    … …  
    499501                                       psPolynomial2D *trans2)
    500502{
    501     //TRACE: printf("multiplyDPoly2D(%d %d: %d %d)\n", trans1->nX, trans1->nY, trans2->nX, trans2->nY);
    502     psS32 orderX = (trans1->nX + trans2->nX) - 1;
    503     psS32 orderY = (trans1->nY + trans2->nY) - 1;
    504 
    505     psPolynomial2D *out = psPolynomial2DAlloc(orderX, orderY, PS_POLYNOMIAL_ORD);
     503    //TRACE: printf("multiplyDPoly2D(%d %d: %d %d)\n", trans1->COOL_2D_nX, trans1->COOL_2D_nY, trans2->COOL_2D_nX, trans2->COOL_2D_nY);
     504    // XXX: OLD:  psS32 orderX = (trans1->nX + trans2->nX) - 1;
     505    // XXX: Verify this (the poly order change)
     506    psS32 orderX = trans1->COOL_2D_nX + trans2->COOL_2D_nX + 1;
     507    psS32 orderY = trans1->COOL_2D_nY + trans2->COOL_2D_nY + 1;
     508
     509    // XXX: Verify this (the poly order change)
     510    psPolynomial2D *out = psPolynomial2DAlloc(orderX-1, orderY-1, PS_POLYNOMIAL_ORD);
    506511    //TRACE: printf("Creating poly (%d, %d)\n", orderX, orderY);
    507     for (psS32 i = 0 ; i < out->nX; i++) {
    508         for (psS32 j = 0 ; j < out->nY; j++) {
     512    for (psS32 i = 0 ; i < (1 + out->COOL_2D_nX); i++) {
     513        for (psS32 j = 0 ; j < (1 + out->COOL_2D_nY); j++) {
    509514            out->coeff[i][j] = 0.0;
    510515            out->mask[i][j] = 0;
    … …  
    512517    }
    513518
    514     for (psS32 t1x = 0 ; t1x < trans1->nX ; t1x++) {
    515         for (psS32 t1y = 0 ; t1y < trans1->nY ; t1y++) {
     519    for (psS32 t1x = 0 ; t1x < (1 + trans1->COOL_2D_nX) ; t1x++) {
     520        for (psS32 t1y = 0 ; t1y < (1 + trans1->COOL_2D_nY) ; t1y++) {
    516521            if (0.0 != trans1->coeff[t1x][t1y]) {
    517                 for (psS32 t2x = 0 ; t2x < trans2->nX ; t2x++) {
    518                     for (psS32 t2y = 0 ; t2y < trans2->nY ; t2y++) {
     522                for (psS32 t2x = 0 ; t2x < (1 + trans2->COOL_2D_nX) ; t2x++) {
     523                    for (psS32 t2y = 0 ; t2y < (1 + trans2->COOL_2D_nY) ; t2y++) {
    519524                        /* Possible debug-only macro which checks these coords?
    520525                        if ((t1x+t2x) >= orderX)
    … …  
    547552    PS_ASSERT_PTR_NON_NULL(trans1, NULL);
    548553    PS_ASSERT_PTR_NON_NULL(trans2, NULL);
    549     //TRACE: printf("psPlaneTransformCombine(%d, %d, %d, %d: %d, %d, %d, %d)\n", trans1->x->nX, trans1->x->nY, trans1->y->nX, trans1->y->nY, trans2->x->nX, trans2->x->nY, trans2->y->nX, trans2->y->nY);
     554    //TRACE: printf("psPlaneTransformCombine(%d, %d, %d, %d: %d, %d, %d, %d)\n", trans1->x->COOL_2D_nX, trans1->x->COOL_2D_nY, trans1->y->COOL_2D_nX, trans1->y->COOL_2D_nY, trans2->x->COOL_2D_nX, trans2->x->COOL_2D_nY, trans2->y->COOL_2D_nX, trans2->y->COOL_2D_nY);
    550555    //
    551556    // Determine the size of the new psPlaneTransform.
    … …  
    553558    // PS_MAX(  Number of x terms in T2->x * number of x terms in T1->x,
    554559    //          Number of y terms in T2->x * number of x terms in T1->y,
    555     psS32 orderXnX = PS_MAX((trans2->x->nX * trans1->x->nX),
    556                             (trans2->x->nY * trans1->y->nX));
    557     psS32 orderXnY = PS_MAX((trans2->x->nX * trans1->x->nY),
    558                             (trans2->x->nY * trans1->y->nY));
    559 
    560     psS32 orderYnX = PS_MAX((trans2->y->nX * trans1->x->nX),
    561                             (trans2->y->nY * trans1->y->nX));
    562     psS32 orderYnY = PS_MAX((trans2->y->nX * trans1->x->nY),
    563                             (trans2->y->nY * trans1->y->nY));
     560    psS32 orderXnX = PS_MAX(((1 + trans2->x->COOL_2D_nX) * (1 + trans1->x->COOL_2D_nX)),
     561                            ((1 + trans2->x->COOL_2D_nY) * (1 + trans1->y->COOL_2D_nX)));
     562    psS32 orderXnY = PS_MAX(((1 + trans2->x->COOL_2D_nX) * (1 + trans1->x->COOL_2D_nY)),
     563                            ((1 + trans2->x->COOL_2D_nY) * (1 + trans1->y->COOL_2D_nY)));
     564
     565    psS32 orderYnX = PS_MAX(((1 + trans2->y->COOL_2D_nX) * (1 + trans1->x->COOL_2D_nX)),
     566                            ((1 + trans2->y->COOL_2D_nY) * (1 + trans1->y->COOL_2D_nX)));
     567    psS32 orderYnY = PS_MAX(((1 + trans2->y->COOL_2D_nX) * (1 + trans1->x->COOL_2D_nY)),
     568                            ((1 + trans2->y->COOL_2D_nY) * (1 + trans1->y->COOL_2D_nY)));
    564569    psS32 orderX = PS_MAX(orderXnX, orderYnX);
    565570    psS32 orderY = PS_MAX(orderXnY, orderYnY);
    … …  
    568573    // Allocate the new psPlaneTransform, if necessary.
    569574    //
     575    // XXX: rename, or verify, or recode, after the poly norder/nterm change.
    570576    psPlaneTransform *myPT = NULL;
    571577    if (out == NULL) {
    572578        myPT = psPlaneTransformAlloc(orderX, orderY);
    573579    } else {
    574         if ((out->x->nX == orderX) && (out->x->nY == orderY) &&
    575                 (out->y->nX == orderX) && (out->y->nY == orderY)) {
     580        if (((1 + out->x->COOL_2D_nX) == orderX) && ((1 + out->x->COOL_2D_nY) == orderY) &&
     581                ((1 + out->y->COOL_2D_nX) == orderX) && ((1 + out->y->COOL_2D_nY) == orderY)) {
    576582            myPT = out;
    577583        } else {
    … …  
    603609    //
    604610
    605     for (psS32 t2x = 0 ; t2x < trans2->x->nX ; t2x++) {
    606         for (psS32 t2y = 0 ; t2y < trans2->x->nY ; t2y++) {
    607             psPolynomial2D *currPoly = psPolynomial2DAlloc(1, 1, PS_POLYNOMIAL_ORD);
     611    for (psS32 t2x = 0 ; t2x < (1 + trans2->x->COOL_2D_nX) ; t2x++) {
     612        for (psS32 t2y = 0 ; t2y < (1 + trans2->x->COOL_2D_nX) ; t2y++) {
     613            // XXX: Verify this (poly order/nterm change)
     614            psPolynomial2D *currPoly = psPolynomial2DAlloc(1-1, 1-1, PS_POLYNOMIAL_ORD);
    608615
    609616            currPoly->coeff[0][0] = 1.0;
    … …  
    627634
    628635                // Set the appropriate coeffs in myPT->x
    629                 for (psS32 i = 0 ; i < currPoly->nX ; i++) {
    630                     for (psS32 j = 0 ; j < currPoly->nY ; j++) {
     636                for (psS32 i = 0 ; i < (1 + currPoly->COOL_2D_nX) ; i++) {
     637                    for (psS32 j = 0 ; j < (1 + currPoly->COOL_2D_nY) ; j++) {
    631638                        myPT->x->coeff[i][j]+= currPoly->coeff[i][j] * trans2->x->coeff[t2x][t2y];
    632639                    }
    … …  
    638645
    639646
    640     for (psS32 t2x = 0 ; t2x < trans2->y->nX ; t2x++) {
    641         for (psS32 t2y = 0 ; t2y < trans2->y->nY ; t2y++) {
    642             psPolynomial2D *currPoly = psPolynomial2DAlloc(1, 1, PS_POLYNOMIAL_ORD);
     647    for (psS32 t2x = 0 ; t2x < (1 + trans2->y->COOL_2D_nX) ; t2x++) {
     648        for (psS32 t2y = 0 ; t2y < (1 + trans2->y->COOL_2D_nY) ; t2y++) {
     649            // XXX: Verify this (poly order/nterm change)
     650            psPolynomial2D *currPoly = psPolynomial2DAlloc(1-1, 1-1, PS_POLYNOMIAL_ORD);
    643651            currPoly->coeff[0][0] = 1.0;
    644652            currPoly->mask[0][0] = 0;
    … …  
    662670
    663671                // Set the appropriate coeffs in myPT->x
    664                 for (psS32 i = 0 ; i < currPoly->nX ; i++) {
    665                     for (psS32 j = 0 ; j < currPoly->nY ; j++) {
     672                for (psS32 i = 0 ; i < (1 + currPoly->COOL_2D_nX) ; i++) {
     673                    for (psS32 j = 0 ; j < (1 + currPoly->COOL_2D_nY) ; j++) {
    666674                        myPT->y->coeff[i][j]+= currPoly->coeff[i][j] * trans2->y->coeff[t2x][t2y];
    667675                    }
    … …  
    694702
    695703    psS32 numCoords = PS_MIN(source->n, dest->n);
    696     psS32 order = PS_MAX(trans->x->nX, trans->x->nY);
    697     order = PS_MAX(order, trans->y->nX);
    698     order = PS_MAX(order, trans->y->nY);
     704    psS32 order = PS_MAX(trans->x->COOL_2D_nX, trans->x->COOL_2D_nY);
     705    order = PS_MAX(order, trans->y->COOL_2D_nX);
     706    order = PS_MAX(order, trans->y->COOL_2D_nY);
     707    // XXX: Verify, or recode, after the poly norder/nterm change.
     708    order++;
    699709
    700710    //
    701711    // Create fake polynomial to use in evaluation
    702712    //
    703     psPolynomial2D *fakePoly = psPolynomial2DAlloc(order, order, PS_POLYNOMIAL_ORD);
     713    // XXX: Verify this (poly order/nterm change)
     714    psPolynomial2D *fakePoly = psPolynomial2DAlloc(order-1, order-1, PS_POLYNOMIAL_ORD);
    704715    for (int i = 0; i < order; i++) {
    705716        for (int j = 0; j < order; j++) {
    … …  
    811822
    812823    // Ensure that the input transformation is symmetrical.
    813     if ((in->x->nX != in->x->nY) ||
    814             (in->y->nX != in->y->nY) ||
    815             (in->x->nX != in->y->nX)) {
     824    if ((in->x->COOL_2D_nX != in->x->COOL_2D_nY) ||
     825            (in->y->COOL_2D_nX != in->y->COOL_2D_nY) ||
     826            (in->x->COOL_2D_nX != in->y->COOL_2D_nX)) {
    816827        psError(PS_ERR_BAD_PARAMETER_TYPE, true, "Input transformation must have same nX==nY.");
    817828    }
    818     psS32 order = in->x->nX;
     829    // XXX: recode or verify after poly changes
     830    psS32 order = 1 + in->x->COOL_2D_nX;
    819831
    820832    psPlaneTransform *myPT = NULL;
    … …  
    826838    //
    827839    if (out == NULL) {
     840        // XXX: Verify this (poly order/nterm change)
    828841        myPT = psPlaneTransformAlloc(order, order);
    829842    } else {
    830         if ((out->x->nX == order) && (out->x->nY == order) &&
    831                 (out->y->nX == order) && (out->y->nY == order)) {
     843        // XXX: recode or verify after poly changes
     844        if (((1 + out->x->COOL_2D_nX) == order) && ((1 + out->x->COOL_2D_nY) == order) &&
     845                ((1 + out->y->COOL_2D_nX) == order) && ((1 + out->y->COOL_2D_nX) == order)) {
    832846            myPT = out;
    833847        } else {
    … …  
    840854    // Copy the input transform to myPT.
    841855    //
    842     for (psS32 i = 0 ; i < in->x->nX ; i++) {
    843         for (psS32 j = 0 ; j < in->x->nY ; j++) {
     856    for (psS32 i = 0 ; i < (1 + in->x->COOL_2D_nX) ; i++) {
     857        for (psS32 j = 0 ; j < (1 + in->x->COOL_2D_nY) ; j++) {
    844858            myPT->x->coeff[i][j] = in->x->coeff[i][j];
    845859        }
    846860    }
    847     for (psS32 i = 0 ; i < in->y->nX ; i++) {
    848         for (psS32 j = 0 ; j < in->y->nY ; j++) {
     861    for (psS32 i = 0 ; i < (1 + in->y->COOL_2D_nX) ; i++) {
     862        for (psS32 j = 0 ; j < (1 + in->y->COOL_2D_nY) ; j++) {
    849863            myPT->y->coeff[i][j] = in->y->coeff[i][j];
    850864        }
    … …  
    926940
    927941    // This loop starts at loop_x=1 since the derivative of the loop_x=0 terms are all 0.0
    928     for (psS32 loop_x = 1; loop_x < xPoly->nX; loop_x++) {
     942    for (psS32 loop_x = 1; loop_x < (1 + xPoly->COOL_2D_nX); loop_x++) {
    929943        ySum = 1.0;
    930         for (psS32 loop_y = 0; loop_y < xPoly->nY; loop_y++) {
     944        for (psS32 loop_y = 0; loop_y < (1 + xPoly->COOL_2D_nY); loop_y++) {
    931945            //
    932946            // For each iteration of the loop, we multiple the (x, y) coefficient
    … …  
    946960
    947961    // This loop starts at loop_y=1 since the derivative of the loop_y=0 terms are all 0.0
    948     for (psS32 loop_x = 0; loop_x < yPoly->nX; loop_x++) {
     962    for (psS32 loop_x = 0; loop_x < (1 + yPoly->COOL_2D_nX); loop_x++) {
    949963        ySum = 1.0;
    950         for (psS32 loop_y = 1; loop_y < yPoly->nY; loop_y++) {
     964        for (psS32 loop_y = 1; loop_y < (1 + yPoly->COOL_2D_nY); loop_y++) {
    951965            //
    952966            // For each iteration of the loop, we multiple the (x, y) coefficient
Note: See TracChangeset for help on using the changeset viewer.