Changeset 42967 for trunk/Ohana/src/opihi/cmd.data/vweave.c

Timestamp:

Mar 9, 2026, 11:26:15 AM (4 months ago)

Author:

eugene

Message:

add spline options to extrapolate or saturate at end points; add vweave function to interpolate a vector with varying bin sizes

File:

• trunk/Ohana/src/opihi/cmd.data/vweave.c (copied) (copied from branches/eam_branches/ipp-ppmerge-20241229/Ohana/src/opihi/cmd.data/vweave.c ) (3 diffs)

trunk/Ohana/src/opihi/cmd.data/vweave.c

@@ -6 +6 @@
   
   int i, j, N;
-  Vector *xin, *yin, *xout, *yout;
+  Vector *xin, *yin, *yout;
+  Vector *xminV = NULL, *xmaxV = NULL, *xout = NULL;
   opihi_flt xmin, xmax, Xmin, Xmax;
 
@@ -22 +23 @@
   }
 
-  if (argc != 5) {
+  if ((argc != 5) && (argc != 6)) {
     gprint (GP_ERR, "USAGE: vweave <xin> <yin> <xout> <yout>\n");
+    gprint (GP_ERR, "   OR: vweave <xin> <yin> <xmin> <xmax> <yout>\n");
     gprint (GP_ERR, " fully interpolate input bin values into output bins, using fractional input bin contributions\n");
     gprint (GP_ERR, " by default, yout has the sum of the associated input values\n");
     gprint (GP_ERR, " use -count to find the fractional input bin contributions\n");
-    gprint (GP_ERR, " output bin boundaries defined as midpoint of xout values\n");
+    gprint (GP_ERR, " with a single <xout> vector, output bin boundaries defined as midpoint of xout values\n");
+    gprint (GP_ERR, " with <xmin> and <xmax>, the output bin boundaries are explicitly provided\n");
     gprint (GP_ERR, " use -binsize to specify a fixed bin width (<xout> will define the bin center)\n");
     gprint (GP_ERR, " input values are treated as uniformly spread between midpoints of xin values\n");
@@ -38 +41 @@
     return FALSE;
   }
-  if ((xout = SelectVector (argv[3], OLDVECTOR, TRUE)) == NULL) return (FALSE);
-  if (xout->Nelements < 2) {
-    gprint (GP_ERR, "output vector of length < 2 is invalid\n");
-    return FALSE;
+  if ((yin  = SelectVector (argv[2], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+
+  int Nout = -1;
+  if (argc == 5) {
+    if ((xout = SelectVector (argv[3], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+    if (xout->Nelements < 2) {
+      gprint (GP_ERR, "output vector of length < 2 is invalid\n");
+      return FALSE;
+    }
+    Nout = xout->Nelements;
+    if ((yout = SelectVector (argv[4], ANYVECTOR, TRUE)) == NULL) return (FALSE);
+  } else {
+    if (!isnan(binsize)) {
+      gprint (GP_ERR, "-binsize is incompatible with <xmin>, <xmax>\n");
+      return FALSE;
+    }
+    if ((xminV = SelectVector (argv[3], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+    if (xminV->Nelements < 2) {
+      gprint (GP_ERR, "output vector of length < 2 is invalid\n");
+      return FALSE;
+    }
+    Nout = xminV->Nelements;
+    if ((xmaxV = SelectVector (argv[4], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+    if (xmaxV->Nelements != Nout) {
+      gprint (GP_ERR, "<xmin> and <xmax> lengths must match\n");
+      return FALSE;
+    }
+    if ((yout = SelectVector (argv[5], ANYVECTOR, TRUE)) == NULL) return (FALSE);
   }
 
-  if ((yin  = SelectVector (argv[2], OLDVECTOR, TRUE)) == NULL) return (FALSE);
-  if ((yout = SelectVector (argv[4], ANYVECTOR, TRUE)) == NULL) return (FALSE);
-
   // re-binning creates a float vector
-  ResetVector (yout, OPIHI_FLT, xout[0].Nelements);
+  ResetVector (yout, OPIHI_FLT, Nout);
 
   // XXX this algorithm is not optimized: full scale of input vector for each output bin
   // could speed up with sorting and pre-defined boundaries
-  for (i = 0; i < xout[0].Nelements; i++) {
+  for (i = 0; i < Nout; i++) {
     if (!isnan(binsize)) {
       xmin = xout[0].elements.Flt[i] - 0.5*binsize;
       xmax = xout[0].elements.Flt[i] + 0.5*binsize;
     } else {
-      if (i > 0) {
-        xmin = 0.5*(xout[0].elements.Flt[i] + xout[0].elements.Flt[i-1]);
+      if (xout) {
+        if (i > 0) {
+          xmin = 0.5*(xout[0].elements.Flt[i] + xout[0].elements.Flt[i-1]);
+        } else {
+          xmin = 1.5*xout[0].elements.Flt[i] - 0.5*xout[0].elements.Flt[i+1];
+        }
+        if (i < xout[0].Nelements - 1) {
+          xmax = 0.5*(xout[0].elements.Flt[i] + xout[0].elements.Flt[i+1]);
+        } else {
+          xmax = 1.5*xout[0].elements.Flt[i] - 0.5*xout[0].elements.Flt[i-1];
+        }
       } else {
-        xmin = 1.5*xout[0].elements.Flt[i] + 0.5*xout[0].elements.Flt[i+1];
-      }
-      if (i < xout[0].Nelements - 1) {
-        xmax = 0.5*(xout[0].elements.Flt[i] + xout[0].elements.Flt[i+1]);
-      } else {
-        xmax = 1.5*xout[0].elements.Flt[i] - 0.5*xout[0].elements.Flt[i-1];
+        xmin = xminV[0].elements.Flt[i];
+        xmax = xmaxV[0].elements.Flt[i];
       }
     }