IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 33095


Ignore:
Timestamp:
Jan 11, 2012, 2:23:22 PM (15 years ago)
Author:
eugene
Message:

tap test to test variance modes

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/eam_branches/ipp-20111122/psphot/test/tap_psphot_varmodel.pro

    r33088 r33095  
    1 #!/usr/bin/env mana 
     1#!/usr/bin/env mana
    22# -*-sh-*-
    33
    4 # macros to support tap tests in mana/dvo
    5 input tap.pro
    6 
    7 # create a fake image (ppSim, all stars are fake)
    8 macro make.fake
    9 
    10   local options
    11 
    12   # basic options for these image (r, 100sec)
    13   $options = -type OBJECT -filter r -exptime 100.0
    14   $options = $options -skymags 20.86 -ra 270.70 -dec -23.70 -pa 0.0 -seeing 1.0
    15   $options = $options -D PSF.MODEL PS_MODEL_GAUSS
     4# $KAPA = kapa -noX
     5
     6# PSF.CONVOLVE : if true, we insert delta functions (and optionally
     7#                galaxies) and smooth the image with the psf model
     8#                (uses a GAUSS regardless of the model). Note that
     9#                PSF.CONVOLVE = T is faster than F, but (a) only
     10#                allows Gauss models and (b) only yields quantized
     11#                locations
     12
     13# config for ppImage to generate chip, mask, weight
     14$ppImageConfig = -recipe PPIMAGE PPIMAGE_N
     15$ppImageConfig = $ppImageConfig -Db BACKGROUND T
     16$ppImageConfig = $ppImageConfig -Db CHIP.FITS T
     17$ppImageConfig = $ppImageConfig -Db CHIP.MASK.FITS T
     18$ppImageConfig = $ppImageConfig -Db CHIP.VARIANCE.FITS T
     19$ppImageConfig = $ppImageConfig -Db BASE.FITS F
     20$ppImageConfig = $ppImageConfig -Db VARIANCE.BUILD T
     21$ppImageConfig = $ppImageConfig -Db PHOTOM F
     22
     23# basic options for the these images (filter, location, obstype)
     24$BaseOptions = -type OBJECT -filter r -skymags 20.86 -ra 270.70 -dec -23.70 -pa 0.0
     25$BaseOptions = $BaseOptions -Df PSASTRO:DVO.GETSTAR.MAX.RHO 50000.0
     26
     27# create an image with fake sources and insert the resulting cmf file into a dvodb
     28$RefConfig = -camera SIMTEST
     29$RefConfig = $RefConfig -recipe PPSIM STACKTEST.MAKE
     30$RefConfig = $RefConfig -D PSASTRO:PSASTRO.CATDIR catdir.ref
     31$RefConfig = $RefConfig -Db PSF.CONVOLVE F
     32
     33# options for the reference image
     34$RefOptions = $BaseOptions -exptime 100.0
     35$RefOptions = $RefOptions -seeing 1.0
     36$RefOptions = $RefOptions -D PSF.MODEL PS_MODEL_GAUSS
     37$RefOptions = $RefOptions -Df STARS.DENSITY 10.0
     38$RefOptions = $RefOptions -Df STARS.SIGMA.LIM 0.5
     39
     40# basic config for ppSim with randomly distributed stars and NO galaxies
     41$RealConfig = -camera SIMTEST
     42$RealConfig = $RealConfig -recipe PPSIM STACKTEST.RUN
     43$RealConfig = $RealConfig -D PSASTRO:PSASTRO.CATDIR catdir.ref
     44$RealConfig = $RealConfig -Db STARS.FAKE F
     45$RealConfig = $RealConfig -Db STARS.REAL T
     46$RealConfig = $RealConfig -Db MATCH.DENSITY F
     47$RealConfig = $RealConfig -Db PSF.CONVOLVE F
     48$RealConfig = $RealConfig -Df STARS.DENSITY 10.0
     49$RealConfig = $RealConfig -Df STARS.SIGMA.LIM 2.5
     50$RealConfig = $RealConfig -Db GALAXY.FAKE F
     51$RealConfig = $RealConfig -Db GALAXY.GRID F
     52
     53# options for the repeated images
     54$RealOptions = $BaseOptions -exptime 30.0
    1655 
    17   $options = $options -Df STARS.DENSITY 20.0
    18 
    19   $options = $options -camera SIMTEST -recipe PPSIM STACKTEST.MAKE
    20   exec ppSim $options testimage
    21 end
    22 
    23 macro run.psphot
    24   exec psphot -file testimage.fits testimage.phot
    25 end
    26 
    27 macro test.init
    28   data testimage.phot.cmf
    29 
    30   # a bit hackish : cannot append to non-existent scalar
    31   $fields = Chip.psf
    32   for i 0 $pairs:n
    33     list word -split $pairs:$i
    34     $fields = $fields $word:0
    35   end
    36 
    37   read -fits $fields
    38 
    39   data testimage.dat
    40   read Xraw 1 Yraw 2 Fraw 3 Galaxy 4 Mraw 5 dMraw 6 SXXraw 7 SYYraw 8 SXYraw 9 S7raw 10
    41   set dFraw = Fraw * dMraw
    42   set dXraw = 4.0 * dMraw / 2.35
    43   set dYraw = 4.0 * dMraw / 2.35
    44 
    45   match2d X_PSF Y_PSF Xraw Yraw 2.0 -index1 index1 -index2 index2
    46 end
    47 
    48 macro test.plots
    49   if ($0 != 2)
    50     echo "USAGE: test.plots (outroot)"
    51     break
    52   end
    53 
    54   clear -s
    55 
    56   resize 1200 1175
    57 
     56$ExtraOptions = -D PSF.MODEL PS_MODEL_GAUSS
     57
     58# sample alternate options:
     59# $ppSimOptions = $FakeOptions -D PSF.MODEL PS_MODEL_PS1_V1
     60# $ppSimOptions = $FakeOptions -Df PSF.ARATIO 1.2
     61# $ppSimOptions = $FakeOptions -Df PSF.THETA +30.0
     62# $ppSimOptions = $FakeOptions -D PSF.MODEL PS_MODEL_GAUSS
     63
     64list fwhm
     65 1.0
     66 1.1
     67 1.2
     68 1.5
     69end
     70
     71macro go
     72
     73  mkdir test
     74
     75  $ExtraOptions = -D PSF.MODEL PS_MODEL_GAUSS
     76  mkexp test/image.00 1.0
     77  $ExtraOptions = -D PSF.MODEL PS_MODEL_PGAUSS
     78  mkexp test/image.01 1.0
     79  $ExtraOptions = -D PSF.MODEL PS_MODEL_PS1_V1
     80  mkexp test/image.02 1.0
     81end
     82
     83# create a reference database of fake stars to be used by ppSim below
     84macro mkref
     85  exec rm -rf catdir.ref
     86  exec rm -f refimage.fits
     87 
     88  exec time ppSim $RefConfig $RefOptions refimage -nx 3000 -ny 3000
     89 
     90  file synth.photcodes found
     91  if (not($found))
     92    echo "making photcodes file"
     93    mkphotcodes synth.photcodes
     94  end
     95
     96  exec addstar -D CAMERA simtest -D CATDIR catdir.ref -accept-astrom -photcode SYNTH.r -D PHOTCODE_FILE synth.photcodes refimage.cmf
     97  exec relphot -averages -D CATDIR catdir.ref -update -region 260 280 -33 -13
     98end
     99
     100# create a realistic distribution of fake stars, GAUSS PSF
     101macro mkexp
     102  if ($0 != 3)
     103    echo "USAGE: mkexp basename fwhm"
     104    break
     105  end
     106
     107  local fwhm basename
     108  $basename = $1
     109  $fwhm = $2
     110
     111  # create the raw image
     112  echo ppSim -seeing $fwhm -nx 3000 -ny 3000 $RealConfig $RealOptions $ExtraOptions $basename
     113  exec ppSim -seeing $fwhm -nx 3000 -ny 3000 $RealConfig $RealOptions $ExtraOptions $basename
     114  exec /bin/mv -f $basename.cmf $basename.in.cmf
     115
     116  # create the chip output
     117  echo ppImage $ppImageConfig -file $basename.fits $basename
     118  exec ppImage $ppImageConfig -file $basename.fits $basename
     119end
     120
     121macro runphot
     122  if ($0 != 4)
     123    echo "USAGE: runphot basename outname options"
     124    break
     125  end
     126
     127  local basename
     128  $basename = $1
     129  $outname = $2
     130  $options = $3
     131
     132  # create the chip output
     133  echo psphot -threads 4 -file $basename.ch.fits -mask $basename.ch.mk.fits -variance $basename.ch.wt.fits $outname $options
     134  exec psphot -threads 4 -file $basename.ch.fits -mask $basename.ch.mk.fits -variance $basename.ch.wt.fits $outname $options
     135end
     136
     137# compare two cmf files with extname Chip.psf
     138# things to compare:
     139# * completeness (which sources in (1) are not detected in (2)
     140# * positions (X_PSF, Y_PSF)
     141# * instrumental psf mags
     142# * position errors (no input errors; use a model?)
     143# * measured FWHM?
     144# * kron mags (fluxes)
     145# * etc, etc
     146macro ckchip
     147  if ($0 != 5)
     148    echo "USAGE: ckchip (raw) (out) (output) (zpt_off)"
     149    break
     150  end
     151
     152  load.cmf $1 Chip.psf raw
     153  load.cmf $2 Chip.psf out
     154
     155  # images generated with convolution will not have the right output positions
     156  set X_raw = int(X_PSF_raw) + 0.5
     157  set Y_raw = int(Y_PSF_raw) + 0.5
     158  set M_raw = PSF_INST_MAG_raw + $4
     159  set K_out = -2.5*log(KRON_FLUX_out)
     160  match2d X_PSF_out Y_PSF_out X_PSF_raw Y_PSF_raw 1.5 -index1 index1 -index2 index2
     161
     162  local i NX NY nx ny N
     163
     164  device -n compare
     165  resize 1000 1000
     166
     167  # plot trends as a function of mag
     168  $NX = 2
     169  $NY = 5
    58170  $nx = 0
    59171  $ny = 0
    60   $np = 0
     172  $N = 0
     173  clear -s
    61174  for i 0 $pairs:n
    62     section a$i {0.5*$nx} {0.8 - 0.2*$ny} 0.5 0.2
     175    section a$nx\$ny {$nx/$NX} {$ny/$NY} {1/$NX} {1/$NY}
    63176    show.pair $i
    64177    $ny ++
    65     if ($ny == 5)
     178    if ($ny == $NY)
     179      $ny = 0
    66180      $nx ++
     181    end
     182    if ($nx == $NX)
     183      png -name $3.$N.png
     184      clear -s
     185      $nx = 0
    67186      $ny = 0
     187      $N ++
    68188    end
    69     if ($nx == 2)
    70       png -name $1.$np.png
     189  end
     190
     191  # plot (input - output) vs mag
     192end
     193
     194macro stchip
     195  if ($0 != 5)
     196    echo "USAGE: stchip (raw) (out) (output) (zpt_off)"
     197    break
     198  end
     199
     200  load.cmf $1 Chip.psf raw
     201  load.cmf $2 Chip.psf out
     202
     203  # images generated with convolution will not have the right output positions
     204  set X_raw = int(X_PSF_raw) + 0.5
     205  set Y_raw = int(Y_PSF_raw) + 0.5
     206  set M_raw = PSF_INST_MAG_raw + $4
     207  set K_out = -2.5*log(KRON_FLUX_out)
     208  match2d X_PSF_out Y_PSF_out X_PSF_raw Y_PSF_raw 1.5 -index1 index1 -index2 index2
     209
     210  local i
     211
     212  for i 0 $spairs:n
     213    stats.pair $i $3
     214  end
     215end
     216
     217# compare chip to warp
     218macro ckwarp
     219  if ($0 != 5)
     220    echo "USAGE: ckwarp (raw) (out) (output) (zpt_off)"
     221    break
     222  end
     223
     224  load.cmf $1 Chip.psf raw
     225  load.cmf $2 SkyChip.psf out
     226  set X_raw = int(X_PSF_raw) + 0.5
     227  set Y_raw = int(Y_PSF_raw) + 0.5
     228  set M_raw = PSF_INST_MAG_raw + $4
     229  set K_out = -2.5*log(KRON_FLUX_out)
     230  match2d X_PSF_out Y_PSF_out X_PSF_raw Y_PSF_raw 1.0 -index1 index1 -index2 index2
     231
     232  local i nx ny NX NY N
     233
     234  device -n compare
     235  resize 1000 1000
     236
     237  # plot trends as a function of mag
     238  $NX = 2
     239  $NY = 5
     240  $nx = 0
     241  $ny = 0
     242  $N = 0
     243  clear -s
     244
     245  for i 0 $pairs:n
     246    section a$nx\$ny {$nx/$NX} {$ny/$NY} {1/$NX} {1/$NY}
     247    show.pair $i
     248    $ny ++
     249    if ($ny == $NY)
     250      $ny = 0
     251      $nx ++
     252    end
     253    if ($nx == $NX)
     254      png -name $3.$N.png
    71255      clear -s
    72       $np ++
    73256      $nx = 0
     257      $ny = 0
     258      $N ++
    74259    end
    75   end 
    76   png -name $1.$np.png
     260  end
     261end
     262
     263macro stwarp
     264  if ($0 != 5)
     265    echo "USAGE: stwarp (raw) (out) (output) (zpt_off)"
     266    break
     267  end
     268
     269  load.cmf $1 Chip.psf raw
     270  load.cmf $2 SkyChip.psf out
     271
     272  # images generated with convolution will not have the right output positions
     273  set X_raw = int(X_PSF_raw) + 0.5
     274  set Y_raw = int(Y_PSF_raw) + 0.5
     275  set M_raw = PSF_INST_MAG_raw + $4
     276  set K_out = -2.5*log(KRON_FLUX_out)
     277  match2d X_PSF_out Y_PSF_out X_PSF_raw Y_PSF_raw 1.5 -index1 index1 -index2 index2
     278
     279  local i
     280
     281  for i 0 $spairs:n
     282    stats.pair $i $3
     283  end
     284end
     285
     286macro stats.pair
     287  if ($0 != 3)
     288    echo "USAGE: stats.pair (N) (output)"
     289    break
     290  end
     291
     292  list word -split $spairs:$1
     293  if ($word:n != 8)
     294    echo "invalid pair $1"
     295    break
     296  end
     297
     298  $Nr = $word:3
     299
     300  reindex v1 = $word:0 using index1
     301  reindex v2 = $word:1 using index2
     302  reindex rv = $word:2 using index$Nr
     303
     304  set delta = v1 - v2
     305  subset d1 = delta if ($word:4 < rv) && (rv < $word:5) && (abs(delta) < $word:7)
     306  subset d2 = delta if ($word:5 < rv) && (rv < $word:6) && (abs(delta) < $word:7)
     307
     308  vstats -q d1 -sigma-clip 3.0
     309  $M1 = $MEAN
     310  $S1 = $SIGMA
     311  vstats -q d2 -sigma-clip 3.0
     312  $M2 = $MEAN
     313  $S2 = $SIGMA
     314
     315  output $2
     316  fprintf "%-18s  %7.4f %7.4f  %7.4f %7.4f" $word:0  $M1 $S1  $M2 $S2
     317  output stdout
    77318end
    78319
     
    106347    lim rv $word:4 $word:5; box; plot rv delta
    107348  end
    108   label -y $word:0 -x $word:3
     349  label -y '$word:0' -x '$word:2'
     350end
     351
     352# This list is used to compare a pair of vectors (sans error) or a
     353# vector and an expected (constant) value.  The last field defines a
     354# vector or constant for the comparison.  It is assumed that the
     355# vector sets have been loaded and matched with match2d to generate
     356# index vectors 'index1' and index2'.  The macro 'show.pair' generates
     357# a plot of the range vector vs (v1 - v2).  The indices for v1, v2 are
     358# index1 and 2 respectively.  The index for the range vector is defined
     359# by the integer following that vector.  the y-limits of the plot are
     360# given by the last two numbers
     361list pairs
     362  X_PSF_out             X_raw                 PSF_INST_MAG_raw 2 -1.01 1.01 V
     363  Y_PSF_out             Y_raw                 PSF_INST_MAG_raw 2 -1.01 1.01 V
     364  X_PSF_out             X_PSF_raw             PSF_INST_MAG_raw 2 -1.01 1.01 V
     365  Y_PSF_out             Y_PSF_raw             PSF_INST_MAG_raw 2 -1.01 1.01 V
     366  X_PSF_SIG_out         X_PSF_SIG_raw         PSF_INST_MAG_raw 2 -1.01 1.01 V
     367  Y_PSF_SIG_out         Y_PSF_SIG_raw         PSF_INST_MAG_raw 2 -1.01 1.01 V
     368  #PSF_INST_MAG_out      PSF_INST_MAG_raw      PSF_INST_MAG_raw 2 -1.01 1.01 V
     369  PSF_INST_MAG_out      M_raw                 PSF_INST_MAG_raw 2 -1.01 1.01 V
     370  PSF_INST_MAG_SIG_out  PSF_INST_MAG_SIG_raw  PSF_INST_MAG_raw 2 -1.01 1.01 V
     371  #PSF_INST_FLUX_out     PSF_INST_FLUX_raw     PSF_INST_MAG_raw 2   def  def V
     372  #PSF_INST_FLUX_SIG_out PSF_INST_FLUX_SIG_raw PSF_INST_MAG_raw 2   def  def V
     373  AP_MAG_out            M_raw                 PSF_INST_MAG_raw 2 -1.01 1.01 V
     374  AP_MAG_RAW_out        M_raw                 PSF_INST_MAG_raw 2 -1.01 1.01 V
     375  AP_MAG_RADIUS_out     0.0                   PSF_INST_MAG_raw 2 -0.01 20.1 S
     376  SKY_out               0.0                   PSF_INST_MAG_raw 2   def  def S
     377  SKY_SIGMA_out         0.0                   PSF_INST_MAG_raw 2   def  def S
     378  PSF_CHISQ_out         1.0                   PSF_INST_MAG_raw 2   def  def S
     379  CR_NSIGMA_out         0.0                   PSF_INST_MAG_raw 2   def  def S
     380  EXT_NSIGMA_out        0.0                   PSF_INST_MAG_raw 2 -5.01 5.01 S
     381  PSF_MAJOR_out         0.0                   PSF_INST_MAG_raw 2 -0.01 5.01 S
     382  PSF_MINOR_out         0.0                   PSF_INST_MAG_raw 2 -0.01 5.01 S
     383  PSF_THETA_out         0.0                   PSF_INST_MAG_raw 2 -1.61 1.61 S
     384  PSF_QF_out            0.0                   PSF_INST_MAG_raw 2 -0.10 1.10 S
     385  PSF_QF_PERFECT_out    0.0                   PSF_INST_MAG_raw 2 -0.10 1.10 S
     386  PSF_NDOF_out          0.0                   PSF_INST_MAG_raw 2  def  def  S
     387  PSF_NPIX_out          0.0                   PSF_INST_MAG_raw 2  def  def  S
     388  MOMENTS_XX_out        0.0                   PSF_INST_MAG_raw 2 -0.01 3.01 S
     389  MOMENTS_XY_out        0.0                   PSF_INST_MAG_raw 2 -3.01 3.01 S
     390  MOMENTS_YY_out        0.0                   PSF_INST_MAG_raw 2 -0.01 3.01 S
     391  MOMENTS_M3C_out       0.0                   PSF_INST_MAG_raw 2 -3.01 3.01 S
     392  MOMENTS_M3S_out       0.0                   PSF_INST_MAG_raw 2 -3.01 3.01 S
     393  MOMENTS_M4C_out       0.0                   PSF_INST_MAG_raw 2 -2.01 2.01 S
     394  MOMENTS_M4S_out       0.0                   PSF_INST_MAG_raw 2 -2.01 2.01 S
     395  MOMENTS_R1_out        0.0                   PSF_INST_MAG_raw 2 -5.01 5.01 S
     396  MOMENTS_RH_out        0.0                   PSF_INST_MAG_raw 2 -5.01 5.01 S
     397  K_out                 M_raw                 PSF_INST_MAG_raw 2  def  def  V
     398  KRON_FLUX_ERR_out     0.0                   PSF_INST_MAG_raw 2  def  def  S
     399  KRON_FLUX_INNER_out   0.0                   PSF_INST_MAG_raw 2  def  def  S
     400  KRON_FLUX_OUTER_out   0.0                   PSF_INST_MAG_raw 2  def  def  S
     401# CAL_PSF_MAG          CAL_PSF_MAG          none Mraw 2 -1.01 1.01 V
     402# CAL_PSF_MAG_SIG      CAL_PSF_MAG_SIG      none Mraw 2 -1.01 1.01 V
     403# RA_PSF               RA_PSF               none Mraw 2 -1.01 1.01 V
     404# DEC_PSF              DEC_PSF              none Mraw 2 -1.01 1.01 V
     405# PEAK_FLUX_AS_MAG     PEAK_FLUX_AS_MAG     none Mraw 2 -1.01 1.01 V
     406# FLAGS                FLAGS                0.0  Mraw 2 -1.01 1.01 S
     407# FLAGS2               FLAGS2               0.0  Mraw 2 -1.01 1.01 S
     408end
     409
     410macro load.cmf
     411  if ($0 != 4)
     412   echo "load.cmf (filename) (ext) (label)"
     413   break
     414  end
     415
     416  data $1
     417
     418  # create the list of fields to load
     419  delete -q myFields
     420  for i 0 $fields:n
     421    if ($?myFields)
     422      $myFields = $myFields $fields:$i
     423    else
     424      $myFields = $fields:$i
     425    end
     426  end
     427
     428  read -fits $2 $myFields
     429
     430  # rename the loaded vectors appending the supplied lable
     431  for i 0 $fields:n
     432    set $fields:$i\_$3 = $fields:$i
     433    delete $fields:$i
     434  end
     435end
     436
     437# this list defines the fields to be loaded from file
     438list fields
     439  X_PSF             
     440  Y_PSF             
     441  X_PSF_SIG         
     442  Y_PSF_SIG         
     443  PSF_INST_MAG       
     444  PSF_INST_MAG_SIG   
     445  PSF_INST_FLUX     
     446  PSF_INST_FLUX_SIG 
     447  AP_MAG             
     448  AP_MAG_RAW         
     449  AP_MAG_RADIUS     
     450  SKY               
     451  SKY_SIGMA         
     452  PSF_CHISQ         
     453  CR_NSIGMA         
     454  EXT_NSIGMA         
     455  PSF_MAJOR         
     456  PSF_MINOR         
     457  PSF_THETA         
     458  PSF_QF             
     459  PSF_QF_PERFECT     
     460  PSF_NDOF           
     461  PSF_NPIX           
     462  MOMENTS_XX         
     463  MOMENTS_XY         
     464  MOMENTS_YY         
     465  MOMENTS_M3C       
     466  MOMENTS_M3S       
     467  MOMENTS_M4C       
     468  MOMENTS_M4S       
     469  MOMENTS_R1         
     470  MOMENTS_RH         
     471  KRON_FLUX         
     472  KRON_FLUX_ERR     
     473  KRON_FLUX_INNER   
     474  KRON_FLUX_OUTER   
     475#   CAL_PSF_MAG     
     476#   CAL_PSF_MAG_SIG 
     477#   RA_PSF           
     478#   DEC_PSF         
     479#   PEAK_FLUX_AS_MAG
     480#  FLAGS           
     481#  FLAGS2           
     482end
     483
     484# use these cmf entries to measure average stats of the given pairs
     485list spairs
     486  X_PSF_out             X_raw                 PSF_INST_MAG_raw 2 -15 -10 -8.5  0.1
     487  Y_PSF_out             Y_raw                 PSF_INST_MAG_raw 2 -15 -10 -8.5  0.1
     488  PSF_INST_MAG_out      M_raw                 PSF_INST_MAG_raw 2 -15 -10 -8.5  0.1
     489  AP_MAG_out            M_raw                 PSF_INST_MAG_raw 2 -15 -10 -8.5  0.1
     490end
     491
     492#  # XXX this is a hack: the cmf file created by ppSim is not
     493#  # compatible with the image for warp (because one is a chip-mosaic
     494#  # and the other is not)
     495#  exec rm -f fix.hdr
     496#  output fix.hdr
     497#  echo "PSMOSAIC= 'CHIP    '           / Mosaicked level"
     498#  output stdout
     499# fix the header to be compatible with the chip file (to avoid running psphot)
     500# exec fits_insert $basename.cmf fix.hdr
     501
     502macro completeness
     503 if ($0 != 4)
     504   echo "USAGE: completeness (raw) (out) (output)"
     505   break
     506 end
     507
     508 load.cmf $1 Chip.psf raw
     509 load.cmf $2 Chip.psf out
     510 set X_raw = int(X_PSF_raw) + 0.5
     511 set Y_raw = int(Y_PSF_raw) + 0.5
     512 match2d X_PSF_raw Y_PSF_raw X_PSF_out Y_PSF_out 1.5 -index1 index1 -index2 index2 -closest
     513
     514 histogram PSF_INST_MAG_raw nMag -16.0 -3.0 0.25 -range dMag
     515 set fMag = zero(dMag)
     516 for i 0 {dMag[]-1}
     517  set inrange = (PSF_INST_MAG_raw > dMag[$i]) && (PSF_INST_MAG_raw <= dMag[$i+1])
     518  subset all = index1 if (inrange)
     519  subset got = index1 if (inrange) && (index1 >= 0)
     520  if (all[] == 0)
     521    fMag[$i] = 0
     522  else
     523   fMag[$i] = got[] / all[]
     524  end
     525 end
     526
     527 device -n complete
     528 resize 1000 600
     529
     530 clear -s
     531
     532 section default 0 0 1 1
     533 lim dMag fMag; clear; box -ypad 5 +ypad 5 -ticks 1110; plot -x 1 dMag fMag
     534 label -x mag_inst -y det_frac
     535
     536 set found = (index1 >= 0)
     537 plot PSF_INST_MAG_raw found -c red
     538
     539 set ldmag = log(PSF_INST_MAG_SIG_raw)
     540 section overlay 0 0 1 1; lim dMag -5 1.2; box -ypad 5 +ypad 5 -ticks 1011 -labels 1001; plot PSF_INST_MAG_raw ldmag
     541 label +y log(S/N)
     542
     543 png -name $3
    109544end
    110545
     
    134569    lim rv $word:5 $word:6; box; plot rv delta
    135570  end
    136   label -y $word:0 -x $word:3
    137 end
    138 
    139 macro show.pair.old
    140   set delta = $1 - $2
    141   lim mag2 delta; clear; box; plot mag2 delta
    142 end
    143 
    144 # list the entries that should be compared and the vector to use as a range
     571  label -y '$word:0' -x '$word:3'
     572end
     573
     574# this list is used to compare a pair of vectors with an error it is
     575# assumed that the vector sets have been loaded and matched with
     576# match2d to generate index vectors 'index1' and index2'.  the macro
     577# show.dpair generates a plot of the range vector vs (v1 - v2) / dv1.
     578# The indices for v1, dv1, and v2 are index1,1, and 2 respectively.  The
     579# index for the range vector is defined by the integer following that
     580# vector.  The y-limits of the plot are given by the last two numbers
     581# (use 'def') for the full default range of the delta vector
    145582list dpairs
    146   X_PSF             X_PSF_SIG         Xraw Mraw 2 -10.0 10.0
    147   Y_PSF             Y_PSF_SIG         Yraw Mraw 2 -10.0 10.0
    148   PSF_INST_MAG      PSF_INST_MAG_SIG  Mraw Mraw 2 -10.0 10.0
    149   PSF_INST_FLUX     PSF_INST_FLUX_SIG Fraw Mraw 2 -10.0 10.0
    150   AP_MAG            PSF_INST_MAG_SIG  Mraw Mraw 2 -10.0 10.0
    151   AP_MAG_RAW        PSF_INST_MAG_SIG  Mraw Mraw 2 -10.0 10.0
    152 end
    153 
    154 # list the entries that should be compared and the vector to use as a range
    155 list pairs
    156   X_PSF              Xraw Mraw 2 -1.01 1.01 V
    157   Y_PSF              Yraw Mraw 2 -1.01 1.01 V
    158   X_PSF_SIG         dXraw Mraw 2 -1.01 1.01 V
    159   Y_PSF_SIG         dYraw Mraw 2 -1.01 1.01 V
    160   PSF_INST_MAG       Mraw Mraw 2 -1.01 1.01 V
    161   PSF_INST_MAG_SIG    0.0 Mraw 2 -0.01 1.01 S
    162   PSF_INST_FLUX      Fraw Mraw 2   def  def V
    163   PSF_INST_FLUX_SIG dFraw Mraw 2   def  def V
    164   AP_MAG             Mraw Mraw 2 -1.01 1.01 V
    165   AP_MAG_RAW         Mraw Mraw 2 -1.01 1.01 V
    166   AP_MAG_RADIUS      0.0  Mraw 2 -0.01 20.1 S
    167   SKY                0.0  Mraw 2   def  def S
    168   SKY_SIGMA          0.0  Mraw 2   def  def S
    169   PSF_CHISQ          1.0  Mraw 2   def  def S
    170   CR_NSIGMA          0.0  Mraw 2   def  def S
    171   EXT_NSIGMA         0.0  Mraw 2 -5.01 5.01 S
    172   PSF_MAJOR          0.0  Mraw 2 -0.01 5.01 S
    173   PSF_MINOR          0.0  Mraw 2 -0.01 5.01 S
    174   PSF_THETA          0.0  Mraw 2 -1.61 1.61 S
    175   PSF_QF             0.0  Mraw 2 -0.10 1.10 S
    176   PSF_QF_PERFECT     0.0  Mraw 2 -0.10 1.10 S
    177   PSF_NDOF           0.0  Mraw 2  def  def  S
    178   PSF_NPIX           0.0  Mraw 2  def  def  S
    179   MOMENTS_XX         0.0  Mraw 2 -0.01 3.01 S
    180   MOMENTS_XY         0.0  Mraw 2 -3.01 3.01 S
    181   MOMENTS_YY         0.0  Mraw 2 -0.01 3.01 S
    182   MOMENTS_M3C        0.0  Mraw 2 -3.01 3.01 S
    183   MOMENTS_M3S        0.0  Mraw 2 -3.01 3.01 S
    184   MOMENTS_M4C        0.0  Mraw 2 -2.01 2.01 S
    185   MOMENTS_M4S        0.0  Mraw 2 -2.01 2.01 S
    186   MOMENTS_R1         0.0  Mraw 2 -5.01 5.01 S
    187   MOMENTS_RH         0.0  Mraw 2 -5.01 5.01 S
    188   KRON_FLUX          0.0  Mraw 2  def  def  S
    189   KRON_FLUX_ERR      0.0  Mraw 2  def  def  S
    190   KRON_FLUX_INNER    0.0  Mraw 2  def  def  S
    191   KRON_FLUX_OUTER    0.0  Mraw 2  def  def  S
    192 #   CAL_PSF_MAG      none Mraw 2 -1.01 1.01 V
    193 #   CAL_PSF_MAG_SIG  none Mraw 2 -1.01 1.01 V
    194 #   RA_PSF           none Mraw 2 -1.01 1.01 V
    195 #   DEC_PSF          none Mraw 2 -1.01 1.01 V
    196 #   PEAK_FLUX_AS_MAG none Mraw 2 -1.01 1.01 V
    197 #  FLAGS             0.0  Mraw 2 -1.01 1.01 S
    198 #  FLAGS2             0.0  Mraw 2 -1.01 1.01 S
    199 end
    200 
    201 if ($SCRIPT)
    202   make.fake
     583  # v1              dv                v2   range
     584  X_PSF             X_PSF_SIG         Xraw Mraw  2 -10.0 10.0
     585  Y_PSF             Y_PSF_SIG         Yraw Mraw  2 -10.0 10.0
     586  PSF_INST_MAG      PSF_INST_MAG_SIG  Mraw Mraw  2 -10.0 10.0
     587  PSF_INST_FLUX     PSF_INST_FLUX_SIG Fraw Mraw  2 -10.0 10.0
     588  AP_MAG            PSF_INST_MAG_SIG  Mraw Mraw  2 -10.0 10.0
     589  AP_MAG_RAW        PSF_INST_MAG_SIG  Mraw Mraw  2 -10.0 10.0
     590end
     591
     592# if we run this test as a stand-alone program somewhere, we may need to create a local copy of the photcode file:
     593macro mkphotcodes
     594  if ($0 != 2)
     595    echo "USAGE: mkphotcodes (filename)"
     596    break
     597  end
     598
     599  exec /bin/rm -f $1
     600  output $1
     601  echo "#                                           airmass      color                         astrometry  mag    photom  astrom mask    photom mask"
     602  echo "# code  name                type    zero  slope offset c1    c2   slope   zero  equiv  sys scale   scale  sys     poor   bad     poor   bad"
     603  echo "  1     g_SYNTH              sec   0.000  0.000 0.000     1     3 0.0000     0    21   0.000 0.000 0.000  0.000   0x0000 0x0000  0x0000 0x0000"
     604  echo "  2     r_SYNTH              sec   0.000  0.000 0.000     2     3 0.0000     0    22   0.000 0.000 0.000  0.000   0x0000 0x0000  0x0000 0x0000"
     605  echo "  3     i_SYNTH              sec   0.000  0.000 0.000     2     3 0.0000     0    23   0.000 0.000 0.000  0.000   0x0000 0x0000  0x0000 0x0000"
     606  echo "  4     z_SYNTH              sec   0.000  0.000 0.000     3     4 0.0000     0    24   0.000 0.000 0.000  0.000   0x0000 0x0000  0x0000 0x0000"
     607  echo "  5     y_SYNTH              sec   0.000  0.000 0.000     4     5 0.0000     0    25   0.000 0.000 0.000  0.000   0x0000 0x0000  0x0000 0x0000"
     608  echo "  3001  SYNTH.g              ref   0.000  0.000 0.000     -     - 0.0000     0     1   0.000 0.000 0.000  0.000   0x0000 0x0000  0x0000 0x0000"
     609  echo "  3002  SYNTH.r              ref   0.000  0.000 0.000     -     - 0.0000     0     2   0.000 0.000 0.000  0.000   0x0000 0x0000  0x0000 0x0000"
     610  echo "  3003  SYNTH.i              ref   0.000  0.000 0.000     -     - 0.0000     0     3   0.000 0.000 0.000  0.000   0x0000 0x0000  0x0000 0x0000"
     611  echo "  3004  SYNTH.z              ref   0.000  0.000 0.000     -     - 0.0000     0     4   0.000 0.000 0.000  0.000   0x0000 0x0000  0x0000 0x0000"
     612  echo "  3005  SYNTH.y              ref   0.000  0.000 0.000     -     - 0.0000     0     5   0.000 0.000 0.000  0.000   0x0000 0x0000  0x0000 0x0000"
     613  output stdout
     614end
     615
     616macro ckradialflux
     617 data $1
     618 read -fits Chip.xrad X_APER Y_APER PSF_FWHM APER_FLUX APER_FLUX_ERR APER_FLUX_STDEV APER_FILL
     619 read -fits Chip.psf X_PSF Y_PSF PSF_INST_MAG
     620 set mag = PSF_INST_MAG
     621 # XXX include nradii in header
     622 for i 0 11
     623   set mA$i = -2.5*log(APER_FLUX:$i)
     624   set dM$i = mA$i - PSF_INST_MAG
     625   vstat dM$i
     626   $DM$i = $MEDIAN
     627 end
     628 for i 1 11
     629   $j = $i - 1
     630   echo {$DM$i - $DM$j}
     631 end
     632
     633 lim mag -2 2; clear; box
     634 for i 0 11
     635   plot mag dM$i -c black
     636 end
     637end
     638
     639if ($SCRIPT)
     640  fulltest 4
    203641  exit 0
    204642end
    205 
    206 
    207   # write tmp.dat X_PSF Y_PSF X_PSF_SIG Y_PSF_SIG PSF_INST_MAG PSF_INST_MAG_SIG
    208   # exec gcompare -m tmp.dat 1 2 testimage.dat 1 2 2.0 > tmp.match.dat
    209   #
    210   # data tmp.match.dat
    211   # read x1 1 y1 2 dx1 3 dy1 4 psfmag1 5 dpsfmag1 6 x2 7 y2 8 flux 9 galaxy 10 mag2 11 dmag2 12 sxx 13 syy 14 sxy 15 s7 16
Note: See TracChangeset for help on using the changeset viewer.