Index: trunk/Ohana/src/opihi/cmd.data/Makefile
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/Makefile	(revision 40165)
+++ trunk/Ohana/src/opihi/cmd.data/Makefile	(revision 40291)
@@ -155,4 +155,5 @@
 $(SRC)/type.$(ARCH).o		   \
 $(SRC)/uniq.$(ARCH).o		   \
+$(SRC)/uniqpair.$(ARCH).o		   \
 $(SRC)/unsign.$(ARCH).o	           \
 $(SRC)/vbin.$(ARCH).o		   \
Index: trunk/Ohana/src/opihi/cmd.data/init.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/init.c	(revision 40165)
+++ trunk/Ohana/src/opihi/cmd.data/init.c	(revision 40291)
@@ -141,6 +141,7 @@
 int tvcontour        PROTO((int, char **));
 int tvgrid           PROTO((int, char **));
-int opihi_type             PROTO((int, char **));
+int opihi_type       PROTO((int, char **));
 int uniq             PROTO((int, char **));
+int uniqpair         PROTO((int, char **));
 int unsign           PROTO((int, char **));
 int vbin             PROTO((int, char **));
@@ -325,4 +326,5 @@
   {1, "ungridify",    ungridify,        "convert image region to vector triplet"},
   {1, "uniq",         uniq,             "create a uniq vector subset from a vector"},
+  {1, "uniqpair",     uniqpair,         "create a uniq vector subset from a pair of vectors, saving duplicates if desired"},
   {1, "unsign",       unsign,           "toggle the UNSIGN status"},
   {1, "vbin",         vbin,             "rebin vector data by a factor of N"},
Index: trunk/Ohana/src/opihi/cmd.data/limits.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/limits.c	(revision 40165)
+++ trunk/Ohana/src/opihi/cmd.data/limits.c	(revision 40291)
@@ -3,7 +3,6 @@
 int limits (int argc, char **argv) {
 
-  int N, APPLY, dX, dY;
+  int N, dX, dY;
   int kapa;
-  char *name;
   Graphdata graphmode;
   Vector *xvec, *yvec;
@@ -11,10 +10,56 @@
   xvec = yvec = NULL;
 
-  APPLY = FALSE;
+  float minLimitX = NAN;
+  float minLimitY = NAN;
+  float maxLimitX = NAN;
+  float maxLimitY = NAN;
+  float delLimitX = NAN;
+  float delLimitY = NAN;
+
+  if ((N = get_argument (argc, argv, "-minX"))) {
+    remove_argument (N, &argc, argv);
+    minLimitX = atof (argv[N]);
+    remove_argument (N, &argc, argv);
+  }
+  if ((N = get_argument (argc, argv, "-maxX"))) {
+    remove_argument (N, &argc, argv);
+    maxLimitX = atof (argv[N]);
+    remove_argument (N, &argc, argv);
+  }
+  if ((N = get_argument (argc, argv, "-delX"))) {
+    if (!isnan(minLimitX) || !isnan(maxLimitX)) {
+      gprint (GP_ERR, "-minX & -maxX cannot be mixed with -delX\n");
+      return (FALSE);
+    }
+    remove_argument (N, &argc, argv);
+    delLimitX = atof (argv[N]);
+    remove_argument (N, &argc, argv);
+  }
+  if ((N = get_argument (argc, argv, "-minY"))) {
+    remove_argument (N, &argc, argv);
+    minLimitY = atof (argv[N]);
+    remove_argument (N, &argc, argv);
+  }
+  if ((N = get_argument (argc, argv, "-maxY"))) {
+    remove_argument (N, &argc, argv);
+    maxLimitY = atof (argv[N]);
+    remove_argument (N, &argc, argv);
+  }
+  if ((N = get_argument (argc, argv, "-delY"))) {
+    if (!isnan(minLimitY) || !isnan(maxLimitY)) {
+      gprint (GP_ERR, "-minY & -maxY cannot be mixed with -delY\n");
+      return (FALSE);
+    }
+    remove_argument (N, &argc, argv);
+    delLimitY = atof (argv[N]);
+    remove_argument (N, &argc, argv);
+  }
+
+  int APPLY = FALSE;
   if ((N = get_argument (argc, argv, "-a"))) {
     remove_argument (N, &argc, argv);
     APPLY = TRUE;
   }
-  name = NULL;
+  char *name = NULL;
   if ((N = get_argument (argc, argv, "-n"))) {
     remove_argument (N, &argc, argv);
@@ -22,4 +67,5 @@
     remove_argument (N, &argc, argv);
   }
+
   if (!GetGraph (&graphmode, &kapa, name)) return (FALSE);
   FREE (name);
@@ -97,4 +143,27 @@
  success:
   SetLimits (xvec, yvec, &graphmode);
+
+  if (!isnan(minLimitX)) graphmode.xmin = MIN (minLimitX, graphmode.xmin);
+  if (!isnan(maxLimitX)) graphmode.xmax = MAX (maxLimitX, graphmode.xmax);
+  if (!isnan(minLimitY)) graphmode.ymin = MIN (minLimitY, graphmode.ymin);
+  if (!isnan(maxLimitY)) graphmode.ymax = MAX (maxLimitY, graphmode.ymax);
+
+  if (!isnan(delLimitX)) {
+    float delta = graphmode.xmax - graphmode.xmin;
+    if (fabs(delLimitX) > fabs(delta)) {
+      float midpt = 0.5*(graphmode.xmax + graphmode.xmin);
+      graphmode.xmax = midpt + 0.5*delLimitX;
+      graphmode.xmin = midpt - 0.5*delLimitX;
+    }
+  }
+  if (!isnan(delLimitY)) {
+    float delta = graphmode.ymax - graphmode.ymin;
+    if (fabs(delLimitY) > fabs(delta)) {
+      float midpt = 0.5*(graphmode.ymax + graphmode.ymin);
+      graphmode.ymax = midpt + 0.5*delLimitY;
+      graphmode.ymin = midpt - 0.5*delLimitY;
+    }
+  }
+
   if (APPLY) KapaSetLimits (kapa, &graphmode);
   return (TRUE);
Index: trunk/Ohana/src/opihi/cmd.data/print_vectors.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/print_vectors.c	(revision 40165)
+++ trunk/Ohana/src/opihi/cmd.data/print_vectors.c	(revision 40291)
@@ -4,5 +4,19 @@
 
   Vector **vec;
-  int i, j;
+  int i, j, N;
+
+  int START_VALUE = 0;
+  if ((N = get_argument (argc, argv, "-s"))) {
+    remove_argument (N, &argc, argv);
+    START_VALUE = atoi (argv[N]);
+    remove_argument (N, &argc, argv);
+  }
+
+  int END_VALUE = -1;
+  if ((N = get_argument (argc, argv, "-e"))) {
+    remove_argument (N, &argc, argv);
+    END_VALUE = atoi (argv[N]);
+    remove_argument (N, &argc, argv);
+  }
 
   if (argc < 2) {
@@ -27,5 +41,13 @@
   }
 
-  for (j = 0; j < MaxLen; j++) {
+  // start and end may be 0 - N (truncated to N) or may be negative, in which case it refers to 
+  // distance from the end (just like vector[-5])
+  START_VALUE = (START_VALUE < 0) ? MaxLen + START_VALUE + 1 : MIN (START_VALUE, MaxLen);
+  START_VALUE = MAX (0, START_VALUE);
+
+  END_VALUE = (END_VALUE < 0) ? MaxLen + END_VALUE + 1 : MIN (END_VALUE, MaxLen);
+  END_VALUE = MAX (0, END_VALUE);
+
+  for (j = START_VALUE; j < END_VALUE; j++) {
     for (i = 0; i < Nvec; i++) {
       if (j >= vec[i][0].Nelements) {
@@ -35,5 +57,5 @@
 	  gprint (GP_LOG, "%f ", vec[i][0].elements.Flt[j]);
 	} else {
-	  gprint (GP_LOG, "%d ", vec[i][0].elements.Int[j]);
+	  gprint (GP_LOG, OPIHI_INT_FMT" ", vec[i][0].elements.Int[j]);
 	}
       }
Index: trunk/Ohana/src/opihi/cmd.data/reindex.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/reindex.c	(revision 40165)
+++ trunk/Ohana/src/opihi/cmd.data/reindex.c	(revision 40291)
@@ -47,5 +47,5 @@
 	continue;
       }
-      if (*vx > Nmax) ESCAPE("unexpected value in index: %d (%d)\n", *vx, i);
+      if (*vx > Nmax) ESCAPE("unexpected value in index: "OPIHI_INT_FMT" (%d)\n", *vx, i);
       ovec[0].elements.Flt[Npts] = vi[*vx];
       Npts++;
@@ -67,5 +67,5 @@
 	continue;
       }
-      if (*vx > Nmax) ESCAPE("unexpected value in index: %d (%d)\n", *vx, i);
+      if (*vx > Nmax) ESCAPE("unexpected value in index: "OPIHI_INT_FMT" (%d)\n", *vx, i);
       ovec[0].elements.Int[Npts] = vi[*vx];
       Npts++;
Index: trunk/Ohana/src/opihi/cmd.data/test/periodogram-fm.sh
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/test/periodogram-fm.sh	(revision 40165)
+++ trunk/Ohana/src/opihi/cmd.data/test/periodogram-fm.sh	(revision 40291)
@@ -55,5 +55,4 @@
 
  periodogram_fm t f df 5 50 period power
-#periodogram t f 5 50 period power
 
  # lim -n 0 t f; clear; box; plot -x 2 -pt 2 t f
@@ -88,5 +87,4 @@
 
  periodogram_fm t f df 1 10 period power
-#periodogram t f 1 10 period power
 
  # lim -n 0 t f; clear; box; plot -x 2 -pt 2 t f
@@ -97,4 +95,8 @@
  if (abs ($peakpos - $P) > 0.05)
    $PASS = 0
+ end
+
+  if ($PLOT)
+  lim period power; clear; box; line -c red70 -lw 3 $P 0 to $P $peakval; plot period power -x line
  end
 end
@@ -117,5 +119,4 @@
 
  periodogram_fm t f df 2 30 period power
-#periodogram t f 2 30 period power
 
 #  lim -n 0 t f; clear; box; plot -x 2 -pt 2 t f
@@ -126,4 +127,8 @@
  if (abs ($peakpos - $P) > 0.05)
    $PASS = 0
+ end
+
+ if ($PLOT)
+  lim period power; clear; box; line -c red70 -lw 3 $P 0 to $P $peakval; plot period power -x line
  end
 end
@@ -145,5 +150,5 @@
  set df = 0.01 + zero(f)
 
- periodogram_fm t f 2 30 period power
+ periodogram_fm t f df 2 30 period power
 
 #  lim -n 0 t f; clear; box; plot -x 2 -pt 2 t f
@@ -155,7 +160,11 @@
    $PASS = 0
  end
-end
-
-# test using random samples, offset start, non-zero DC, some noise
+
+ if ($PLOT)
+  lim period power; clear; box; line -c red70 -lw 3 $P 0 to $P $peakval; plot period power -x line
+ end
+end
+
+# test using 300 random samples, offset start, non-zero DC, some noise
 macro test6
  $PASS = 1
@@ -178,5 +187,5 @@
  set df = 0.01 + zero(f)
 
- periodogram_fm t f 2 30 period power
+ periodogram_fm t f df 2 30 period power
 
 #  lim -n 0 t f; clear; box; plot -x 2 -pt 2 t f
@@ -188,7 +197,11 @@
    $PASS = 0
  end
-end
-
-# test using fewer random samples, offset start, non-zero DC, some noise
+
+ if ($PLOT)
+  lim period power; clear; box; line -c red70 -lw 3 $P 0 to $P $peakval; plot period power -x line
+ end
+end
+
+# test using 100 fewer random samples, offset start, non-zero DC, some noise
 macro test7
  $PASS = 1
@@ -211,5 +224,5 @@
  set df = 0.01 + zero(f)
 
- periodogram_fm t f 2 30 period power
+ periodogram_fm t f df 2 30 period power
 
 #  lim -n 0 t f; clear; box; plot -x 2 -pt 2 t f
@@ -221,22 +234,28 @@
    $PASS = 0
  end
-end
-
-# test using fewer random samples, high frequency, non-zero DC, some noise
+
+ if ($PLOT)
+  lim period power; clear; box; line -c red70 -lw 3 $P 0 to $P $peakval; plot period power -x line
+ end
+end
+
+# test using Ndays random samples, RR Lyrae-sized light curves (0.7 mag),
+# optional noise level 
 macro test8
- if ($0 != 2)
-   echo "USAGE: test8: Ndays")
+ if ($0 != 4)
+   echo "USAGE: test8: Period Ndays (df)"
    break
  end
  
  local Ndays 
- $Ndays = $1
-
- $PASS = 1
- break -auto off
-
- local P PI
- $PI = 3.14159265359
- $P  = 0.8*rnd(0) + 0.2
+ $P = $1
+ $Ndays = $2
+ $dM = $3
+
+ $PASS = 1
+ break -auto off
+
+ local PI
+ $PI = 3.14159265359
  $trueP = $P
 
@@ -246,7 +265,7 @@
 
  # t is a time in days, but we always have 4 within 1 hour:
- set t0 = int(100 * rnd(x))
- set dtx = (3/24) * rnd(x)
- set t0 = t0 + dtx
+ set tday = int(100 * rnd(x)); # choose Ndays random days between 0 and 100
+ set dtx = (3/24) * rnd(x);  # choose a starting time within that night
+ set t0 = tday + dtx
 
  set dt1 = (15.0 / 1440) * rnd(x) + ( 0 + 7.5) / 1440
@@ -260,13 +279,13 @@
  set tmp = t0 + dt3; concat tmp t
 
- set fraw = sin(2*$PI*t/$P) + 0.5
+ set fraw = 0.75*sin(2*$PI*t/$P)
 
  # 0.05 : peakpos = 14.95
  # 0.10 : peakpos = 15.04 (
- gaussdev df t[] 0.0 0.25
+ gaussdev df t[] 0.0 $dM
  set f = fraw + df
- set df = 0.01 + zero(f)
-
- periodogram_fm t f 0.1 2.0 period power
+ set df = $dM + zero(f)
+
+ periodogram_fm t f df 0.1 20.0 period power
 
 #  lim -n 0 t f; clear; box; plot -x 2 -pt 2 t f
@@ -278,5 +297,233 @@
    $PASS = 0
  end
-end
+ if ($PLOT)
+  lim period power; clear; box; line -c red70 -lw 3 $P 0 to $P $peakval; plot period power -x line
+ end
+end
+
+# test using Ndays random samples, RR Lyrae-sized light curves (0.7 mag),
+# optional noise level 
+# compare periodogram and periodogram_fm
+macro test9
+ if ($0 != 4)
+   echo "USAGE: test8: Period Ndays (df)"
+   break
+ end
+ 
+ local Ndays 
+ $P = $1
+ $Ndays = $2
+ $dM = $3
+
+ $PASS = 1
+ break -auto off
+
+ local PI
+ $PI = 3.14159265359
+ $trueP = $P
+
+ delete -q x t f period power
+
+ create x 0 $Ndays
+
+ # t is a time in days, but we always have 4 within 1 hour:
+ set tday = int(100 * rnd(x)); # choose Ndays random days between 0 and 100
+ set dtx = (3/24) * rnd(x);  # choose a starting time within that night
+ set t0 = tday + dtx
+
+ set dt1 = (15.0 / 1440) * rnd(x) + ( 0 + 7.5) / 1440
+ set dt2 = (15.0 / 1440) * rnd(x) + (15 + 7.5) / 1440
+ set dt3 = (15.0 / 1440) * rnd(x) + (30 + 7.5) / 1440
+
+ delete -q t
+ concat t0 t
+ set tmp = t0 + dt1; concat tmp t
+ set tmp = t0 + dt2; concat tmp t
+ set tmp = t0 + dt3; concat tmp t
+
+ set fraw = 0.75*sin(2*$PI*t/$P)
+
+ # 0.05 : peakpos = 14.95
+ # 0.10 : peakpos = 15.04 (
+ gaussdev df t[] 0.0 $dM
+ set f = fraw + df
+ set df = $dM + zero(f)
+
+ periodogram_fm t f df 0.1 20.0 period_fm power_fm
+ periodogram t f 0.1 20.0 period power
+
+#  lim -n 0 t f; clear; box; plot -x 2 -pt 2 t f
+#  lim -n 1 period power; clear; box; plot period power
+
+ peak -q period_fm power_fm
+ $peakval_fm = $peakval
+
+ peak -q period power
+ vstat -q power
+ set power = power / $MAX
+
+# if (abs ($peakpos - $P) > 0.05)
+#   $PASS = 0
+# end
+
+ set freq = 1 / period
+ set freq_fm = 1 / period_fm
+ $Freq = 1 / $P
+
+ if ($PLOT)
+  if (1)
+    lim period power; clear; box
+    line -c red70 -lw 3 $P 0 to $P $peakval_fm;
+    plot period power -x line -c grey70 -lw 2
+    plot period_fm power_fm -x line -c black
+  else
+    lim freq power; clear; box
+    line -c red70 -lw 3 $Freq 0 to $Freq 1.0
+    plot freq power -x line -c grey70 -lw 2
+    plot freq_fm power_fm -x line -c black
+  end
+ end
+end
+
+
+# test using Ndays random samples, RR Lyrae-sized light curves (0.7 mag),
+# optional noise level 
+# compare periodogram and periodogram_fm
+macro test10
+ if ($0 != 4)
+   echo "USAGE: test8: Period Ndays (df)"
+   break
+ end
+ 
+ local Ndays 
+ $P = $1
+ $Ndays = $2
+ $dM = $3
+
+ $PASS = 1
+ break -auto off
+
+ local PI
+ $PI = 3.14159265359
+ $trueP = $P
+
+ delete -q x t f period power
+
+ create x 0 $Ndays
+
+ # t is a time in days, but we always have 4 within 1 hour:
+ set tday = int(100 * rnd(x)); # choose Ndays random days between 0 and 100
+ set dtx = (3/24) * rnd(x);  # choose a starting time within that night
+ set t0 = tday + dtx
+
+ set dt1 = (15.0 / 1440) * rnd(x) + ( 0 + 7.5) / 1440
+ set dt2 = (15.0 / 1440) * rnd(x) + (15 + 7.5) / 1440
+ set dt3 = (15.0 / 1440) * rnd(x) + (30 + 7.5) / 1440
+
+ delete -q t
+ concat t0 t
+ set tmp = t0 + dt1; concat tmp t
+ set tmp = t0 + dt2; concat tmp t
+ set tmp = t0 + dt3; concat tmp t
+
+ set fraw = 0.75*sin(2*$PI*t/$P)
+
+ # 0.05 : peakpos = 14.95
+ # 0.10 : peakpos = 15.04 (
+ gaussdev df t[] 0.0 $dM
+ set f = fraw + df
+ set df = $dM + zero(f)
+
+ periodogram_fm t f df 0.05 20.0 period_fm power_fm
+
+ gaussdev df t[] 0.0 $dM
+ set Fo = df
+ periodogram_fm t Fo df 0.05 20.0 period power
+
+#  lim -n 0 t f; clear; box; plot -x 2 -pt 2 t f
+#  lim -n 1 period power; clear; box; plot period power
+
+ peak -q period_fm power_fm
+ $peakval_fm = $peakval
+
+ peak -q period power
+
+# if (abs ($peakpos - $P) > 0.05)
+#   $PASS = 0
+# end
+
+ set freq = 1 / period
+ set freq_fm = 1 / period_fm
+ $Freq = 1 / $P
+
+ if ($PLOT)
+  if (1)
+    lim period_fm power_fm; clear; box
+    line -c red70 -lw 3 $P 0 to $P $peakval_fm;
+    plot period power -x line -c grey70 -lw 2
+    plot period_fm power_fm -x line -c black
+  else
+    lim freq power; clear; box
+    line -c red70 -lw 3 $Freq 0 to $Freq 1.0
+    plot freq power -x line -c grey70 -lw 2
+    plot freq_fm power_fm -x line -c black
+  end
+ end
+end
+
+# we have time (MJD) and mag
+# we generate the folded lightcure and measure sigma relative to the smoothed version (bins of 0.1 period)
+macro fold.one.period
+  if ($0 != 5)
+    echo "USAGE: fold.one.period (time) (mag) (magErr) (period)"
+    break
+  end
+
+  local myTime myMag myMagErr myPeriod
+  $myTime = $1
+  $myMag  = $2
+  $myMagErr  = $3
+  $myPeriod = $4
+
+  set phi = $myTime / $myPeriod - int($myTime / $myPeriod)
+
+  if ($PLOT_FOLD)
+    lim -n phi phi $myMag; clear; box; 
+  end
+
+  delete -q magResid
+
+  $dPhi = 0.05; # half of bin size
+  create nphi $dPhi {1 + $dPhi} {2*$dPhi}
+  set magR = zero(nphi)
+  set magS = zero(nphi)
+  for i 0 nphi[]
+    subset tmp_mag_sub = $myMag where (phi >= nphi[$i] - $dPhi) && (phi < nphi[$i] + $dPhi)
+    vstat -q tmp_mag_sub
+    magR[$i] = $MEDIAN
+    magS[$i] = $SIGMA
+
+    set magDelta = tmp_mag_sub - $MEDIAN
+    concat magDelta magResid 
+
+    if ($PLOT_FOLD)
+      subset tmp_phi_sub = phi where (phi >= nphi[$i] - $dPhi) && (phi < nphi[$i] + $dPhi)
+      if ($i % 2)
+        plot tmp_phi_sub tmp_mag_sub -pt 7 -sz 3 -c blue -lw 2
+      else
+        plot tmp_phi_sub tmp_mag_sub -pt 7 -sz 3 -c red -lw 2
+      end
+    end  
+  end
+
+  if ($PLOT_FOLD)  
+    plot -pt 10 -sz 1.5 phi $myMag -dy $myMagErr
+    plot -pt 2 -sz 2.0 -c red nphi magR -dy magS
+  end
+
+  vstat -q magResid
+end
+
+
 
 # Memory test
@@ -298,5 +545,5 @@
 
  for i 0 100
-  periodogram_fm t f 2 30 period power
+  periodogram_fm t f df 2 30 period power
  end
   
Index: trunk/Ohana/src/opihi/cmd.data/test/periodogram.sh
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/test/periodogram.sh	(revision 40165)
+++ trunk/Ohana/src/opihi/cmd.data/test/periodogram.sh	(revision 40291)
@@ -209,18 +209,19 @@
 # test using fewer random samples, high frequency, non-zero DC, some noise
 macro test8
- if ($0 != 2)
-   echo "USAGE: test8: Ndays")
+ if ($0 != 4)
+   echo "USAGE: test8: Period Ndays df"
    break
  end
  
  local Ndays 
- $Ndays = $1
-
- $PASS = 1
- break -auto off
-
- local P PI
- $PI = 3.14159265359
- $P  = 0.8*rnd(0) + 0.2
+ $P = $1
+ $Ndays = $2
+ $dM = $3
+
+ $PASS = 1
+ break -auto off
+
+ local PI
+ $PI = 3.14159265359
  $trueP = $P
 
@@ -230,7 +231,7 @@
 
  # t is a time in days, but we always have 4 within 1 hour:
- set t0 = int(100 * rnd(x))
- set dtx = (3/24) * rnd(x)
- set t0 = t0 + dtx
+ set tday = int(100 * rnd(x)); # choose Ndays random days between 0 and 100
+ set dtx = (3/24) * rnd(x);  # choose a starting time within that night
+ set t0 = tday + dtx
 
  set dt1 = (15.0 / 1440) * rnd(x) + ( 0 + 7.5) / 1440
@@ -244,9 +245,9 @@
  set tmp = t0 + dt3; concat tmp t
 
- set fraw = sin(2*$PI*t/$P) + 0.5
+ set fraw = 0.75*sin(2*$PI*t/$P)
 
  # 0.05 : peakpos = 14.95
  # 0.10 : peakpos = 15.04 (
- gaussdev df t[] 0.0 0.25
+ gaussdev df t[] 0.0 $dM
  set f = fraw + df
 
@@ -260,4 +261,7 @@
  if (abs ($peakpos - $P) > 0.05)
    $PASS = 0
+ end
+ if ($PLOT)
+  lim period power; clear; box; line -c red70 -lw 3 $P 0 to $P $peakval; plot period power -x line
  end
 end
Index: trunk/Ohana/src/opihi/cmd.data/test/uniqpair.sh
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/test/uniqpair.sh	(revision 40291)
+++ trunk/Ohana/src/opihi/cmd.data/test/uniqpair.sh	(revision 40291)
@@ -0,0 +1,83 @@
+
+
+list tests
+ test1
+ memtest1
+end
+
+# Test if uniqpair works
+macro test1
+
+ $PASS = 1
+
+ local i
+
+ delete ID1 ID2
+
+ vlist -int ID1 1   1   2   3   3   4   5   6
+ vlist -int ID2 2   3   4   5   5   6   7   8
+
+ uniqpair ID1 ID2 IDu -c IDn -d index
+
+ # if ((xvec[1024] != 17) || (yvec[1024] != 100))
+ #  $PASS = 0
+ #  echo "Value mismatch: xvec[1024] yvec[1024] (should be 17,100)"
+ # end
+ # 
+ # imhist -q buff xvec yvec -region 40 0 25 10 -range 0 10
+ # 
+ # if ((xvec[1024] != 10) || (yvec[1024] != 100))
+ #  $PASS = 0
+ #  echo "Value mismatch: xvec[1024] yvec[1024] (should be 10,100)"
+ # end
+
+end
+
+# Test if uniqpair works
+macro test2
+
+ $PASS = 1
+
+ delete ID1 ID2
+
+ vlist -int ID1 {2^16 + 2} {2^17}     {2^18 + 3} {2^16 + 2} {2^17 + 2} {2^18 + 3} 
+ vlist -int ID2 {2^16 + 0} {2^17 + 2} {2^18 + 5} {2^16 + 0} {2^17 + 2} {2^18 + 5} 
+
+ uniqpair ID1 ID2 IDu -c IDn -d index
+
+ reindex ID2s = ID2 using index
+ reindex ID1s = ID1 using index
+
+ vectors
+
+ echo "unique"
+ print_v IDu IDn
+
+ echo "dups"
+ print_v ID1s ID2s index
+end
+
+
+# Memory test
+macro memtest1
+
+ local i
+
+ list word -x "ps -p $PID -o rss"
+ $startmem = $word:1
+
+ for i 0 1000
+  imhist -q buff xvec yvec -region 40 0 25 10 -range 0 10
+ end
+
+ list word -x "ps -p $PID -o rss"
+ $endmem = $word:1
+
+ $PASS = 1
+
+ if ($endmem - $startmem > 10)
+   $PASS = 0
+   echo "growth: {$endmem-$startmem}"
+   echo "kB/loop: {($endmem-$startmem)/1000}"
+ end
+end
Index: trunk/Ohana/src/opihi/cmd.data/uniq.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/uniq.c	(revision 40165)
+++ trunk/Ohana/src/opihi/cmd.data/uniq.c	(revision 40291)
@@ -80,5 +80,5 @@
     memcpy (indata, ivec->elements.Int, ivec[0].Nelements*sizeof(opihi_int));
 
-    isort (indata, ivec->Nelements);
+    llsort (indata, ivec->Nelements);
 
     Nnew = 0;
Index: trunk/Ohana/src/opihi/cmd.data/uniqpair.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/uniqpair.c	(revision 40291)
+++ trunk/Ohana/src/opihi/cmd.data/uniqpair.c	(revision 40291)
@@ -0,0 +1,146 @@
+# include "data.h"
+
+/* given two 32bit IDs, a 64bit joint ID and sort it, along with an index (sequence) vector */ 
+int uniqpair (int argc, char **argv) {
+
+  int N;
+  Vector *ID1vec = NULL, *ID2vec = NULL, *OUTvec = NULL;
+
+  int VERBOSE = FALSE;
+  if ((N = get_argument (argc, argv, "-v"))) {
+    remove_argument (N, &argc, argv);
+    VERBOSE = TRUE;
+  }
+
+  Vector *SEQdup = NULL;
+  if ((N = get_argument (argc, argv, "-d"))) {
+    remove_argument (N, &argc, argv);
+    if ((SEQdup = SelectVector (argv[N], ANYVECTOR, TRUE)) == NULL) {
+      gprint (GP_ERR, "invalid vector %s\n", "SEQdups");
+      return FALSE;
+    }
+    remove_argument (N, &argc, argv);
+  }
+
+  Vector *CNTvec = NULL;
+  if ((N = get_argument (argc, argv, "-c"))) {
+    remove_argument (N, &argc, argv);
+    if ((CNTvec = SelectVector (argv[N], ANYVECTOR, TRUE)) == NULL) {
+      gprint (GP_ERR, "invalid vector %s\n", argv[N]);
+      return FALSE;
+    }
+    remove_argument (N, &argc, argv);
+  }
+
+  if (argc != 4) {
+    gprint (GP_ERR, "USAGE: uniqpair (ID1) (ID2) (IDout) [-c count] [-d dupindex]\n");
+    gprint (GP_ERR, "  merge ID1 and ID2 (32 bit int values) into a single index and find unique entries\n");
+    gprint (GP_ERR, "  -c count: save the number of each unique entry in the count vector\n");
+    gprint (GP_ERR, "  -d dupindex: save the sequence number of duplicate entries in the dupindex vector\n");
+    return (FALSE);
+  }
+
+  if ((ID1vec = SelectVector (argv[1], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((ID2vec = SelectVector (argv[2], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((OUTvec = SelectVector (argv[3], ANYVECTOR, TRUE)) == NULL) return (FALSE);
+
+  if ((ID1vec->type == OPIHI_FLT) || (ID2vec->type == OPIHI_FLT)) {
+    gprint (GP_ERR, "ERROR: ID vectors much be int\n");
+    return (FALSE);
+  }
+  if (ID1vec->Nelements != ID2vec->Nelements) {
+    gprint (GP_ERR, "ERROR: ID vectors not all the same length\n");
+    return (FALSE);
+  }
+
+  // storage for the output and count vectors
+  int D_NOUT = 0.02*ID1vec->Nelements;
+
+  int Nout = 0;
+  int NOUT = D_NOUT;
+  ResetVector (OUTvec, OPIHI_INT, NOUT);
+  if (CNTvec) {
+    ResetVector (CNTvec, OPIHI_INT, NOUT);
+  }
+
+  // storage for the duplicate sequences
+  int NseqDup = 0;
+  int NSEQDUP = 10000;
+  if (SEQdup) {
+    ResetVector (SEQdup, OPIHI_INT, NSEQDUP);
+  }
+  
+  // incrementing pointers to the input IDs
+  opihi_int *ID1 = ID1vec->elements.Int;
+  opihi_int *ID2 = ID2vec->elements.Int;
+
+  // generate the joined and sequence vectors
+  opihi_int *IDfull = NULL;
+  opihi_int *SEQdata = NULL;
+  ALLOCATE (IDfull, opihi_int, ID1vec->Nelements);
+  ALLOCATE (SEQdata, opihi_int, ID1vec->Nelements);
+  for (int i = 0; i < ID1vec->Nelements; i++) {
+    IDfull[i] = ID1[i] + (ID2[i] << 32);
+    SEQdata[i] = i;
+  }
+
+  // XXX watch out: this all needs to get rationalized...
+  llsortpair ((off_t *)IDfull, (off_t *)SEQdata, ID1vec->Nelements);
+
+  opihi_int *vtgt = OUTvec->elements.Int;
+  opihi_int *vsrc = IDfull;
+
+  int onePercent = ID1vec->Nelements / 100;
+
+  struct sigaction *old_sigaction = SetInterrupt();
+  for (int i = 0; (i < ID1vec->Nelements) && !interrupt; Nout++) {
+    if (Nout >= NOUT) {
+      NOUT += D_NOUT;
+      REALLOCATE (vtgt, opihi_int, NOUT);
+      if (CNTvec) {
+	REALLOCATE (CNTvec->elements.Int, opihi_int, NOUT);
+      }
+    }
+    vtgt[Nout] = *vsrc;
+    int Ndup = 0;
+    opihi_int lastValue = *vsrc;
+    while ((i < ID1vec->Nelements) && (*vsrc == lastValue)) {
+      i++;
+      vsrc ++;
+      Ndup ++;
+      if (VERBOSE && (i % onePercent == 0)) gprint (GP_ERR, ".");
+    }
+    if (CNTvec) {
+      CNTvec->elements.Int[Nout] = Ndup;
+    }
+    if (SEQdup && (Ndup > 1)) {
+      // we want to save SEQdata values for the duplicates
+      // SEQdata[i-Ndup] .. SEQdata[i-N
+      for (int j = i - Ndup; j < i; j++) {
+	SEQdup->elements.Int[NseqDup] = SEQdata[j];
+	NseqDup ++;
+	if (NseqDup >= NSEQDUP) {
+	  NSEQDUP += 10000;
+	  REALLOCATE (SEQdup->elements.Int, opihi_int, NSEQDUP);
+	}
+      }
+    }
+  }
+  ClearInterrupt (old_sigaction);
+  if (VERBOSE) gprint (GP_ERR, "\n");
+
+  free (SEQdata);
+  free (IDfull);
+
+  if (SEQdup) {
+    SEQdup->Nelements = NseqDup;
+    REALLOCATE (SEQdup->elements.Int, opihi_int, NseqDup);
+  }
+
+  // fix references and free up extra memory:
+  OUTvec->elements.Int = vtgt;
+  ResetVector (OUTvec, OPIHI_INT, Nout);
+  if (CNTvec) ResetVector (CNTvec, OPIHI_INT, Nout);
+
+  return (TRUE);
+}
Index: trunk/Ohana/src/opihi/cmd.data/write_vectors.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/write_vectors.c	(revision 40165)
+++ trunk/Ohana/src/opihi/cmd.data/write_vectors.c	(revision 40291)
@@ -165,7 +165,7 @@
 	} else {
 	  if (CSV) {
-	    fprintf (f, "%d,", vec[j][0].elements.Int[i]);
+	    fprintf (f, OPIHI_INT_FMT",", vec[j][0].elements.Int[i]);
 	  } else {
-	    fprintf (f, "%d ", vec[j][0].elements.Int[i]);
+	    fprintf (f, OPIHI_INT_FMT" ", vec[j][0].elements.Int[i]);
 	  }
 	}