IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 7535


Ignore:
Timestamp:
Jun 13, 2006, 10:02:05 AM (20 years ago)
Author:
Paul Price
Message:

x refers to the column, y refers to the row (bug 765).

File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/psLib/src/imageops/psImageStats.c

    r6959 r7535  
    99 *  @author GLG, MHPCC
    1010 *
    11  *  @version $Revision: 1.95 $ $Name: not supported by cvs2svn $
    12  *  @date $Date: 2006-04-22 20:52:30 $
     11 *  @version $Revision: 1.96 $ $Name: not supported by cvs2svn $
     12 *  @date $Date: 2006-06-13 20:02:05 $
    1313 *
    1414 *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
    … …  
    225225    chebPolys[i][j]
    226226    sums[i][j]: This will contain the sum of
    227                 input->data.F32[x][y] *
     227                input->data.F32[y][x] *
    228228                psPolynomial1DEval(
    229229chebPolys[i],
    … …  
    289289            //            for (x = 0; x < input->numRows; x++) {
    290290            //                for (y = 0; y < input->numCols; y++) {
    291             for (x = input->row0; x < (input->row0 + input->numRows); x++) {
    292                 for (y = input->col0; y < (input->col0 + input->numCols); y++) {
     291            for (y = input->row0; y < (input->row0 + input->numRows); y++) {
     292                for (x = input->col0; x < (input->col0 + input->numCols); x++) {
    293293                    double pixel = 0.0;
    294294                    if (input->type.type == PS_TYPE_S8) {
    295                         pixel = (double) input->data.S8[x][y];
     295                        pixel = (double) input->data.S8[y][x];
    296296                    } else if (input->type.type == PS_TYPE_U16) {
    297                         pixel = (double) input->data.U16[x][y];
     297                        pixel = (double) input->data.U16[y][x];
    298298                    } else if (input->type.type == PS_TYPE_F32) {
    299                         pixel = (double) input->data.F32[x][y];
     299                        pixel = (double) input->data.F32[y][x];
    300300                    } else if (input->type.type == PS_TYPE_F64) {
    301                         pixel = input->data.F64[x][y];
     301                        pixel = input->data.F64[y][x];
    302302                    }
    303303                    sums[i][j] += pixel * psPolynomial1DEval(chebPolys[i],rScalingFactors[x]) *
    … …  
    341341}
    342342
    343 
    344 
    345 
    346 /*****************************************************************************
    347 psImageFitPolynomial(): This routine takes as input a 2-D image and produces
    348 as output the coefficients of the Chebyshev polynomials which match that input
    349 image.  This is a TEST version of the code.  It is not used by anything.
    350   Input:
    351   Output:
    352   Internal Data Structures:
    353     chebPolys[i][j]
    354     sums[i][j]: This will contain the sum of
    355                 input->data.F32[x][y] *
    356                 psPolynomial1DEval(
    357 chebPolys[i],
    358 (float) x) *
    359                 psPolynomial1DEval(
    360 chebPolys[j],
    361 (float) y,
    362 );
    363         over all pixels (x,y) in the image.
    364   *****************************************************************************/
    365 psPolynomial2D* psImageFitPolynomialTest(psPolynomial2D* coeffs,
    366         const psImage* input)
    367 {
    368     PS_ASSERT_IMAGE_NON_NULL(input, NULL);
    369     PS_ASSERT_IMAGE_NON_EMPTY(input, NULL);
    370     if ((input->type.type != PS_TYPE_S8) &&
    371             (input->type.type != PS_TYPE_U16) &&
    372             (input->type.type != PS_TYPE_F32) &&
    373             (input->type.type != PS_TYPE_F64)) {
    374         psError(PS_ERR_BAD_PARAMETER_TYPE, true, "Unallowable image type.\n");
    375     }
    376     PS_ASSERT_POLY_NON_NULL(coeffs, NULL);
    377     PS_ASSERT_POLY_TYPE(coeffs, PS_POLYNOMIAL_CHEB, NULL);
    378     psS32 x = 0;
    379     psS32 y = 0;
    380     psS32 i = 0;
    381     psS32 j = 0;
    382     double **sums = NULL;
    383     psPolynomial1D* *chebPolys = NULL;
    384     psS32 maxChebyPoly = 0;
    385     double *cScalingFactors = NULL;
    386     double *rScalingFactors = NULL;
    387     psImage *nodes = psImageAlloc(input->numCols, input->numRows, PS_TYPE_F64);
    388 
    389     double min = -1.0;
    390     double max = 1.0;
    391     double bma = 0.5 * (max-min);  // 1
    392     double bpa = 0.5 * (max+min);  // 0
    393     // We must calculate the value of the image at the nodes where the
    394     // Chebyshev polynomials are 0.
    395     for (x = 0; x < input->numRows; x++) {
    396         double xTmp = cos(M_PI * (0.5 + ((float) x)) / ((float) input->numRows));
    397         double xNode = - ((xTmp + bma + bpa) - 1.0);
    398         double xOrig = ((float) input->numRows) * (xNode - min) / (max - min);
    399 
    400         for (y = 0; y < input->numCols; y++) {
    401             double yTmp = cos(M_PI * (0.5 + ((float) y)) / ((float) input->numCols));
    402             double yNode = - ((yTmp + bma + bpa) - 1.0);
    403             double yOrig = ((float) input->numCols) * (yNode - min) / (max - min);
    404 
    405             //            nodes->data.F64[x][y] = psImagePixelInterpolate(input, yNode, xNode, NULL, 0, 0.0, PS_INTERPOLATE_BILINEAR);
    406             //            nodes->data.F64[x][y] = psImagePixelInterpolate(input, yTmp, xTmp, NULL, 0, 0.0, PS_INTERPOLATE_BILINEAR);
    407             nodes->data.F64[x][y] = psImagePixelInterpolate(input, yOrig, xOrig, NULL,
    408                                     0, 0.0, PS_INTERPOLATE_BILINEAR);
    409         }
    410     }
    411 
    412     // Create the sums[][] data structure.  This
    413     // will hold the LHS of
    414     // equation
    415     // 29 in the ADD: sums[k][l] = SUM {
    416     // image(x,y) * Tk(x) * Tl(y) }
    417     sums = (double **)psAlloc((1 + coeffs->nX) * sizeof(double *));
    418     for (i = 0; i < (1 + coeffs->nX); i++) {
    419         sums[i] = (double *)psAlloc((1 + coeffs->nY) * sizeof(double));
    420     }
    421     // We scale the pixel positions to values
    422     // between -1.0 and 1.0
    423     rScalingFactors = calcScaleFactors(input->numRows);
    424     cScalingFactors = calcScaleFactors(input->numCols);
    425 
    426     // Determine how many Chebyshev polynomials
    427     // are needed, then create them.
    428     maxChebyPoly = coeffs->nX;
    429     if (coeffs->nY > coeffs->nX) {
    430         maxChebyPoly = coeffs->nY;
    431     }
    432     chebPolys = p_psCreateChebyshevPolys(maxChebyPoly + 1);
    433 
    434     // Compute the sums[][] data structure.
    435     for (i = 0; i < (1 + coeffs->nX); i++) {
    436         for (j = 0; j < (1 + coeffs->nY); j++) {
    437             sums[i][j] = 0.0;
    438             for (x = 0; x < input->numRows; x++) {
    439                 for (y = 0; y < input->numCols; y++) {
    440                     double pixel;
    441                     if (0) {
    442                         if (input->type.type == PS_TYPE_S8) {
    443                             pixel = (double) input->data.S8[x][y];
    444                         } else if (input->type.type == PS_TYPE_U16) {
    445                             pixel = (double) input->data.U16[x][y];
    446                         } else if (input->type.type == PS_TYPE_F32) {
    447                             pixel = (double) input->data.F32[x][y];
    448                         } else if (input->type.type == PS_TYPE_F64) {
    449                             pixel = input->data.F64[x][y];
    450                         }
    451                     }
    452                     pixel = nodes->data.F64[x][y];
    453                     sums[i][j] += pixel * psPolynomial1DEval(chebPolys[i], rScalingFactors[x]) *
    454                                   psPolynomial1DEval(chebPolys[j], cScalingFactors[y]);
    455                 }
    456             }
    457         }
    458     }
    459 
    460     for (i = 0; i < (1 + coeffs->nX); i++) {
    461         for (j = 0; j < (1 + coeffs->nY); j++) {
    462             coeffs->coeff[i][j] = sums[i][j];
    463             coeffs->coeff[i][j] /= (double)(input->numRows * input->numCols);
    464 
    465             if ((i != 0) && (j != 0)) {
    466                 coeffs->coeff[i][j] *= 4.0;
    467             } else if ((i == 0) && (j == 0)) {
    468                 coeffs->coeff[i][j] *= 1.0;
    469             } else {
    470                 coeffs->coeff[i][j] *= 2.0;
    471             }
    472         }
    473     }
    474 
    475     // Free the Chebyshev polynomials that were
    476     // created in this routine.
    477     for (i = 0; i < maxChebyPoly + 1; i++) {
    478         psFree(chebPolys[i]);
    479     }
    480     psFree(chebPolys);
    481 
    482     // Free some data
    483     for (i = 0; i < (1 + coeffs->nX); i++) {
    484         psFree(sums[i]);
    485     }
    486     psFree(sums);
    487     psFree(cScalingFactors);
    488     psFree(rScalingFactors);
    489     psFree(nodes);
    490 
    491     return (coeffs);
    492 }
    493 
    494343/*****************************************************************************
    495344XXX: Use static variables for Chebyshev polynomials and scaling factors.
    … …  
    500349    PS_ASSERT_POLY_TYPE(coeffs, PS_POLYNOMIAL_CHEB, NULL);
    501350
    502     psS32 x = 0;
    503     psS32 y = 0;
    504     psS32 i = 0;
    505     psS32 j = 0;
    506351    psPolynomial1D* *chebPolys = NULL;
    507     psS32 maxChebyPoly = 0;
     352    long maxChebyPoly = 0;
    508353    double *cScalingFactors = NULL;
    509354    double *rScalingFactors = NULL;
    510     float polySum = 0.0;
    511355
    512356    // We scale the pixel positions to values between -1.0 and 1.0
    … …  
    526370    chebPolys = p_psCreateChebyshevPolys(maxChebyPoly + 1);
    527371
    528     //    for (x = 0; x < input->numRows; x++) {
    529     //        for (y = 0; y < input->numCols; y++) {
    530     for (x = input->row0; x < (input->row0 + input->numRows); x++) {
    531         for (y = input->col0; y < (input->col0 + input->numCols); y++) {
    532             polySum = 0.0;
    533             for (i = 0; i < (1 + coeffs->nX); i++) {
    534                 for (j = 0; j < (1 + coeffs->nY); j++) {
     372    for (long y = input->row0; y < (input->row0 + input->numRows); y++) {
     373        for (long x = input->col0; x < (input->col0 + input->numCols); x++) {
     374            double polySum = 0.0;
     375            for (long i = 0; i < (1 + coeffs->nX); i++) {
     376                for (long j = 0; j < (1 + coeffs->nY); j++) {
    535377                    polySum +=
    536378                        psPolynomial1DEval(chebPolys[i], rScalingFactors[x]) *
    … …  
    541383
    542384            if (input->type.type == PS_TYPE_S8) {
    543                 input->data.S8[x][y] = (char) polySum;
     385                input->data.S8[y][x] = (char) polySum;
    544386            } else if (input->type.type == PS_TYPE_U16) {
    545                 input->data.U16[x][y] = (short int) polySum;
     387                input->data.U16[y][x] = (short int) polySum;
    546388            } else if (input->type.type == PS_TYPE_F32) {
    547                 input->data.F32[x][y] = (float) polySum;
     389                input->data.F32[y][x] = (float) polySum;
    548390            } else if (input->type.type == PS_TYPE_F64) {
    549                 input->data.F64[x][y] = polySum;
     391                input->data.F64[y][x] = polySum;
    550392            }
    551393        }
    … …  
    555397    // created in this routine.
    556398    // XXX: Use static data structures here.
    557     for (i = 0; i < maxChebyPoly + 1; i++) {
     399    for (long i = 0; i < maxChebyPoly + 1; i++) {
    558400        psFree(chebPolys[i]);
    559401    }
    … …  
    750592        return -1;
    751593    }/* else if (col0 == col1 && row0 == row1) {
    752                             psError(PS_ERR_BAD_PARAMETER_VALUE, true,
    753                                     "Invalid psRegion specified.  Region contains only 1 pixel.\n");
    754                             return -1;
    755                         }
    756                     */
     594                                psError(PS_ERR_BAD_PARAMETER_VALUE, true,
     595                                        "Invalid psRegion specified.  Region contains only 1 pixel.\n");
     596                                return -1;
     597                            }
     598                        */
    757599    x0 = col0;
    758600    x1 = col1;
Note: See TracChangeset for help on using the changeset viewer.