Index: /branches/eam_branches/ipp-20130711/psphot/test/tap_psphot_galaxygrid.pro
===================================================================
--- /branches/eam_branches/ipp-20130711/psphot/test/tap_psphot_galaxygrid.pro	(revision 36028)
+++ /branches/eam_branches/ipp-20130711/psphot/test/tap_psphot_galaxygrid.pro	(revision 36029)
@@ -155,15 +155,4 @@
   $BaseConfig = $FakeConfig
 
-  # $FakeConfig = $FakeConfig -D GALAXY.MODEL PS_MODEL_GAUSS
-  # $FakeConfig = $FakeConfig -D GALAXY.MODEL PS_MODEL_EXP
-  # $FakeConfig = $FakeConfig -D GALAXY.MODEL PS_MODEL_SERSIC
-  # $FakeConfig = $FakeConfig -D GALAXY.MODEL PS_MODEL_DEV
-  # $FakeConfig = $FakeConfig -Df GALAXY.RMAJOR.MIN 10.0
-  # $FakeConfig = $FakeConfig -Df GALAXY.RMAJOR.MAX 10.0
-  # $FakeConfig = $FakeConfig -Df GALAXY.ARATIO.MIN 0.25
-  # $FakeConfig = $FakeConfig -Df GALAXY.ARATIO.MAX 0.25
-  # $FakeConfig = $FakeConfig -Df GALAXY.INDEX.MIN 1.66
-  # $FakeConfig = $FakeConfig -Df GALAXY.INDEX.MAX 1.66
-
   $Nseq = 0
   foreach type EXP DEV
@@ -229,5 +218,5 @@
         foreach fwhm 0.8 1.0 1.5
           sprint name "sample.%02d" $Nseq
-          cmf.load.concat $1/$name.dat $2/$name.fit.cmf DEVEXP
+          cmf.load.concat $1/$name.dat $2/$name.fit.cmf $type
 	  $Nseq ++
         end
@@ -280,4 +269,11 @@
 
 macro grid.mkexp.sersic
+  if ($0 != 2)
+    echo "USAGE: grid.mkexp.devexp (dir)"
+    break
+  end
+
+  mkdir $1
+
   $FakeConfig = -camera SIMTEST
   $FakeConfig = $FakeConfig -recipe PPSIM STACKTEST.RUN
@@ -288,4 +284,5 @@
   $FakeConfig = $FakeConfig -Db GALAXY.FAKE T                        ; # generate a "realistic" distribution of galaxies
   $FakeConfig = $FakeConfig -Df GALAXY.MAG 17.0
+  $FakeConfig = $FakeConfig -Df GALAXY.GRID.MAG 14.5
   $FakeConfig = $FakeConfig -Db GALAXY.GRID T                        ; # generate a grid of galaxies (constant mag)
   $FakeConfig = $FakeConfig -Df GALAXY.THETA.MIN 0 
@@ -295,13 +292,4 @@
   $FakeConfig = $FakeConfig -D GALAXY.MODEL PS_MODEL_SERSIC
   $BaseConfig = $FakeConfig
-
-  # $FakeConfig = $FakeConfig -Df GALAXY.INDEX.MIN 1.0
-  # $FakeConfig = $FakeConfig -Df GALAXY.INDEX.MAX 1.0
-  # $FakeConfig = $FakeConfig -Df GALAXY.RMAJOR.MIN 10.0
-  # $FakeConfig = $FakeConfig -Df GALAXY.RMAJOR.MAX 10.0
-  # $FakeConfig = $FakeConfig -Df GALAXY.ARATIO.MIN 0.25
-  # $FakeConfig = $FakeConfig -Df GALAXY.ARATIO.MAX 0.25
-  # $FakeConfig = $FakeConfig -Df GALAXY.INDEX.MIN 1.66
-  # $FakeConfig = $FakeConfig -Df GALAXY.INDEX.MAX 1.66
 
   $Nseq = 0
@@ -318,6 +306,28 @@
   	  $FakeConfig = $FakeConfig -Df GALAXY.INDEX.MAX $index
 	  
+          sprint name "$1/sersic.%02d" $Nseq
+          mkexp $name $fwhm SERSIC
+	  $Nseq ++
+        end
+      end
+    end
+  end
+end
+
+macro grid.fitexp.sersic.devexp
+  if ($0 != 3)
+    echo "USAGE: grid.fitexp.sersic.devexp (srcdir) (outdir)"
+    break
+  end
+
+  mkdir $2
+
+  $Nseq = 0
+  foreach index 1 2 3 4
+    foreach Rmajor 3 10 30
+      foreach Aratio 0.25 0.5 1.0
+        foreach fwhm 0.8 1.0 1.5
           sprint name "sersic.%02d" $Nseq
-          mkexp $name $fwhm SERSIC
+          fitexp $1/$name $2/$name.fit EXP_CONV,DEV_CONV
 	  $Nseq ++
         end
@@ -343,6 +353,11 @@
 
 macro grid.load.sersic
-  delete -q Xin_s Yin_s Min_s Tin_s Rin_s rin_s Iin_s
-  delete -q Xot_s Yot_s Mot_s Tot_s Rot_s rot_s Iot_s
+  if ($0 != 3)
+    echo "USAGE: grid.load.devexp (srcdir) (fitdir)"
+    break
+  end
+
+  delete -q Xin_s Yin_s Min_s Tin_s Rin_s rin_s MTin_s Iin_s
+  delete -q Xot_s Yot_s Mot_s Tot_s Rot_s rot_s MTot_s Iot_s
 
   $Nseq = 0
@@ -352,5 +367,5 @@
         foreach fwhm 0.8 1.0 1.5
           sprint name "sersic.%02d" $Nseq
-          cmf.load.concat $name.dat $name.fit.cmf SERSIC
+          cmf.load.concat $1/$name.dat $2/$name.fit.cmf SERSIC
 	  $Nseq ++
         end
@@ -359,4 +374,10 @@
   end
 end
+
+# I want to make plots of Iin_s vs Mkron, Mxx, and similar things
+# this means I need to be able to join Chip.xfit things against Chip.psf
+# and to join Chip.xfit(DEV) to Chip.xfit(EXP)
+
+# I think I need a generic 'JOIN' function
 
 macro grid.plots.sersic
@@ -482,10 +503,28 @@
   subset IndexIn = IndexIn_all if (Type == 1)
 
+  $TYPE_S = 83
+  $TYPE_D = 68
+  $TYPE_E = 69
+  if ("$3" == "SERSIC")
+    $InType = $TYPE_S
+  end
+  if ("$3" == "DEV")
+    $InType = $TYPE_D
+  end
+  if ("$3" == "EXP")
+    $InType = $TYPE_E
+  end
+
   data $2
-  if ("$3" == "SERSIC")
-    read -fits Chip.xfit X_EXT Y_EXT EXT_INST_MAG EXT_WIDTH_MAJ EXT_WIDTH_MIN EXT_THETA EXT_PAR_07
-  else
-    read -fits Chip.xfit X_EXT Y_EXT EXT_INST_MAG EXT_WIDTH_MAJ EXT_WIDTH_MIN EXT_THETA
-  end
+
+  break -auto off
+  read -fits Chip.xfit X_EXT Y_EXT EXT_INST_MAG EXT_WIDTH_MAJ EXT_WIDTH_MIN EXT_THETA MODEL_TYPE EXT_PAR_07
+  $reread = not($STATUS)
+  break -auto on
+  if ($reread)
+    read -fits Chip.xfit X_EXT Y_EXT EXT_INST_MAG EXT_WIDTH_MAJ EXT_WIDTH_MIN EXT_THETA MODEL_TYPE 
+    set EXT_PAR_07 = (MODEL_TYPE:9 == 68)*4 + (MODEL_TYPE:9 == 69)
+  end
+
   set EXT_THETA_ALT = EXT_THETA * (EXT_THETA >= 0.0) + (EXT_THETA + 3.14159265) * (EXT_THETA < 0.0)
   set EXT_THETA = EXT_THETA_ALT * 180 / 3.14159265
@@ -498,4 +537,7 @@
   reindex Xin_m = Xin using index2
   reindex Yin_m = Yin using index2
+
+  set MTin_m = $InType + zero(Xin_m)
+  reindex MTot_m = MODEL_TYPE:9 using index1
 
   reindex Mot_m = EXT_INST_MAG using index1
@@ -514,7 +556,7 @@
     reindex Iot_m = EXT_PAR_07 using index1
     reindex Iin_m = IndexIn using index2
-   $fields = X Y M T R r I
+   $fields = X Y M T R r MT I
   else
-   $fields = X Y M T R r
+   $fields = X Y M T R r MT
   end
   
@@ -523,4 +565,167 @@
       concat $field\$set\_m $field\$set\_s
     end
+  end
+
+  concat min min_S
+  concat Min Min_S
+end
+
+macro grid.load.sersic.test
+  if ($0 != 3)
+    echo "USAGE: grid.load.devexp.test (srcdir) (fitdir)"
+    break
+  end
+
+  $fields_bt = X Y M T R r I
+  $fields_ot = Pmag Kmag Amag IDx MTot
+
+  foreach field $fields_bt
+    foreach set in ot
+      delete -q $field\$set\_dev_s
+      delete -q $field\$set\_exp_s
+    end
+  end
+  foreach field $fields_ot
+    delete -q $field\_dev_s
+    delete -q $field\_exp_s
+  end
+
+  $Nseq = 0
+  foreach index 1 2 3 4
+    foreach Rmajor 3 10 30
+      foreach Aratio 0.25 0.5 1.0
+        foreach fwhm 0.8 1.0 1.5
+          sprint name "sersic.%02d" $Nseq
+          cmf.load.sersic.test $1/$name.dat $2/$name.fit.cmf
+	  $Nseq ++
+        end
+      end
+    end
+  end
+end
+
+# I have run DEV and EXP against input models of type SERSIC
+macro cmf.load.sersic.test
+  if ($0 != 3)
+    echo "USAGE: cmf.load.sersic.test (dat) (cmf)"
+    break
+  end
+
+  # input parameters
+  data $1
+  read Xin_all 1 Yin_all 2 Fin_all 3 Type 4 Min_all 5 RmajIn_all 7 RminIn_all 8 ThetaIn_all 9 IndexIn_all 10
+
+  # galaxies only
+  subset Xin = Xin_all if (Type == 1)
+  subset Yin = Yin_all if (Type == 1)
+  subset Min = Min_all if (Type == 1)
+  subset Fin = Fin_all if (Type == 1)
+  set min = -2.5*log(Fin)
+
+  subset Tin_rad = ThetaIn_all if (Type == 1)
+  set Tin = Tin_rad * 180 / 3.14159265
+
+  subset RmajIn = RmajIn_all if (Type == 1)
+  subset RminIn = RminIn_all if (Type == 1)
+
+  subset IndexIn = IndexIn_all if (Type == 1)
+
+  $TYPE_S = 83
+  $TYPE_D = 68
+  $TYPE_E = 69
+
+  data $2
+
+  # load measured values from xfit
+  read -fits Chip.xfit IPP_IDET X_EXT Y_EXT EXT_INST_MAG EXT_WIDTH_MAJ EXT_WIDTH_MIN EXT_THETA MODEL_TYPE PSF_INST_MAG AP_MAG KRON_MAG
+  set EXT_PAR_07 = (MODEL_TYPE:9 == $TYPE_D)*4 + (MODEL_TYPE:9 == $TYPE_E)
+
+  set EXT_THETA_ALT = EXT_THETA * (EXT_THETA >= 0.0) + (EXT_THETA + 3.14159265) * (EXT_THETA < 0.0)
+  set EXT_THETA = EXT_THETA_ALT * 180 / 3.14159265
+  set IPP_IDET_EXT = IPP_IDET
+  
+  match2d X_EXT Y_EXT Xin Yin 1.0 -index1 index1 -index2 index2
+
+  reindex Xot_m = X_EXT using index1
+  reindex Yot_m = Y_EXT using index1
+
+  reindex Xin_m = Xin using index2
+  reindex Yin_m = Yin using index2
+
+  reindex Mot_m = EXT_INST_MAG using index1
+  reindex Tot_m = EXT_THETA using index1
+
+  reindex Min_m = Min using index2
+  reindex Tin_m = Tin using index2
+
+  reindex Rot_m = EXT_WIDTH_MAJ using index1
+  reindex rot_m = EXT_WIDTH_MIN using index1
+
+  reindex Rin_m = RmajIn using index2
+  reindex rin_m = RminIn using index2
+  
+  reindex Pmag_m = PSF_INST_MAG using index1
+  reindex Kmag_m = KRON_MAG using index1
+  reindex Amag_m = AP_MAG using index1
+
+  reindex IDx_m = IPP_IDET_EXT using index1
+
+  reindex MTot_m = MODEL_TYPE:9 using index1
+
+  reindex Iot_m = EXT_PAR_07 using index1
+  reindex Iin_m = IndexIn using index2
+  
+  # load moments and other kron values from Chip.psf
+  read -fits Chip.psf IPP_IDET X_PSF Y_PSF MOMENTS_XX MOMENTS_XY MOMENTS_YY KRON_FLUX_INNER MOMENTS_R1 MOMENTS_RH
+  set IPP_IDET_PSF = IPP_IDET
+
+  join -outer IPP_IDET_PSF IDx_m 
+  reindex Xp = X_PSF using index2
+  reindex Yp = Y_PSF using index2
+
+  reindex Mxx_m = MOMENTS_XX using index2
+  reindex Mxy_m = MOMENTS_XY using index2
+  reindex Myy_m = MOMENTS_YY using index2
+  reindex Mr1_m = MOMENTS_R1 using index2
+  reindex Mrh_m = MOMENTS_RH using index2
+  reindex Kfi_m = KRON_FLUX_INNER using index2
+  set Kmi_m = -2.5*log(Kfi_m)
+
+  $fields_bt = X Y M T R r I
+  $fields_ot = Pmag Kmag Amag IDx MTot Mxx Mxy Myy Mr1 Mrh Kmi
+
+  foreach field $fields_ot
+    subset $field\_exp_m = $field\_m where (MTot_m == $TYPE_E)
+    subset $field\_dev_m = $field\_m where (MTot_m == $TYPE_D)
+  end
+  foreach field $fields_bt
+    foreach set in ot
+      subset $field\$set\_exp_m = $field\$set\_m where (MTot_m == $TYPE_E)
+      subset $field\$set\_dev_m = $field\$set\_m where (MTot_m == $TYPE_D)
+    end
+  end
+
+  join IDx_exp_m IDx_dev_m
+  foreach field $fields_ot
+    reindex $field\_exp_mr = $field\_exp_m using index1
+    reindex $field\_dev_mr = $field\_dev_m using index2
+  end
+  foreach field $fields_bt
+    foreach set in ot
+      reindex $field\$set\_exp_mr = $field\$set\_exp_m using index1
+      reindex $field\$set\_dev_mr = $field\$set\_dev_m using index2
+    end
+  end
+
+  # concat
+  foreach field $fields_bt
+    foreach set in ot
+      concat $field\$set\_dev_mr $field\$set\_dev_s
+      concat $field\$set\_exp_mr $field\$set\_exp_s
+    end
+  end
+  foreach field $fields_ot
+    concat $field\_dev_mr $field\_dev_s
+    concat $field\_exp_mr $field\_exp_s
   end