Changeset 42506 for branches/eam_branches/ipp-20230313/psLib/src/math/psPolynomial.c

Timestamp:

Aug 16, 2023, 10:30:15 AM (3 years ago)

Author:

eugene

Message:

add 2D, 3D, 4D versions of psVectorIRLSFitPolynomial, psPolynomialCopy, psPolynomialRecycle

File:

• branches/eam_branches/ipp-20230313/psLib/src/math/psPolynomial.c (modified) (7 diffs)

branches/eam_branches/ipp-20230313/psLib/src/math/psPolynomial.c

@@ -944 +944 @@
     bool match = true;
     match &= (poly->type == type);
-    match &= (poly->nX == type);
+    match &= (poly->nX   == nX);
 
     if (!match) {
-        psFree (poly->coeff);
-        psFree (poly->coeffErr);
-        psFree (poly->coeffMask);
+        polynomial1DFree (poly); // frees the coeffs
 
         poly->type = type;
-        poly->nX = nX;
+        poly->nX   = nX;
 
         poly->coeff = psAlloc((1 + nX) * sizeof(psF64));
@@ -976 +974 @@
     bool match = true;
     match &= (poly->type == type);
-    match &= (poly->nX == type);
-    match &= (poly->nY == type);
+    match &= (poly->nX ==   nX);
+    match &= (poly->nY ==   nY);
 
     if (!match) {
-        for (int i = 0; i < poly->nX + 1; i++) {
-            psFree (poly->coeff[i]);
-            psFree (poly->coeffErr[i]);
-            psFree (poly->coeffMask[i]);
-        }
-        psFree (poly->coeff);
-        psFree (poly->coeffErr);
-        psFree (poly->coeffMask);
+        polynomial2DFree (poly); // frees the coeffs
 
         poly->type = type;
-        poly->nX = nX;
-        poly->nY = nY;
+        poly->nX   = nX;
+        poly->nY   = nY;
 
         poly->coeff = psAlloc((1 + nX) * sizeof(psF64 *));
@@ -1012 +1003 @@
 }
 
-// XXX 3D, 4D versions
+bool psPolynomial3DRecycle(psPolynomial3D *poly,
+                           psPolynomialType type,
+                           unsigned int nX,
+                           unsigned int nY,
+                           unsigned int nZ)
+{
+    PS_ASSERT_INT_NONNEGATIVE(nX, NULL);
+    PS_ASSERT_INT_NONNEGATIVE(nY, NULL);
+    PS_ASSERT_INT_NONNEGATIVE(nZ, NULL);
+
+    bool match = true;
+    match &= (poly->type == type);
+    match &= (poly->nX   == nX);
+    match &= (poly->nY   == nY);
+    match &= (poly->nZ   == nZ);
+
+    if (!match) {
+        polynomial3DFree (poly); // frees the coeffs
+
+        poly->type = type;
+        poly->nX = nX;
+        poly->nY = nY;
+        poly->nZ = nZ;
+
+        poly->coeff = psAlloc((nX + 1) * sizeof(psF64 **));
+        poly->coeffErr = psAlloc((nX + 1) * sizeof(psF64 **));
+        poly->coeffMask = (psMaskType ***)psAlloc((nX + 1) * sizeof(psMaskType **));
+        for (int ix = 0; ix < (1 + nX); ix++) {
+            poly->coeff[ix] = psAlloc((nY + 1) * sizeof(psF64 *));
+            poly->coeffErr[ix] = psAlloc((nY + 1) * sizeof(psF64 *));
+            poly->coeffMask[ix] = (psMaskType **)psAlloc((nY + 1) * sizeof(psMaskType *));
+            for (int iy = 0; iy < (nY + 1); iy++) {
+                poly->coeff[ix][iy] = psAlloc((nZ + 1) * sizeof(psF64));
+                poly->coeffErr[ix][iy] = psAlloc((nZ + 1) * sizeof(psF64));
+                poly->coeffMask[ix][iy] = (psMaskType *)psAlloc((nZ + 1) * sizeof(psMaskType));
+            }
+        }
+    }
+    for (int ix = 0; ix < (1 + nX); ix++) {
+        for (int iy = 0; iy < (1 + nY); iy++) {
+            for (int iz = 0; iz < (1 + nZ); iz++) {
+                poly->coeff[ix][iy][iz]     = 0.0;
+                poly->coeffErr[ix][iy][iz]  = 0.0;
+                poly->coeffMask[ix][iy][iz] = PS_POLY_MASK_NONE;
+            }
+        }
+    }
+    return(true);
+}
+
+bool psPolynomial4DRecycle(psPolynomial4D *poly,
+                           psPolynomialType type,
+                           unsigned int nX,
+                           unsigned int nY,
+                           unsigned int nZ,
+                           unsigned int nT)
+{
+    PS_ASSERT_INT_NONNEGATIVE(nX, NULL);
+    PS_ASSERT_INT_NONNEGATIVE(nY, NULL);
+    PS_ASSERT_INT_NONNEGATIVE(nZ, NULL);
+    PS_ASSERT_INT_NONNEGATIVE(nT, NULL);
+
+    bool match = true;
+    match &= (poly->type == type);
+    match &= (poly->nX   == nX);
+    match &= (poly->nY   == nY);
+    match &= (poly->nZ   == nZ);
+    match &= (poly->nT   == nT);
+
+    if (!match) {
+        polynomial4DFree (poly); // frees the coeffs
+
+        poly->type = type;
+        poly->nX = nX;
+        poly->nY = nY;
+        poly->nZ = nZ;
+        poly->nT = nT;
+
+        poly->coeff = psAlloc((nX + 1) * sizeof(psF64 ***));
+        poly->coeffErr = psAlloc((nX + 1) * sizeof(psF64 ***));
+        poly->coeffMask = (psMaskType ****)psAlloc((nX + 1) * sizeof(psMaskType ***));
+        for (int ix = 0; ix < (nX + 1); ix++) {
+            poly->coeff[ix] = psAlloc((nY + 1) * sizeof(psF64 **));
+            poly->coeffErr[ix] = psAlloc((nY + 1) * sizeof(psF64 **));
+            poly->coeffMask[ix] = (psMaskType ***)psAlloc((nY + 1) * sizeof(psMaskType **));
+            for (int iy = 0; iy < (nY + 1); iy++) {
+                poly->coeff[ix][iy] = psAlloc((nZ + 1) * sizeof(psF64 *));
+                poly->coeffErr[ix][iy] = psAlloc((nZ + 1) * sizeof(psF64 *));
+                poly->coeffMask[ix][iy] = (psMaskType **)psAlloc((nZ + 1) * sizeof(psMaskType *));
+                for (int iz = 0; iz < (nZ + 1); iz++) {
+                    poly->coeff[ix][iy][iz] = psAlloc((nT + 1) * sizeof(psF64));
+                    poly->coeffErr[ix][iy][iz] = psAlloc((nT + 1) * sizeof(psF64));
+                    poly->coeffMask[ix][iy][iz] = (psMaskType *)psAlloc((nT + 1) * sizeof(psMaskType));
+                }
+            }
+        }
+    }
+    for (int ix = 0; ix < (1 + nX); ix++) {
+        for (int iy = 0; iy < (1 + nY); iy++) {
+            for (int iz = 0; iz < (1 + nZ); iz++) {
+                for (int it = 0; it < (1 + nT); it++) {
+                    poly->coeff[ix][iy][iz][it]     = 0.0;
+                    poly->coeffErr[ix][iy][iz][it]  = 0.0;
+                    poly->coeffMask[ix][iy][iz][it] = PS_POLY_MASK_NONE;
+                }
+            }
+        }
+    }
+    return(true);
+}
+
+// ######## Copy polynomials ########
 psPolynomial1D *psPolynomial1DCopy(psPolynomial1D *out,
-                                   psPolynomial1D *poly)
+                                   psPolynomial1D *poly)
 {
     if (out == NULL) {
-        out = psPolynomial1DAlloc (poly->type, poly->nX);
+        out = psPolynomial1DAlloc (poly->type, poly->nX);
     } else {
-        psPolynomial1DRecycle (out, poly->type, poly->nX);
+        psPolynomial1DRecycle (out, poly->type, poly->nX);
     }
 
@@ -1030 +1132 @@
 }
 psPolynomial2D *psPolynomial2DCopy(psPolynomial2D *out,
-                                   psPolynomial2D *poly)
+                                   psPolynomial2D *poly)
 {
     if (out == NULL) {
-        out = psPolynomial2DAlloc (poly->type, poly->nX, poly->nY);
+        out = psPolynomial2DAlloc (poly->type, poly->nX, poly->nY);
     } else {
-        psPolynomial2DRecycle (out, poly->type, poly->nX, poly->nY);
+        psPolynomial2DRecycle (out, poly->type, poly->nX, poly->nY);
     }
 
     for (int i = 0; i < (1 + poly->nX); i++) {
-        for (int j = 0; j < (1 + poly->nY); j++) {
-            out->coeff[i][j] = poly->coeff[i][j];
-            out->coeffErr[i][j] = poly->coeffErr[i][j];
-            out->coeffMask[i][j] = poly->coeffMask[i][j];
-        }
+        for (int j = 0; j < (1 + poly->nY); j++) {
+            out->coeff[i][j] = poly->coeff[i][j];
+            out->coeffErr[i][j] = poly->coeffErr[i][j];
+            out->coeffMask[i][j] = poly->coeffMask[i][j];
+        }
+    }
+    return(out);
+}
+psPolynomial3D *psPolynomial3DCopy(psPolynomial3D *out,
+                                   psPolynomial3D *poly)
+{
+    if (out == NULL) {
+        out = psPolynomial3DAlloc (poly->type, poly->nX, poly->nY, poly->nZ);
+    } else {
+        psPolynomial3DRecycle (out, poly->type, poly->nX, poly->nY, poly->nZ);
+    }
+
+    for (int ix = 0; ix < (1 + poly->nX); ix++) {
+        for (int iy = 0; iy < (1 + poly->nY); iy++) {
+            for (int iz = 0; iz < (1 + poly->nZ); iz++) {
+                out->coeff[ix][iy][iz] = poly->coeff[ix][iy][iz];
+                out->coeffErr[ix][iy][iz] = poly->coeffErr[ix][iy][iz];
+                out->coeffMask[ix][iy][iz] = poly->coeffMask[ix][iy][iz];
+            }
+        }
+    }
+    return(out);
+}
+psPolynomial4D *psPolynomial4DCopy(psPolynomial4D *out,
+                                   psPolynomial4D *poly)
+{
+    if (out == NULL) {
+        out = psPolynomial4DAlloc (poly->type, poly->nX, poly->nY, poly->nZ, poly->nT);
+    } else {
+        psPolynomial4DRecycle (out, poly->type, poly->nX, poly->nY, poly->nZ, poly->nT);
+    }
+
+    for (int ix = 0; ix < (1 + poly->nX); ix++) {
+        for (int iy = 0; iy < (1 + poly->nY); iy++) {
+            for (int iz = 0; iz < (1 + poly->nZ); iz++) {
+                for (int it = 0; it < (1 + poly->nT); it++) {
+                    out->coeff[ix][iy][iz][it] = poly->coeff[ix][iy][iz][it];
+                    out->coeffErr[ix][iy][iz][it] = poly->coeffErr[ix][iy][iz][it];
+                    out->coeffMask[ix][iy][iz][it] = poly->coeffMask[ix][iy][iz][it];
+                }
+            }
+        }
     }
     return(out);
@@ -1080 +1224 @@
 
     switch (x->type.type) {
-    case PS_TYPE_F64:
+      case PS_TYPE_F64:
         tmp = psVectorAlloc(x->n, PS_TYPE_F64);
         for (unsigned int i=0;i<x->n;i++) {
@@ -1086 +1230 @@
         }
         break;
-    case PS_TYPE_F32:
+      case PS_TYPE_F32:
         tmp = psVectorAlloc(x->n, PS_TYPE_F32);
         for (unsigned int i=0;i<x->n;i++) {
@@ -1092 +1236 @@
         }
         break;
-    default:
+      default:
         psError(PS_ERR_UNKNOWN, false, "invalid input data type.\n");
         return (NULL);