Index: /trunk/psLib/src/astro/psCoord.c =================================================================== --- /trunk/psLib/src/astro/psCoord.c (revision 15253) +++ /trunk/psLib/src/astro/psCoord.c (revision 15254) @@ -10,6 +10,6 @@ * @author GLG, MHPCC * - * @version $Revision: 1.140 $ $Name: not supported by cvs2svn $ - * @date $Date: 2007-09-28 00:37:27 $ + * @version $Revision: 1.141 $ $Name: not supported by cvs2svn $ + * @date $Date: 2007-10-09 19:25:44 $ * * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii @@ -765,7 +765,7 @@ for (psS32 j = 0 ; j < orderY+1 ; j++) { myPT->x->coeff[i][j] = 0.0; - myPT->x->mask[i][j] = 0; myPT->y->coeff[i][j] = 0.0; - myPT->y->mask[i][j] = 0; + myPT->x->coeffMask[i][j] = PS_POLY_MASK_NONE; + myPT->y->coeffMask[i][j] = PS_POLY_MASK_NONE; } } @@ -807,24 +807,26 @@ for (psS32 t2y = 0 ; t2y < (trans2->x->nY + 1) ; t2y++) { psTrace("psLib.astro", 6, "X: -------------------- (t2x, t2y) (%d, %d) --------------------\n", t2x, t2y); - if (trans2->x->mask[t2x][t2y] == 0) { - psTrace("psLib.astro", 6, "In this iteration, we raise trans1->x to the %d power and trans1->y to the %d-power.\n", t2x, t2y); - psPolynomial2D *newPoly = multiplyDPoly2D(trans1XPolys[t2x], trans1YPolys[t2y]); - - if (psTraceGetLevel("psLib.astro") >= 6) { - PS_POLY_PRINT_2D(newPoly); - } - - // Set the appropriate coeffs in myPT->x - for (psS32 i = 0 ; i < (1 + newPoly->nX) ; i++) { - for (psS32 j = 0 ; j < (1 + newPoly->nY) ; j++) { - myPT->x->coeff[i][j]+= newPoly->coeff[i][j] * trans2->x->coeff[t2x][t2y]; - } - } - - if (psTraceGetLevel("psLib.astro") >= 6) { - PS_POLY_PRINT_2D(myPT->x); - } - psFree(newPoly); - } + if (trans2->x->coeffMask[t2x][t2y] & PS_POLY_MASK_SET) { + continue; + } + + psTrace("psLib.astro", 6, "In this iteration, we raise trans1->x to the %d power and trans1->y to the %d-power.\n", t2x, t2y); + psPolynomial2D *newPoly = multiplyDPoly2D(trans1XPolys[t2x], trans1YPolys[t2y]); + + if (psTraceGetLevel("psLib.astro") >= 6) { + PS_POLY_PRINT_2D(newPoly); + } + + // Set the appropriate coeffs in myPT->x + for (psS32 i = 0 ; i < (1 + newPoly->nX) ; i++) { + for (psS32 j = 0 ; j < (1 + newPoly->nY) ; j++) { + myPT->x->coeff[i][j]+= newPoly->coeff[i][j] * trans2->x->coeff[t2x][t2y]; + } + } + + if (psTraceGetLevel("psLib.astro") >= 6) { + PS_POLY_PRINT_2D(myPT->x); + } + psFree(newPoly); } } @@ -841,23 +843,24 @@ for (psS32 t2y = 0 ; t2y < (trans2->y->nY + 1) ; t2y++) { psTrace("psLib.astro", 5, "Y: -------------------- (t2x, t2y) (%d, %d) --------------------\n", t2x, t2y); - if (trans2->y->mask[t2x][t2y] == 0) { - psTrace("psLib.astro", 5, "In this iteration, we raise trans1->x to the %d power and trans1->y to the %d-power.\n", t2x, t2y); - psPolynomial2D *newPoly = multiplyDPoly2D(trans1XPolys[t2x], trans1YPolys[t2y]); - - if (psTraceGetLevel("psLib.astro") >= 6) { - PS_POLY_PRINT_2D(newPoly); - } - - // Set the appropriate coeffs in myPT->x - for (psS32 i = 0 ; i < (1 + newPoly->nX) ; i++) { - for (psS32 j = 0 ; j < (1 + newPoly->nY) ; j++) { - myPT->y->coeff[i][j]+= newPoly->coeff[i][j] * trans2->y->coeff[t2x][t2y]; - } - } - if (psTraceGetLevel("psLib.astro") >= 6) { - PS_POLY_PRINT_2D(myPT->x); - } - psFree(newPoly); - } + if (trans2->y->coeffMask[t2x][t2y] & PS_POLY_MASK_SET) { + continue; + } + psTrace("psLib.astro", 5, "In this iteration, we raise trans1->x to the %d power and trans1->y to the %d-power.\n", t2x, t2y); + psPolynomial2D *newPoly = multiplyDPoly2D(trans1XPolys[t2x], trans1YPolys[t2y]); + + if (psTraceGetLevel("psLib.astro") >= 6) { + PS_POLY_PRINT_2D(newPoly); + } + + // Set the appropriate coeffs in myPT->x + for (psS32 i = 0 ; i < (1 + newPoly->nX) ; i++) { + for (psS32 j = 0 ; j < (1 + newPoly->nY) ; j++) { + myPT->y->coeff[i][j]+= newPoly->coeff[i][j] * trans2->y->coeff[t2x][t2y]; + } + } + if (psTraceGetLevel("psLib.astro") >= 6) { + PS_POLY_PRINT_2D(myPT->x); + } + psFree(newPoly); } } Index: /trunk/psLib/src/imageops/psImagePixelInterpolate.c =================================================================== --- /trunk/psLib/src/imageops/psImagePixelInterpolate.c (revision 15253) +++ /trunk/psLib/src/imageops/psImagePixelInterpolate.c (revision 15254) @@ -11,6 +11,6 @@ * @author Eugene Magnier, IfA * - * @version $Revision: 1.2 $ $Name: not supported by cvs2svn $ - * @date $Date: 2007-09-20 23:54:25 $ + * @version $Revision: 1.3 $ $Name: not supported by cvs2svn $ + * @date $Date: 2007-10-09 19:25:44 $ * * Copyright 2007 Institute for Astronomy, University of Hawaii @@ -153,7 +153,7 @@ // allocate a 2D polynomial to fit a quadratic to the valid neighbor pixels. psPolynomial2D *poly = psPolynomial2DAlloc (PS_POLYNOMIAL_ORD, 2, 2); - poly->mask[2][2] = 1; - poly->mask[2][1] = 1; - poly->mask[1][2] = 1; + poly->coeffMask[2][2] = PS_POLY_MASK_SET; + poly->coeffMask[2][1] = PS_POLY_MASK_SET; + poly->coeffMask[1][2] = PS_POLY_MASK_SET; for (int iy = 0; iy < state->numRows; iy++) { @@ -275,11 +275,11 @@ // allocate a 2D polynomial to fit a quadratic to the valid neighbor pixels. psPolynomial2D *poly2o = psPolynomial2DAlloc (PS_POLYNOMIAL_ORD, 2, 2); - poly2o->mask[2][2] = 1; - poly2o->mask[2][1] = 1; - poly2o->mask[1][2] = 1; + poly2o->coeffMask[2][2] = PS_POLY_MASK_SET; + poly2o->coeffMask[2][1] = PS_POLY_MASK_SET; + poly2o->coeffMask[1][2] = PS_POLY_MASK_SET; // allocate a 2D polynomial to fit a plane to the valid neighbor pixels. psPolynomial2D *poly1o = psPolynomial2DAlloc (PS_POLYNOMIAL_ORD, 1, 1); - poly2o->mask[1][1] = 1; + poly2o->coeffMask[1][1] = PS_POLY_MASK_SET; for (int iy = 0; iy < state->numRows; iy++) { Index: /trunk/psLib/src/math/psMinimizePolyFit.c =================================================================== --- /trunk/psLib/src/math/psMinimizePolyFit.c (revision 15253) +++ /trunk/psLib/src/math/psMinimizePolyFit.c (revision 15254) @@ -10,6 +10,6 @@ * @author EAM, IfA * - * @version $Revision: 1.31 $ $Name: not supported by cvs2svn $ - * @date $Date: 2007-05-16 16:16:48 $ + * @version $Revision: 1.32 $ $Name: not supported by cvs2svn $ + * @date $Date: 2007-10-09 19:25:45 $ * * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii @@ -340,5 +340,5 @@ psImage *sums = psImageAlloc(numData, numTerms, PS_TYPE_F64); for (int i = 0; i < numTerms; i++) { - if (myPoly->mask[i]) { + if (myPoly->coeffMask[i] & PS_POLY_MASK_BOTH) { continue; } @@ -358,5 +358,5 @@ dataMask = mask->data.U8; } - psU8 *termMask = myPoly->mask; // Mask for polynomial terms + psU8 *coeffMask = myPoly->coeffMask; // Mask for polynomial terms psF64 *yData = y->data.F64; // Coordinate data psF64 *yErrData = NULL; // Errors in the coordinate @@ -380,5 +380,5 @@ for (int i = 0; i < numTerms; i++) { - if (termMask[i]) { + if (coeffMask[i] & PS_POLY_MASK_BOTH) { matrix[i][i] = 1.0; continue; @@ -387,5 +387,5 @@ matrix[i][i] += sumsData[i][k] * sumsData[i][k] * wt; // The diagonal entry for (int j = i + 1; j < numTerms; j++) { // The upper diagonal only: we will use symmetry - if (termMask[j]) { + if (coeffMask[j] & PS_POLY_MASK_BOTH) { continue; } @@ -533,5 +533,5 @@ dataMask = mask->data.U8; } - psU8 *termMask = myPoly->mask; // Dereferenced version of mask for polynomial terms + psU8 *coeffMask = myPoly->coeffMask; // Dereferenced version of mask for polynomial terms psF64 *ordinates = NULL; // Dereferenced version of ordinate data if (x) { @@ -567,5 +567,5 @@ for (int i = 0; i < nTerm; i++) { - if (termMask[i]) { + if (coeffMask[i] & PS_POLY_MASK_SET) { matrix[i][i] = 1.0; continue; @@ -574,5 +574,5 @@ matrix[i][i] += sums[2 * i] * wt; // The diagonal entry for (int j = i + 1; j < nTerm; j++) { // The upper diagonal only: we will use symmetry - if (termMask[j]) { + if (coeffMask[j] & PS_POLY_MASK_SET) { continue; } @@ -584,4 +584,32 @@ } psFree(xSums); + + // elements which are masked for fitting need to be subtracted from the vector + for (int i = 0; i < nTerm; i++) { + if (coeffMask[i] & PS_POLY_MASK_BOTH) { + continue; + } + for (int j = 0; j < nTerm; j++) { // The upper diagonal only: we will use symmetry + if (coeffMask[j] & PS_POLY_MASK_SET) { + continue; + } + if (!(coeffMask[j] & PS_POLY_MASK_FIT)) { + continue; + } + vector[i] -= matrix[i][j]*myPoly->coeff[j]; + } + } + + // set the un-fitted and un-set elements to 0 or 1 for pivots + for (int i = 0; i < nTerm; i++) { + if (coeffMask[i] & PS_POLY_MASK_BOTH) { + for (int j = 0; j < nTerm; j++) { // The upper diagonal only: we will use symmetry + matrix[i][j] = 0.0; + matrix[j][i] = 0.0; + } + matrix[i][i] = 1.0; + continue; + } + } if (psTraceGetLevel("psLib.math") >= 4) { @@ -610,9 +638,9 @@ // polynomial coefficients. this is only true for the 1D case for (psS32 k = 0; k < nTerm; k++) { - myPoly->coeff[k] = B->data.F64[k]; - myPoly->coeffErr[k] = sqrt(A->data.F64[k][k]); - } - } - + if (coeffMask[k] & PS_POLY_MASK_FIT) continue; + myPoly->coeff[k] = B->data.F64[k]; + myPoly->coeffErr[k] = sqrt(A->data.F64[k][k]); + } + } } else { // LUD version of the fit @@ -633,4 +661,5 @@ } else { for (psS32 k = 0; k < nTerm; k++) { + if (coeffMask[k] & PS_POLY_MASK_FIT) continue; myPoly->coeff[k] = coeffs->data.F64[k]; // XXX LUD does not give inverse of A @@ -990,5 +1019,5 @@ psFree(A); psFree(B); - psTrace("psLib.math", 4, "---- %s() End ----\n", __func__); + psTrace("psLib.math", 6, "---- %s() End ----\n", __func__); return false; } @@ -997,5 +1026,5 @@ psF64 **matrix = A->data.F64; // Dereference the least-squares matrix psF64 *vector = B->data.F64; // Dereference the least-squares vector - psU8 **termMask = myPoly->mask; // Dereference mask for polynomial terms + psU8 **coeffMask = myPoly->coeffMask; // Dereference mask for polynomial terms psU8 *dataMask = NULL; // Dereference mask for data if (mask) { @@ -1031,5 +1060,5 @@ int l = i / nYterm; // x index int m = i % nYterm; // y index - if (termMask[l][m]) { + if (coeffMask[l][m] & PS_POLY_MASK_SET) { matrix[i][i] = 1.0; continue; @@ -1040,5 +1069,5 @@ int p = j / nYterm; // x index int q = j % nYterm; // y index - if (termMask[p][q]) { + if (coeffMask[p][q] & PS_POLY_MASK_SET) { continue; } @@ -1050,4 +1079,53 @@ } psFree(xySums); + + // elements which are masked for fitting need to be subtracted from the vector + for (int i = 0; i < nTerm; i++) { + int ix = i / nYterm; // x index + int iy = i % nYterm; // y index + if (coeffMask[ix][iy] & PS_POLY_MASK_BOTH) { + continue; + } + for (int j = 0; j < nTerm; j++) { // The upper diagonal only: we will use symmetry + int jx = j / nYterm; // x index + int jy = j % nYterm; // y index + if (coeffMask[jx][jy] & PS_POLY_MASK_SET) { + continue; + } + if (!(coeffMask[jx][jy] & PS_POLY_MASK_FIT)) { + continue; + } + vector[i] -= matrix[i][j]*myPoly->coeff[jx][jy]; + } + } + + // set the un-fitted and un-set elements to 0 or 1 for pivots + for (int i = 0; i < nTerm; i++) { + int ix = i / nYterm; // x index + int iy = i % nYterm; // y index + if (coeffMask[ix][iy] & PS_POLY_MASK_BOTH) { + for (int j = 0; j < nTerm; j++) { // The upper diagonal only: we will use symmetry + matrix[i][j] = 0.0; + matrix[j][i] = 0.0; + } + matrix[i][i] = 1.0; + continue; + } + } + + if (psTraceGetLevel("psLib.math") >= 4) { + printf("Least-squares vector:\n"); + for (int i = 0; i < nTerm; i++) { + printf("%f ", B->data.F64[i]); + } + printf("\n"); + printf("Least-squares matrix:\n"); + for (int i = 0; i < nTerm; i++) { + for (int j = 0; j < nTerm; j++) { + printf("%f ", A->data.F64[i][j]); + } + printf("\n"); + } + } if (!psMatrixGJSolve(A, B)) { @@ -1062,11 +1140,14 @@ int l = i / nYterm; // x index int m = i % nYterm; // y index - myPoly->coeff[l][m] = B->data.F64[i]; - myPoly->coeffErr[l][m] = sqrt(A->data.F64[i][i]); + + // retain the incoming values if masked on the fit + if (coeffMask[l][m] & PS_POLY_MASK_FIT) continue; + myPoly->coeff[l][m] = B->data.F64[i]; + myPoly->coeffErr[l][m] = sqrt(A->data.F64[i][i]); } psFree(A); psFree(B); - psTrace("psLib.math", 4, "---- %s() end ----\n", __func__); + psTrace("psLib.math", 6, "---- %s() end ----\n", __func__); return true; } @@ -1387,5 +1468,5 @@ dataMask = mask->data.U8; } - psU8 ***termMask = myPoly->mask; // Mask for polynomial terms + psU8 ***coeffMask = myPoly->coeffMask; // Mask for polynomial terms int nYZterm = nYterm * nZterm; // Multiplication of the numbers, to calculate the index @@ -1411,5 +1492,5 @@ int iy = (i % nYZterm) / nZterm; // y index int iz = (i % nYZterm) % nZterm; // z index - if (termMask[ix][iy][iz]) { + if (coeffMask[ix][iy][iz] & PS_POLY_MASK_BOTH) { matrix[i][i] = 1.0; continue; @@ -1422,5 +1503,5 @@ int jy = (j % nYZterm) / nZterm; // y index int jz = (j % nYZterm) % nZterm; // z index - if (termMask[jx][jy][jz]) { + if (coeffMask[jx][jy][jz] & PS_POLY_MASK_BOTH) { continue; } @@ -1455,4 +1536,5 @@ int iy = (i % nYZterm) / nZterm; // y index int iz = (i % nYZterm) % nZterm; // z index + if (coeffMask[ix][iy][iz] & PS_POLY_MASK_FIT) continue; myPoly->coeff[ix][iy][iz] = B->data.F64[i]; myPoly->coeffErr[ix][iy][iz] = sqrt(A->data.F64[i][i]); @@ -1480,4 +1562,5 @@ for (psS32 iz = 0; iz < nZterm; iz++) { psS32 nx = ix+iy*nXterm+iz*nXterm*nYterm; + if (coeffMask[ix][iy][iz] & PS_POLY_MASK_FIT) continue; myPoly->coeff[ix][iy][iz] = coeffs->data.F64[nx]; // XXX myPoly->coeffErr[ix][iy][iz] = sqrt(A->data.F64[nx][nx]); @@ -1839,5 +1922,5 @@ dataMask = mask->data.U8; } - psU8 ****termMask = myPoly->mask; // Mask for polynomial terms + psU8 ****coeffMask = myPoly->coeffMask; // Mask for polynomial terms int nYZTterm = nYterm * nZterm * nTterm; // Multiplication of the numbers, for calculating the index int nZTterm = nZterm * nTterm; // Multiplication of the numbers, for calculating the index @@ -1865,5 +1948,5 @@ int iz = ((i % (nYZTterm)) % (nZTterm)) / nTterm; // z index int it = ((i % (nYZTterm)) % (nZTterm)) % nTterm; // t index - if (termMask[ix][iy][iz][it]) { + if (coeffMask[ix][iy][iz][it] & PS_POLY_MASK_BOTH) { matrix[i][i] = 1.0; continue; @@ -1877,5 +1960,5 @@ int jz = ((j % nYZTterm) % nZTterm) / nTterm; // z index int jt = ((j % nYZTterm) % nZTterm) % nTterm; // t index - if (termMask[jx][jy][jz][jt]) { + if (coeffMask[jx][jy][jz][jt] & PS_POLY_MASK_BOTH) { continue; } @@ -1918,4 +2001,5 @@ int iz = ((i % nYZTterm) % nZTterm) / nTterm; // z index int it = ((i % nYZTterm) % nZTterm) % nTterm; // t index + if (coeffMask[ix][iy][iz][it] & PS_POLY_MASK_FIT) continue; myPoly->coeff[ix][iy][iz][it] = B->data.F64[i]; myPoly->coeffErr[ix][iy][iz][it] = sqrt(A->data.F64[i][i]); @@ -1944,4 +2028,5 @@ for (psS32 it = 0; it < nTterm; it++) { psS32 nx = ix+iy*nXterm+iz*nXterm*nYterm+it*nXterm*nYterm*nZterm; + if (coeffMask[ix][iy][iz][it] & PS_POLY_MASK_FIT) continue; myPoly->coeff[ix][iy][iz][it] = coeffs->data.F64[nx]; // myPoly->coeffErr[ix][iy][iz][it] = sqrt(A->data.F64[nx][nx]); Index: /trunk/psLib/src/math/psPolynomialMetadata.c =================================================================== --- /trunk/psLib/src/math/psPolynomialMetadata.c (revision 15253) +++ /trunk/psLib/src/math/psPolynomialMetadata.c (revision 15254) @@ -12,6 +12,6 @@ * @author Ross Harman, MHPCC * - * @version $Revision: 1.1 $ $Name: not supported by cvs2svn $ - * @date $Date: 2007-02-06 21:55:28 $ + * @version $Revision: 1.2 $ $Name: not supported by cvs2svn $ + * @date $Date: 2007-10-09 19:25:45 $ * * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii @@ -62,9 +62,9 @@ // an undefined component implies the component was masked // this is symmetrical with the 1DtoMD function - poly->mask[nx] = 1; poly->coeff[nx] = 0; poly->coeffErr[nx] = 0; + poly->coeffMask[nx] = PS_POLY_MASK_SET; } else { - poly->mask[nx] = 0; + poly->coeffMask[nx] = PS_POLY_MASK_NONE; nElements ++; } @@ -123,5 +123,5 @@ // place polynomial entries on folder for (int nx = 0; nx < poly->nX + 1; nx++) { - if (poly->mask[nx] == 0) { + if (!(poly->coeffMask[nx] & PS_POLY_MASK_SET)) { sprintf(namespace, "VAL_X%02d", nx); sprintf(namespace_err, "ERR_X%02d", nx); @@ -174,9 +174,9 @@ // an undefined component implies the component was masked // this is symmetrical with the 2DtoMD function - poly->mask[nx][ny] = 1; poly->coeff[nx][ny] = 0; poly->coeffErr[nx][ny] = 0; + poly->coeffMask[nx][ny] = PS_POLY_MASK_SET; } else { - poly->mask[nx][ny] = 0; + poly->coeffMask[nx][ny] = PS_POLY_MASK_NONE; nElements ++; } @@ -203,14 +203,8 @@ PS_ASSERT_PTR_NON_NULL(md, false); PS_ASSERT_PTR_NON_NULL(poly, false); - //XXX: Current implementation only supports ordinary polynomials. + + // XXX Current implementation only supports ordinary polynomials. if (poly->type != PS_POLYNOMIAL_ORD) return false; - - // XXX I'm puzzled by this test. a polynomial of 0 order with a value of 0 is a - // perfectly valid polynomial and can be written out. a polynomial with all elements - // masked still carries information. - //Make sure polynomial isn't 0, completely empty - //if (poly->nX == 0 && poly->nY == 0 && poly->coeff[0][0] == 0) - //return false; int Nbyte; @@ -245,5 +239,5 @@ for (int nx = 0; nx < poly->nX + 1; nx++) { for (int ny = 0; ny < poly->nY + 1; ny++) { - if (poly->mask[nx][ny] == 0) { + if (!(poly->coeffMask[nx][ny] & PS_POLY_MASK_SET)) { sprintf(namespace, "VAL_X%02d_Y%02d", nx, ny); sprintf(namespace_err, "ERR_X%02d_Y%02d", nx, ny); @@ -304,9 +298,9 @@ // an undefined component implies the component was masked // this is symmetrical with the 3DtoMD function - poly->mask[nx][ny][nz] = 1; poly->coeff[nx][ny][nz] = 0; poly->coeffErr[nx][ny][nz] = 0; + poly->coeffMask[nx][ny][nz] = PS_POLY_MASK_SET; } else { - poly->mask[nx][ny][nz] = 0; + poly->coeffMask[nx][ny][nz] = PS_POLY_MASK_NONE; nElements ++; } @@ -370,5 +364,5 @@ for (int ny = 0; ny < poly->nY + 1; ny++) { for (int nz = 0; nz < poly->nZ + 1; nz++) { - if (poly->mask[nx][ny][nz] == 0) { + if (!(poly->coeffMask[nx][ny][nz] & PS_POLY_MASK_SET)) { sprintf(namespace, "VAL_X%02d_Y%02d_Z%02d", nx, ny, nz); sprintf(namespace_err, "ERR_X%02d_Y%02d_Z%02d", nx, ny, nz); @@ -437,9 +431,9 @@ // an undefined component implies the component was masked // this is symmetrical with the 4DtoMD function - poly->mask[nx][ny][nz][nt] = 1; poly->coeff[nx][ny][nz][nt] = 0; poly->coeffErr[nx][ny][nz][nt] = 0; + poly->coeffMask[nx][ny][nz][nt] = PS_POLY_MASK_SET; } else { - poly->mask[nx][ny][nz][nt] = 0; + poly->coeffMask[nx][ny][nz][nt] = PS_POLY_MASK_NONE; nElements ++; } @@ -506,5 +500,5 @@ for (int nz = 0; nz < poly->nZ + 1; nz++) { for (int nt = 0; nt < poly->nT + 1; nt++) { - if (poly->mask[nx][ny][nz][nt] == 0) { + if (!(poly->coeffMask[nx][ny][nz][nt] & PS_POLY_MASK_SET)) { sprintf(namespace, "VAL_X%02d_Y%02d_Z%02d_T%02d", nx, ny, nz, nt); sprintf(namespace_err, "ERR_X%02d_Y%02d_Z%02d_T%02d", nx, ny, nz, nt); Index: /trunk/psLib/src/math/psPolynomialUtils.c =================================================================== --- /trunk/psLib/src/math/psPolynomialUtils.c (revision 15253) +++ /trunk/psLib/src/math/psPolynomialUtils.c (revision 15254) @@ -159,7 +159,7 @@ psPolynomial2D *poly = psPolynomial2DAlloc (PS_POLYNOMIAL_ORD, 2, 2); - poly->mask[2][2] = 1; - poly->mask[1][2] = 1; - poly->mask[2][1] = 1; + poly->coeffMask[2][2] = PS_POLY_MASK_SET; + poly->coeffMask[1][2] = PS_POLY_MASK_SET; + poly->coeffMask[2][1] = PS_POLY_MASK_SET; poly->coeff[0][0] = Foo*(5.0/9.0) - (Fxp + Fxm)/3.0 - (Fyp + Fym)/3.0 ; @@ -212,5 +212,5 @@ for (int i = 0; i < nXout + 1; i++) { for (int j = 0; j < nYout + 1; j++) { - out->mask[i][j] = poly->mask[i+1][j]; + out->coeffMask[i][j] = poly->coeffMask[i+1][j]; out->coeff[i][j] = poly->coeff[i+1][j] * (i+1); } @@ -240,5 +240,5 @@ for (int i = 0; i < nXout + 1; i++) { for (int j = 0; j < nYout + 1; j++) { - out->mask[i][j] = poly->mask[i][j+1]; + out->coeffMask[i][j] = poly->coeffMask[i][j+1]; out->coeff[i][j] = poly->coeff[i][j+1] * (j+1); } Index: /trunk/psModules/src/astrom/pmAstrometryDistortion.c =================================================================== --- /trunk/psModules/src/astrom/pmAstrometryDistortion.c (revision 15253) +++ /trunk/psModules/src/astrom/pmAstrometryDistortion.c (revision 15254) @@ -7,6 +7,6 @@ * @author EAM, IfA * -* @version $Revision: 1.19 $ $Name: not supported by cvs2svn $ -* @date $Date: 2007-06-06 00:43:19 $ +* @version $Revision: 1.20 $ $Name: not supported by cvs2svn $ +* @date $Date: 2007-10-09 19:27:04 $ * * Copyright 2006 Institute for Astronomy, University of Hawaii @@ -62,5 +62,5 @@ psPolynomial2D *local = psPolynomial2DAlloc (PS_POLYNOMIAL_ORD, 1, 1); - local->mask[1][1] = 1; + local->coeffMask[1][1] = PS_POLY_MASK_SET; // measure gradient for fractional chip regions @@ -211,12 +211,10 @@ for (int i = 0; i <= fpa->toTPA->x->nX; i++) { for (int j = 0; j <= fpa->toTPA->x->nY; j++) { - if (fpa->toTPA->x->mask[i][j]) { - if ((i > 0) && (i <= fpa->toTPA->x->nX)) { - localX->mask[i-1][j] = 1; - } - if ((j > 0) && (j <= fpa->toTPA->x->nY)) { - localY->mask[i][j-1] = 1; - } - } + if ((i > 0) && (i <= fpa->toTPA->x->nX)) { + localX->coeffMask[i-1][j] = fpa->toTPA->x->coeffMask[i][j]; + } + if ((j > 0) && (j <= fpa->toTPA->x->nY)) { + localY->coeffMask[i][j-1] = fpa->toTPA->x->coeffMask[i][j]; + } } } @@ -240,18 +238,21 @@ fpa->toTPA->x->coeff[0][0] = 0; for (int i = 1; i <= fpa->toTPA->x->nX; i++) { - if (!fpa->toTPA->x->mask[i][0]) { - fpa->toTPA->x->coeff[i][0] = localX->coeff[i-1][0] / i; - } + if (fpa->toTPA->x->coeffMask[i][0] & PS_POLY_MASK_SET) { + continue; + } + fpa->toTPA->x->coeff[i][0] = localX->coeff[i-1][0] / i; } for (int j = 1; j <= fpa->toTPA->x->nY; j++) { - if (!fpa->toTPA->x->mask[0][j]) { - fpa->toTPA->x->coeff[0][j] = localY->coeff[0][j-1] / j; - } + if (fpa->toTPA->x->coeffMask[0][j] & PS_POLY_MASK_SET) { + continue; + } + fpa->toTPA->x->coeff[0][j] = localY->coeff[0][j-1] / j; } for (int i = 1; i <= fpa->toTPA->x->nX; i++) { for (int j = 1; j <= fpa->toTPA->x->nY; j++) { - if (!fpa->toTPA->x->mask[i][j]) { - fpa->toTPA->x->coeff[i][j] = 0.5*(localX->coeff[i-1][j] / i + localY->coeff[i][j-1] / j); - } + if (fpa->toTPA->x->coeffMask[i][j] & PS_POLY_MASK_SET) { + continue; + } + fpa->toTPA->x->coeff[i][j] = 0.5*(localX->coeff[i-1][j] / i + localY->coeff[i][j-1] / j); } } @@ -268,12 +269,10 @@ for (int i = 0; i < fpa->toTPA->y->nX; i++) { for (int j = 0; j < fpa->toTPA->y->nY; j++) { - if (fpa->toTPA->y->mask[i][j]) { - if ((i > 0) && (i <= fpa->toTPA->y->nX)) { - localX->mask[i-1][j] = 1; - } - if ((j > 0) && (j <= fpa->toTPA->y->nY)) { - localY->mask[i][j-1] = 1; - } - } + if ((i > 0) && (i <= fpa->toTPA->y->nX)) { + localX->coeffMask[i-1][j] = fpa->toTPA->y->coeffMask[i][j]; + } + if ((j > 0) && (j <= fpa->toTPA->y->nY)) { + localY->coeffMask[i][j-1] = fpa->toTPA->y->coeffMask[i][j]; + } } } @@ -286,18 +285,21 @@ fpa->toTPA->y->coeff[0][0] = 0; for (int i = 1; i <= fpa->toTPA->y->nX; i++) { - if (!fpa->toTPA->y->mask[i][0]) { - fpa->toTPA->y->coeff[i][0] = localX->coeff[i-1][0] / i; - } + if (fpa->toTPA->y->coeffMask[i][0] & PS_POLY_MASK_SET) { + continue; + } + fpa->toTPA->y->coeff[i][0] = localX->coeff[i-1][0] / i; } for (int j = 1; j <= fpa->toTPA->y->nY; j++) { - if (!fpa->toTPA->y->mask[0][j]) { - fpa->toTPA->y->coeff[0][j] = localY->coeff[0][j-1] / j; - } + if (fpa->toTPA->y->coeffMask[0][j] & PS_POLY_MASK_SET) { + continue; + } + fpa->toTPA->y->coeff[0][j] = localY->coeff[0][j-1] / j; } for (int i = 1; i <= fpa->toTPA->y->nX; i++) { for (int j = 1; j <= fpa->toTPA->y->nY; j++) { - if (!fpa->toTPA->y->mask[i][j]) { - fpa->toTPA->y->coeff[i][j] = 0.5*(localX->coeff[i-1][j] / i + localY->coeff[i][j-1] / j); + if (fpa->toTPA->y->coeffMask[i][j] & PS_POLY_MASK_SET) { + continue; } + fpa->toTPA->y->coeff[i][j] = 0.5*(localX->coeff[i-1][j] / i + localY->coeff[i][j-1] / j); } } Index: /trunk/psModules/src/astrom/pmAstrometryUtils.c =================================================================== --- /trunk/psModules/src/astrom/pmAstrometryUtils.c (revision 15253) +++ /trunk/psModules/src/astrom/pmAstrometryUtils.c (revision 15254) @@ -7,6 +7,6 @@ * @author EAM, IfA * - * @version $Revision: 1.4 $ $Name: not supported by cvs2svn $ - * @date $Date: 2007-03-30 21:12:56 $ + * @version $Revision: 1.5 $ $Name: not supported by cvs2svn $ + * @date $Date: 2007-10-09 19:27:04 $ * * Copyright 2006 Institute for Astronomy, University of Hawaii @@ -105,8 +105,8 @@ psPolynomial2D *yPy = psPolynomial2DCopy (NULL, yPx); for (int j = 0; j <= input->x->nY; j++) { - output->x->mask[i][j] = input->x->mask[i][j]; - output->y->mask[i][j] = input->y->mask[i][j]; - output->x->coeff[i][j] = (output->x->mask[i][j]) ? 0 : psPolynomial2DEval (xPy, Xo, Yo) / tgamma(i+1) / tgamma(j+1); - output->y->coeff[i][j] = (output->y->mask[i][j]) ? 0 : psPolynomial2DEval (yPy, Xo, Yo) / tgamma(i+1) / tgamma(j+1); + output->x->coeffMask[i][j] = input->x->coeffMask[i][j]; + output->y->coeffMask[i][j] = input->y->coeffMask[i][j]; + output->x->coeff[i][j] = (output->x->coeffMask[i][j] & PS_POLY_MASK_SET) ? 0 : psPolynomial2DEval (xPy, Xo, Yo) / tgamma(i+1) / tgamma(j+1); + output->y->coeff[i][j] = (output->y->coeffMask[i][j] & PS_POLY_MASK_SET) ? 0 : psPolynomial2DEval (yPy, Xo, Yo) / tgamma(i+1) / tgamma(j+1); // take the next derivative wrt y, catch output (is NULL on last pass) @@ -146,6 +146,6 @@ for (int j = 0; j <= order; j++) { if (i + j > order) { - transform->x->mask [i][j] = 1; - transform->y->mask [i][j] = 1; + transform->x->coeffMask [i][j] = PS_POLY_MASK_SET; + transform->y->coeffMask [i][j] = PS_POLY_MASK_SET; } } @@ -153,6 +153,6 @@ transform->x->coeff[1][0] = 1; transform->y->coeff[0][1] = 1; - transform->x->mask[1][0] = 0; - transform->y->mask[0][1] = 0; + transform->x->coeffMask[1][0] = PS_POLY_MASK_NONE; + transform->y->coeffMask[0][1] = PS_POLY_MASK_NONE; return transform; @@ -197,7 +197,7 @@ if (i + j > order) { // high-order cross terms must be masked (eg, x^3 y^2) - status &= transform->x->mask[i][j]; + status &= (transform->x->coeffMask[i][j] & PS_POLY_MASK_SET); } else { - status &= !transform->x->mask[i][j]; + status &= !(transform->x->coeffMask[i][j] & PS_POLY_MASK_SET); if ((i == 1) && (i + j == 1)) { // linear, diagonal terms must be non-zero @@ -216,7 +216,7 @@ if (i + j > order) { // high-order cross terms must be masked (eg, x^3 y^2) - status &= transform->y->mask[i][j]; + status &= (transform->y->coeffMask[i][j] & PS_POLY_MASK_SET); } else { - status &= !transform->y->mask[i][j]; + status &= !(transform->y->coeffMask[i][j] & PS_POLY_MASK_SET); if ((j == 1) && (i + j == 1)) { // linear, diagonal terms must be 1.0 @@ -249,6 +249,6 @@ for (int j = 0; j <= order; j++) { if (i + j > order) { - distort->x->mask [i][j][0][0] = 1; - distort->y->mask [i][j][0][0] = 1; + distort->x->coeffMask [i][j][0][0] = PS_POLY_MASK_SET; + distort->y->coeffMask [i][j][0][0] = PS_POLY_MASK_SET; } } @@ -306,7 +306,7 @@ if (i + j > order) { // high-order cross terms must be masked (eg, x^3 y^2) - status &= distort->x->mask[i][j][0][0]; + status &= (distort->x->coeffMask[i][j][0][0] & PS_POLY_MASK_SET); } else { - status &= !distort->x->mask[i][j][0][0]; + status &= !(distort->x->coeffMask[i][j][0][0] & PS_POLY_MASK_SET); if ((i == 1) && (i + j == 1)) { // linear, diagonal terms must be 1.0 @@ -325,7 +325,7 @@ if (i + j > order) { // high-order cross terms must be masked (eg, x^3 y^2) - status &= distort->y->mask[i][j][0][0]; + status &= (distort->y->coeffMask[i][j][0][0] & PS_POLY_MASK_SET); } else { - status &= !distort->y->mask[i][j][0][0]; + status &= !(distort->y->coeffMask[i][j][0][0] & PS_POLY_MASK_SET); if ((j == 1) && (i + j == 1)) { // linear, diagonal terms must be 1.0 Index: /trunk/psModules/src/astrom/pmAstrometryWCS.c =================================================================== --- /trunk/psModules/src/astrom/pmAstrometryWCS.c (revision 15253) +++ /trunk/psModules/src/astrom/pmAstrometryWCS.c (revision 15254) @@ -7,6 +7,6 @@ * @author EAM, IfA * - * @version $Revision: 1.22 $ $Name: not supported by cvs2svn $ - * @date $Date: 2007-04-17 19:39:04 $ + * @version $Revision: 1.23 $ $Name: not supported by cvs2svn $ + * @date $Date: 2007-10-09 19:27:04 $ * * Copyright 2006 Institute for Astronomy, University of Hawaii @@ -298,6 +298,6 @@ continue; if (i + j > fitOrder) { - wcs->trans->x->mask[i][j] = 1; - wcs->trans->y->mask[i][j] = 1; + wcs->trans->x->coeffMask[i][j] = PS_POLY_MASK_SET; + wcs->trans->y->coeffMask[i][j] = PS_POLY_MASK_SET; continue; } Index: /trunk/psModules/src/detrend/pmShutterCorrection.c =================================================================== --- /trunk/psModules/src/detrend/pmShutterCorrection.c (revision 15253) +++ /trunk/psModules/src/detrend/pmShutterCorrection.c (revision 15254) @@ -230,6 +230,6 @@ // mask out the terms we will not fit - line->mask[0][0] = 1; - line->mask[1][1] = 1; + line->coeffMask[0][0] = PS_POLY_MASK_SET; + line->coeffMask[1][1] = PS_POLY_MASK_SET; line->coeff[0][0] = 0; line->coeff[1][1] = 0; Index: /trunk/psModules/src/objects/pmPSF_IO.c =================================================================== --- /trunk/psModules/src/objects/pmPSF_IO.c (revision 15253) +++ /trunk/psModules/src/objects/pmPSF_IO.c (revision 15254) @@ -6,6 +6,6 @@ * @author EAM, IfA * - * @version $Revision: 1.26 $ $Name: not supported by cvs2svn $ - * @date $Date: 2007-10-06 01:02:49 $ + * @version $Revision: 1.27 $ $Name: not supported by cvs2svn $ + * @date $Date: 2007-10-09 19:26:25 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -358,5 +358,5 @@ psMetadataAddF32 (row, PS_LIST_TAIL, "VALUE", 0, "", poly->coeff[ix][iy]); psMetadataAddF32 (row, PS_LIST_TAIL, "ERROR", 0, "", poly->coeffErr[ix][iy]); - psMetadataAddU8 (row, PS_LIST_TAIL, "MASK", 0, "", poly->mask[ix][iy]); + psMetadataAddU8 (row, PS_LIST_TAIL, "MASK", 0, "", poly->coeffMask[ix][iy]); psArrayAdd (psfTable, 100, row); @@ -710,7 +710,7 @@ assert (xPow <= poly->nX); assert (yPow <= poly->nY); - poly->coeff[xPow][yPow] = psMetadataLookupF32 (&status, row, "VALUE"); - poly->coeffErr[xPow][yPow] = psMetadataLookupF32 (&status, row, "ERROR"); - poly->mask[xPow][yPow] = psMetadataLookupU8 (&status, row, "MASK"); + poly->coeff[xPow][yPow] = psMetadataLookupF32 (&status, row, "VALUE"); + poly->coeffErr[xPow][yPow] = psMetadataLookupF32 (&status, row, "ERROR"); + poly->coeffMask[xPow][yPow] = psMetadataLookupU8 (&status, row, "MASK"); } } Index: /trunk/psModules/src/objects/pmPSFtry.c =================================================================== --- /trunk/psModules/src/objects/pmPSFtry.c (revision 15253) +++ /trunk/psModules/src/objects/pmPSFtry.c (revision 15254) @@ -5,6 +5,6 @@ * @author EAM, IfA * - * @version $Revision: 1.48 $ $Name: not supported by cvs2svn $ - * @date $Date: 2007-09-25 22:05:05 $ + * @version $Revision: 1.49 $ $Name: not supported by cvs2svn $ + * @date $Date: 2007-10-09 19:26:25 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -280,5 +280,5 @@ // linear clipped fit of ApResid to r2rflux psPolynomial1D *poly = psPolynomial1DAlloc (PS_POLYNOMIAL_ORD, 1); - poly->mask[1] = 1; // fit only a constant offset (no SKYBIAS) + poly->coeffMask[1] = PS_POLY_MASK_SET; // fit only a constant offset (no SKYBIAS) // XXX replace this with a psVectorStats call? since we are not fitting the trend Index: /trunk/psModules/src/objects/pmTrend2D.c =================================================================== --- /trunk/psModules/src/objects/pmTrend2D.c (revision 15253) +++ /trunk/psModules/src/objects/pmTrend2D.c (revision 15254) @@ -3,6 +3,6 @@ * @author EAM, IfA * - * @version $Revision: 1.3 $ $Name: not supported by cvs2svn $ - * @date $Date: 2007-09-24 21:27:58 $ + * @version $Revision: 1.4 $ $Name: not supported by cvs2svn $ + * @date $Date: 2007-10-09 19:26:25 $ * * Copyright 2004 Institute for Astronomy, University of Hawaii @@ -48,7 +48,7 @@ for (int ny = 0; ny < trend->poly->nY + 1; ny++) { if (nx + ny >= PS_MAX (trend->poly->nX, trend->poly->nY) + 1) { - trend->poly->mask[nx][ny] = 1; + trend->poly->coeffMask[nx][ny] = PS_POLY_MASK_SET; } else { - trend->poly->mask[nx][ny] = 0; + trend->poly->coeffMask[nx][ny] = PS_POLY_MASK_NONE; } } @@ -96,7 +96,7 @@ for (int ny = 0; ny < trend->poly->nY + 1; ny++) { if (nx + ny >= PS_MAX (trend->poly->nX, trend->poly->nY) + 1) { - trend->poly->mask[nx][ny] = 1; + trend->poly->coeffMask[nx][ny] = PS_POLY_MASK_SET; } else { - trend->poly->mask[nx][ny] = 0; + trend->poly->coeffMask[nx][ny] = PS_POLY_MASK_NONE; } } @@ -137,7 +137,7 @@ for (int ny = 0; ny < trend->poly->nY + 1; ny++) { if (nx + ny >= PS_MAX (trend->poly->nX, trend->poly->nY) + 1) { - trend->poly->mask[nx][ny] = 1; + trend->poly->coeffMask[nx][ny] = PS_POLY_MASK_SET; } else { - trend->poly->mask[nx][ny] = 0; + trend->poly->coeffMask[nx][ny] = PS_POLY_MASK_NONE; } }