Index: branches/czw_branch/20160809/Ohana/src/opihi/cmd.data/test/fit1d_irls.sh
===================================================================
--- branches/czw_branch/20160809/Ohana/src/opihi/cmd.data/test/fit1d_irls.sh	(revision 39753)
+++ branches/czw_branch/20160809/Ohana/src/opihi/cmd.data/test/fit1d_irls.sh	(revision 39753)
@@ -0,0 +1,347 @@
+
+list tests
+ test1
+ test2
+ test3
+ test4
+ test5
+ test6
+ test7
+ test8
+end
+
+## NOTE: this is not a TAP test
+macro mean.test
+  if ($0 != 2)
+    echo "USAGE: mean.test (errorbar)"
+    break
+  end
+
+  create x 0 100
+
+  delete C0fit dC0fit Cnfit Mfit dMfit
+  for i 0 500
+    gaussdev y x[] 0 1.0   ; # mean of 0, sigma of 1
+    
+    set dy = 1.0 + zero(x) ; # set error to be correct
+    fit1d_irls -q x y 0 -dy dy -wt wt -mask mask -use-median -Nboot 1000
+
+    concat $C0    C0fit
+    concat $dC0  dC0fit
+
+    vstat -q y
+    concat $MEAN Mfit
+    concat {$SIGMA / sqrt($NPTS)} dMfit
+  end
+
+  vstat -q C0fit;  set MeanC  = $MEAN; set StdevC = $SIGMA
+  vstat -q dC0fit; set MeanDC = $MEAN
+  vstat -q Mfit;   set MeanM  = $MEAN; set StdevM = $SIGMA
+  vstat -q dMfit;  set MeanDM = $MEAN
+
+  fprintf "%7.3f %7.3f" $MeanC $MeanM
+  fprintf "%7.3f %7.3f" $StdevC $StdevM
+end
+
+macro cliptest.plot
+  if ($0 != 3)
+    echo "USAGE: cliptest.plot (errorbar) (output)"
+    break
+  end
+
+  set dy = $1 + zero(x) ; # set error to be correct
+  fit1d_irls x y 0 -dy dy -wt wt -mask mask -use-median
+
+  subset ygood = y where not(mask)
+
+  vstat y
+  vstat ygood
+
+  histogram wt Nwt 0 1 0.005 -range dw
+  vstat wt
+
+  subset wts = wt where (wt < 0.3*$MEDIAN)
+  echo "0.3 of MEDIAN : wts[]"
+
+  subset wts = wt where (wt < 0.3*$MEAN)
+  echo "0.3 of MEAN : wts[]"
+
+  dev -n 0; resize 1200 600
+  label -fn courier 24; 
+  section a 0.4 0.0 0.6 1.0; lim x y; clear; box -ypad 0.1 -xpad 3.2 -labelpadx 2.5 -labels 1001; label -x sequence +y value +x "red: clipped, blue: unclipped"; 
+  plot x y; plot -c blue -pt 7 x y where not(mask); plot -c red -pt 7 -lw 2 x y where mask
+
+  section b 0.0 0.0 0.4 1.0; lim dw Nwt; box +ypad 0.1 -xpad 3.2 -labelpadx 2.5; label -x "weight" -y "Npoints"; plot -x 1 dw Nwt
+  png -name $2
+end
+
+## NOTE: this is not a TAP test
+macro cliptest
+
+  create x 0 100
+  gaussdev y x[] 0 1.0   ; # mean of 0, sigma of 1
+
+  if (1)
+    y[10] = 5.0
+    y[20] = 7.0
+    y[30] = 4.0
+  end
+
+  cliptest.plot 1.0 irls.sample.v0.png
+
+  cliptest.plot 0.1 irls.sample.v1.png
+
+  cliptest.plot 0.01 irls.sample.v2.png
+end
+
+## NOTE: this is not a TAP test
+macro outlier.test
+  if ($0 != 2)
+    echo "USAGE: outlier.test (errorbar)"
+    break
+  end
+
+  create x 0 100
+
+  delete C0fit dC0fit Cnfit Mfit dMfit
+  for i 0 300
+    gaussdev y x[] 0 1.0   ; # mean of 0, sigma of 1
+    
+    if (0)
+      y[10] = 5.0
+      y[20] = 7.0
+      y[30] = 4.0
+    end
+    
+    set dy = $1 + zero(x) ; # set error to be correct
+    fit1d_irls -q x y 0 -dy dy -wt wt -mask mask -use-median -Nboot 1000
+    concat $C0 C0fit
+    concat $dC0 dC0fit
+    concat $Cnfit Cnfit
+
+    vstat -q y
+    concat $MEAN Mfit
+    concat {$SIGMA / sqrt($NPTS)} dMfit
+  end
+
+  vstat C0fit
+  vstat dC0fit
+  vstat Cnfit
+  vstat Mfit
+  vstat dMfit
+end
+
+# fit a line without errors
+macro test1
+ $PASS = 1
+ break -auto off
+
+ delete -q x y
+
+ create x 0 100
+ set y = 3 + 5*x
+ fit -q x y 1
+
+ if ($Cn != 1)
+   $PASS = 0
+ end
+ if (abs($C0 - 3) > 1e-5)
+   $PASS = 0
+ end
+ if (abs($C1 - 5) > 1e-5)
+   $PASS = 0
+ end
+end
+
+# fit a line with errors
+macro test2
+ $PASS = 1
+ break -auto off
+
+ delete -q x y dy
+
+ create x 0 100
+ set dy = 0.1*rnd(x) - 0.05
+ set y = 3 + 5*x + dy
+ fit -q x y 1
+
+ if ($Cn != 1)
+   $PASS = 0
+   echo "wrong number of elements : $Cn"
+ end
+ if (abs($C0 - 3) > 0.06)
+   $PASS = 0
+   echo "wrong value for C0 : $C0"
+ end
+ if (abs($C1 - 5) > 0.06)
+   $PASS = 0
+   echo "wrong value for C1 : $C1"
+ end
+end
+
+# fit a line with errors and weights
+macro test3
+ $PASS = 1
+ break -auto off
+
+ delete -q x y dy
+
+ create x 0 100
+ set dy = 0.1*rnd(x) - 0.05
+ set y = 3 + 5*x + dy
+ set dy = 0.1 + zero(x)
+ fit -q x y 1 -dy dy
+
+ if ($Cn != 1)
+   $PASS = 0
+ end
+ if (abs($C0 - 3) > 0.06)
+   $PASS = 0
+ end
+ if (abs($C1 - 5) > 0.06)
+   $PASS = 0
+ end
+end
+
+# fit a line with errors, weights, and outliers 
+macro test4
+ $PASS = 1
+ break -auto off
+
+ delete -q x y dy
+
+ create x 0 100
+ set dy = 0.1*rnd(x) - 0.05
+ set y = 3 + 5*x + dy
+ set dy = 0.1 + zero(x)
+ y[5] = 23
+ y[20] = -10
+ y[50] = 0.0
+ fit -q x y 1 -dy dy -clip 3 3
+
+ if ($Cn != 1)
+   $PASS = 0
+ end
+ if ($Cnv != 97)
+   $PASS = 0
+ end
+ if (abs($C0 - 3) > 0.06)
+   $PASS = 0
+ end
+ if (abs($C1 - 5) > 0.06)
+   $PASS = 0
+ end
+end
+
+# fit a quadratic without errors
+macro test5
+ $PASS = 1
+ break -auto off
+
+ delete -q x y
+
+ create x 0 100
+ set y = 3 + 5*x - 4*x^2
+ fit -q x y 2
+
+ if ($Cn != 2)
+   $PASS = 0
+ end
+ if (abs($C0 - 3) > 1e-5)
+   $PASS = 0
+ end
+ if (abs($C1 - 5) > 1e-5)
+   $PASS = 0
+ end
+ if (abs($C2 + 4) > 1e-5)
+   $PASS = 0
+ end
+end
+
+# fit a quadratic with errors
+macro test6
+ $PASS = 1
+ break -auto off
+
+ delete -q x y dy
+
+ create x 0 100
+ set dy = 0.1*rnd(x) - 0.05
+ set y = 3 + 5*x - 4*x^2 + dy
+ fit -q x y 2
+
+ if ($Cn != 2)
+   $PASS = 0
+ end
+ if (abs($C0 - 3) > 0.06)
+   $PASS = 0
+ end
+ if (abs($C1 - 5) > 0.06)
+   $PASS = 0
+ end
+ if (abs($C2 + 4) > 0.06)
+   $PASS = 0
+ end
+end
+
+# fit a quadratic with errors and weights
+macro test7
+ $PASS = 1
+ break -auto off
+
+ delete -q x y dy
+
+ create x 0 100
+ set dy = 0.1*rnd(x) - 0.05
+ set y = 3 + 5*x - 4*x^2 + dy
+ set dy = 0.1 + zero(x)
+ fit -q x y 2 -dy dy
+
+ if ($Cn != 2)
+   $PASS = 0
+ end
+ if (abs($C0 - 3) > 0.06)
+   $PASS = 0
+ end
+ if (abs($C1 - 5) > 0.06)
+   $PASS = 0
+ end
+ if (abs($C2 + 4) > 0.06)
+   $PASS = 0
+ end
+end
+
+# fit a quadratic with errors, weights, and outliers 
+macro test8
+ $PASS = 1
+ break -auto off
+
+ delete -q x y dy
+
+ create x 0 100
+ set dy = 0.1*rnd(x) - 0.05
+ set y = 3 + 5*x - 4*x^2 + dy
+ set dy = 0.1 + zero(x)
+ y[5] = 23
+ y[20] = -10
+ y[50] = 0.0
+
+ # it takes 4 iterations to successfully reject the outliers above...
+ fit -q x y 2 -dy dy -clip 3 4
+
+ if ($Cn != 2)
+   $PASS = 0
+ end
+ if ($Cnv != 97)
+   $PASS = 0
+ end
+ if (abs($C0 - 3) > 0.06)
+   $PASS = 0
+ end
+ if (abs($C1 - 5) > 0.06)
+   $PASS = 0
+ end
+ if (abs($C2 + 4) > 0.06)
+   $PASS = 0
+ end
+end