IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 1192


Ignore:
Timestamp:
Jul 7, 2004, 1:25:08 PM (22 years ago)
Author:
gusciora
Message:

...

Location:
trunk/psLib/src
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • trunk/psLib/src/dataManip/psMinimize.c

    r1191 r1192  
    588588
    589589
    590 
    591 
    592 
    593 
    594 
    595 
    596 
    597 
    598 
    599 
    600 
    601 
    602 
    603590/******************************************************************************
    604 p_psBuildSums1D(x, sums, polyOrder): this routine calculates the powers of
     591p_psBuildSums1D(x, polyOrder, sums): this routine calculates the powers of
    605592input parameter "x" between 0 and input parameter polyOrder.  The result is
    606 returned as a psVector.
     593returned as a psVector sums.
    607594 *****************************************************************************/
    608 psVector *p_psBuildSums1D(double x,
    609                           int polyOrder)
     595void p_psBuildSums1D(double x,
     596                     int polyOrder,
     597                     psVector *sums)
    610598{
    611599    int       i = 0;
    612600    double    xSum = 0.0;
    613     psVector *sums = NULL;
    614 
    615     sums = psVectorAlloc(polyOrder+1, PS_TYPE_F32);
     601
    616602    xSum = 1.0;
    617603    for(i=0;i<=polyOrder;i++) {
    618         sums->data.F32[i] = xSum;
     604        sums->data.F64[i] = xSum;
    619605        xSum*= x;
    620606    }
    621     return(sums);
    622607}
    623608
    … …  
    639624        imageScalingFactors->data.F32[i] = (((float) 2*i) / ((float) x)) - 1.0;
    640625    }
     626
    641627    return(imageScalingFactors);
    642628}
    … …  
    761747    int j = 0;
    762748    int k = 0;
     749    psVector *xSums = NULL;
    763750
    764751    if (x->n != y->n) {
    … …  
    769756    }
    770757
    771     A = psImageAlloc(polyOrder, polyOrder, PS_TYPE_F64);
    772     ALUD = psImageAlloc(polyOrder, polyOrder, PS_TYPE_F64);
    773     B = psVectorAlloc(polyOrder, PS_TYPE_F64);
    774     coeffs = psVectorAlloc(polyOrder, PS_TYPE_F64);
    775     X = psVectorAlloc(x->n, PS_TYPE_F64);
     758    A       = psImageAlloc(polyOrder, polyOrder, PS_TYPE_F64);
     759    ALUD    = psImageAlloc(polyOrder, polyOrder, PS_TYPE_F64);
     760    B       = psVectorAlloc(polyOrder, PS_TYPE_F64);
     761    coeffs  = psVectorAlloc(polyOrder, PS_TYPE_F64);
     762    X       = psVectorAlloc(x->n, PS_TYPE_F64);
    776763    outPerm = psVectorAlloc(polyOrder, PS_TYPE_F64);
     764    xSums   = psVectorAlloc(1+2*polyOrder, PS_TYPE_F64);
    777765
    778766    // Initialize data structures.
    … …  
    791779
    792780    // Build the B and A data structs.
    793     for(k=0;k<(polyOrder);k++) {
    794         for (i=0;i<X->n;i++) {
    795             B->data.F64[k]+= y->data.F64[i] * (pow(X->data.F64[i], (float) k));
    796         }
    797     }
    798 
    799     for(k=0;k<(polyOrder);k++) {
    800         for(j=0;j<(polyOrder);j++) {
    801             for (i=0;i<X->n;i++) {
    802                 A->data.F64[k][j]+= pow(X->data.F64[i], (float) (k + j));
     781    for (i=0;i<X->n;i++) {
     782        p_psBuildSums1D(X->data.F64[i], 2*polyOrder, xSums);
     783
     784        for(k=0;k<(polyOrder);k++) {
     785            B->data.F64[k]+= y->data.F64[i] * xSums->data.F64[k];
     786        }
     787
     788        for(k=0;k<(polyOrder);k++) {
     789            for(j=0;j<(polyOrder);j++) {
     790                A->data.F64[k][j]+= xSums->data.F64[k+j];
    803791            }
    804792        }
    … …  
    818806    psFree(X);
    819807    psFree(outPerm);
     808    psFree(xSums);
    820809
    821810    return(NULL);

  • trunk/psLib/src/math/psMinimize.c

    r1191 r1192  
    588588
    589589
    590 
    591 
    592 
    593 
    594 
    595 
    596 
    597 
    598 
    599 
    600 
    601 
    602 
    603590/******************************************************************************
    604 p_psBuildSums1D(x, sums, polyOrder): this routine calculates the powers of
     591p_psBuildSums1D(x, polyOrder, sums): this routine calculates the powers of
    605592input parameter "x" between 0 and input parameter polyOrder.  The result is
    606 returned as a psVector.
     593returned as a psVector sums.
    607594 *****************************************************************************/
    608 psVector *p_psBuildSums1D(double x,
    609                           int polyOrder)
     595void p_psBuildSums1D(double x,
     596                     int polyOrder,
     597                     psVector *sums)
    610598{
    611599    int       i = 0;
    612600    double    xSum = 0.0;
    613     psVector *sums = NULL;
    614 
    615     sums = psVectorAlloc(polyOrder+1, PS_TYPE_F32);
     601
    616602    xSum = 1.0;
    617603    for(i=0;i<=polyOrder;i++) {
    618         sums->data.F32[i] = xSum;
     604        sums->data.F64[i] = xSum;
    619605        xSum*= x;
    620606    }
    621     return(sums);
    622607}
    623608
    … …  
    639624        imageScalingFactors->data.F32[i] = (((float) 2*i) / ((float) x)) - 1.0;
    640625    }
     626
    641627    return(imageScalingFactors);
    642628}
    … …  
    761747    int j = 0;
    762748    int k = 0;
     749    psVector *xSums = NULL;
    763750
    764751    if (x->n != y->n) {
    … …  
    769756    }
    770757
    771     A = psImageAlloc(polyOrder, polyOrder, PS_TYPE_F64);
    772     ALUD = psImageAlloc(polyOrder, polyOrder, PS_TYPE_F64);
    773     B = psVectorAlloc(polyOrder, PS_TYPE_F64);
    774     coeffs = psVectorAlloc(polyOrder, PS_TYPE_F64);
    775     X = psVectorAlloc(x->n, PS_TYPE_F64);
     758    A       = psImageAlloc(polyOrder, polyOrder, PS_TYPE_F64);
     759    ALUD    = psImageAlloc(polyOrder, polyOrder, PS_TYPE_F64);
     760    B       = psVectorAlloc(polyOrder, PS_TYPE_F64);
     761    coeffs  = psVectorAlloc(polyOrder, PS_TYPE_F64);
     762    X       = psVectorAlloc(x->n, PS_TYPE_F64);
    776763    outPerm = psVectorAlloc(polyOrder, PS_TYPE_F64);
     764    xSums   = psVectorAlloc(1+2*polyOrder, PS_TYPE_F64);
    777765
    778766    // Initialize data structures.
    … …  
    791779
    792780    // Build the B and A data structs.
    793     for(k=0;k<(polyOrder);k++) {
    794         for (i=0;i<X->n;i++) {
    795             B->data.F64[k]+= y->data.F64[i] * (pow(X->data.F64[i], (float) k));
    796         }
    797     }
    798 
    799     for(k=0;k<(polyOrder);k++) {
    800         for(j=0;j<(polyOrder);j++) {
    801             for (i=0;i<X->n;i++) {
    802                 A->data.F64[k][j]+= pow(X->data.F64[i], (float) (k + j));
     781    for (i=0;i<X->n;i++) {
     782        p_psBuildSums1D(X->data.F64[i], 2*polyOrder, xSums);
     783
     784        for(k=0;k<(polyOrder);k++) {
     785            B->data.F64[k]+= y->data.F64[i] * xSums->data.F64[k];
     786        }
     787
     788        for(k=0;k<(polyOrder);k++) {
     789            for(j=0;j<(polyOrder);j++) {
     790                A->data.F64[k][j]+= xSums->data.F64[k+j];
    803791            }
    804792        }
    … …  
    818806    psFree(X);
    819807    psFree(outPerm);
     808    psFree(xSums);
    820809
    821810    return(NULL);
Note: See TracChangeset for help on using the changeset viewer.