IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 33169


Ignore:
Timestamp:
Jan 27, 2012, 10:06:11 AM (14 years ago)
Author:
eugene
Message:

setphot / relphot now seem to work well with the ubercal input zero table

Location:
branches/eam_branches/ipp-20111122/Ohana/src
Files:
8 edited

Legend:

Unmodified
Added
Removed
  • branches/eam_branches/ipp-20111122/Ohana/src/addstar/src/mkcmf.c

    r33138 r33169  
    8383    remove_argument (N, &argc, argv);
    8484  }
     85
     86  float exptime = 1.0;
     87  if ((N = get_argument (argc, argv, "-exptime"))) {
     88    remove_argument (N, &argc, argv);
     89    exptime = atof (argv[N]);
     90    remove_argument (N, &argc, argv);
     91  } 
     92
     93  // XXX note that the airmass and ra,dec,mjd can be inconsistent (for a given observatory location)
     94  float airmass = 1.0;
     95  if ((N = get_argument (argc, argv, "-airmass"))) {
     96    remove_argument (N, &argc, argv);
     97    airmass = atof (argv[N]);
     98    remove_argument (N, &argc, argv);
     99  } 
    85100
    86101  int NX = 0;
     
    202217  }
    203218  gfits_modify (&header, "ZERO_PT", "%lf", 1, 25.0);
    204   gfits_modify (&header, "EXPTIME", "%lf", 1, 1.0);
    205   gfits_modify (&header, "AIRMASS", "%lf", 1, 1.0);
     219  gfits_modify (&header, "EXPTIME", "%lf", 1, exptime);
     220  gfits_modify (&header, "AIRMASS", "%lf", 1, airmass);
    206221  gfits_modify (&header, "NASTRO",   "%d", 1, 10);
    207222
  • branches/eam_branches/ipp-20111122/Ohana/src/addstar/test/relphot.flatcorr.dvo

    r33159 r33169  
    3434  # we have two sets of images: ubercaled and not-ubercaled
    3535
     36  # sequence to count images (only used in this function)
     37  create tmpseq 0 9
     38
    3639  # mjd and zpt values for ubercal'ed imaged
    37   create tmpseq 0 9
    38 
    39   set zpt_uc = 25.0 + tmpseq*0.0050 - 0.0025
     40  set exptime_uc = 10.0 + zero(tmpseq)
     41  set airmass_uc = 1.3 + zero(tmpseq)
     42
     43  # ubercal zero points are defined as ZP_nominal + 2.5log(exptime) + K*(airmass - 1.0)
     44  set zpt_uc = 25.0 + 2.5*log(exptime_uc) - 0.15*(airmass_uc - 1.0) + tmpseq*0.0050 - 0.0025
     45
    4046  set mjd_uc = zero(zpt_uc)
    4147  mjd_uc[0] = 55000.01
    4248  mjd_uc[1] = 55000.02
    4349  mjd_uc[2] = 55000.03
    44 
    4550  mjd_uc[3] = 55015.01
    4651  mjd_uc[4] = 55015.02
    4752  mjd_uc[5] = 55015.03
    48 
    4953  mjd_uc[6] = 55025.01
    5054  mjd_uc[7] = 55025.02
     
    5357  # mjd and zpt values for not-ubercal'ed images
    5458  # place this within a valid season (55000.0 - 55010.0 - 55020.0 - 55030.0)
    55   set zpt_nc = 25.0 - tmpseq*0.0050 + 0.0025
     59  set exptime_nc = 15.0 + zero(tmpseq)
     60  set airmass_nc = 1.6 + zero(tmpseq)
     61
     62  # ubercal zero points are defined as ZP_nominal + 2.5log(exptime) + K*(airmass - 1.0)
     63  set zpt_nc = 25.0 + 2.5*log(exptime_nc) - 0.15*(airmass_nc - 1.0) - tmpseq*0.0050 + 0.0025
     64
    5665  set mjd_nc = zero(zpt_nc)
    5766  mjd_nc[0] = 55000.11
     
    8796
    8897  mksequence $fileroot $catdir
    89 end
    90 
    91 # go
    92 # ckexposure catdir.test mjd_uc[0] zpt_uc[0] raw
    93 # ckexposure catdir.test mjd_nc[0] zpt_nc[0] raw
     98  ckexposure catdir.test mjd_uc[0] zpt_uc[0] exptime_uc[0] airmass_uc[0] raw
     99  ckexposure catdir.test mjd_nc[0] zpt_nc[0] exptime_nc[0] airmass_nc[0] raw
     100
     101  exec setphot -update -ubercal testzpt.fits -D CATDIR $catdir
     102  ckexposure catdir.test mjd_uc[0] zpt_uc[0] exptime_uc[0] airmass_uc[0] corr
     103  ckexposure catdir.test mjd_nc[0] zpt_nc[0] exptime_nc[0] airmass_nc[0] corr
     104
     105  exec relphot g -v -region 9.5 10.5 19.5 20.5 -D CATDIR catdir.test -D STAR_TOOFEW 1 -D SIGMA_LIM 0.07 -statmode WT_MEAN -cloud-limit 0.5 -keep-ubercal -D IMAGE_OFFSET 0.5 -update
     106  ckexposure catdir.test mjd_uc[0] zpt_uc[0] exptime_uc[0] airmass_uc[0] corr
     107  ckexposure catdir.test mjd_nc[0] zpt_nc[0] exptime_nc[0] airmass_nc[0] corr
     108end
     109
     110# go testdata/t1 catdir.test
     111# ckexposure catdir.test mjd_uc[0] zpt_uc[0] raw : result is mean of 0.0 (since test corrects for ZPT_NOMINAL vs REAL)
     112# ckexposure catdir.test mjd_nc[0] zpt_nc[0] raw : result is mean of 0.0
    94113# setphot -update -ubercal testzpt.fits -D CATDIR catdir.test
    95 # ckexposure catdir.test mjd_uc[0] zpt_uc[0] corr
    96 # ckexposure catdir.test mjd_nc[0] zpt_nc[0] corr
     114# ckexposure catdir.test mjd_uc[0] zpt_uc[0] corr : result is mean of 0.0 (since image is corrected to ZPT_REAL)
     115# ckexposure catdir.test mjd_nc[0] zpt_nc[0] corr : result is mean of -0.42 (since image is NOT corrected to ZPT_REAL)
    97116# relphot g -region 9.5 10.5 19.5 20.5 -update -D CATDIR catdir.test -D STAR_TOOFEW 1 -D SIGMA_LIM 0.075 -statmode MEAN
    98117# dvo: ckexposure catdir.test mjd_uc[0] zpt_uc[0] corr
     
    104123# re-run setphot, then relphot with -keep-ubercal
    105124
    106 
     125# the test case has only 2 exposures, so the STAR_TOOFEW of 1 is necessary.  Also, the mag and sigma distributions are somewhat artificial
     126# relphot g -v -region 9.5 10.5 19.5 20.5 -D CATDIR catdir.test -D STAR_TOOFEW 1 -D SIGMA_LIM 0.07 -statmode WT_MEAN -cloud-limit 0.5 -keep-ubercal -D IMAGE_OFFSET 0.5 -update
     127# dvo: ckexposure catdir.test mjd_uc[0] zpt_uc[0] corr
     128#  mean is now 0.0
     129# dvo: ckexposure catdir.test mjd_nc[0] zpt_nc[0] corr
     130#  mean is now 0.0
     131#  ** setphot correctly assigned the zero point of UC & relphot adjusted NC to match
    107132
    108133macro ckexposure
    109   if ($0 != 5)
    110     echo "ckexposure (catdir) (mjd) (zpt) (mode)"
     134  if ($0 != 7)
     135    echo "ckexposure (catdir) (mjd) (zpt) (exptime) (airmass) (mode)"
    111136    echo "  mode == raw or corr"
    112137    break
    113138  end
    114139
    115   local CATDIR MJD_IMAGE ZPT_REAL ZPT_NOMINAL
     140  local CATDIR MJD_IMAGE ZPT_REAL ZPT_NOMINAL EXPTIME AIRMASS
    116141 
    117142  $CATDIR      = $1
     
    119144  $ZPT_REAL    = $3
    120145  $ZPT_NOMINAL = 24.58
    121   $MODE        = $4
     146  $EXPTIME     = $4
     147  $AIRMASS     = $5
     148  $MODE        = $6
    122149  # XXX need a function to extract the nominal zpt for a given filter / photcode
    123150
     
    153180
    154181  # uncorrected values behave like this:
    155   #   mag_DVO  = m_inst + zpt_nominal
     182  #   mag_DVO  = m_inst + zpt_nominal + 2.5*log(exptime) + K*(airmass - 1.0)
    156183  #   mag_real = m_inst + zpt_real + cell_offset
    157184  #   dm = mag_DVO - mag_real = zpt_nominal - zpt_real - cell_offset
    158185  #   <dm> - zpt_nominal + zpt_real + cell_offset ~ 0.0
    159   #   zpt_real (in this case) = 25.0 (actually, it is the value in the vector 'zpt' for this entry
     186  #   zpt_real (in this case) = 25.0  + 2.5*log(exptime) + K*(airmass - 1.0)
     187  #     (actually, it is the value in the vector 'zpt' for this entry
    160188  #   zpt_nominal = 24.58
    161189
     
    170198        # echo cell_off[$ix][$iy] {$MEDIAN - $ZPT_NOMINAL + $ZPT_REAL + cell_off[$ix][$iy]} $MEDIAN $MEAN $SIGMA
    171199 
    172         set dm_adjust = dm$ix\$iy - $ZPT_NOMINAL + $ZPT_REAL + cell_off[$ix][$iy]
     200        $ZPT_REAL_NORM = $ZPT_REAL - 2.5*log($EXPTIME) + 0.15*($AIRMASS - 1.0)
     201        set dm_adjust = dm$ix\$iy - $ZPT_NOMINAL + $ZPT_REAL_NORM + cell_off[$ix][$iy]
    173202        vstat dm_adjust
    174203      end
     
    190219  # for i 0 mjd_uc[]
    191220  for i 0 1
    192     mkexposure $1.uc.$i $RA_CENTER $DEC_CENTER zpt_uc[$i] mjd_uc[$i] g $2
     221    mkexposure $1.uc.$i $RA_CENTER $DEC_CENTER zpt_uc[$i] exptime_uc[$i] airmass_uc[$i] mjd_uc[$i] g $2
    193222  end
    194223
     
    196225  # for i 0 mjd_nc[]
    197226  for i 0 1
    198     mkexposure $1.nc.$i $RA_CENTER $DEC_CENTER zpt_nc[$i] mjd_nc[$i] g $2
     227    mkexposure $1.nc.$i $RA_CENTER $DEC_CENTER zpt_nc[$i] exptime_nc[$i] airmass_nc[$i] mjd_nc[$i] g $2
    199228  end
    200229end
    201230
    202231macro mkexposure
    203  if ($0 != 8)
    204    echo "mkexposure (fileroot) (ra) (dec) (zpt) (mjd) (filter) (catdir)"
     232 if ($0 != 10)
     233   echo "mkexposure (fileroot) (ra) (dec) (zpt) (exptime) (airmass) (mjd) (filter) (catdir)"
    205234   break
    206235 end
    207236
    208  local ix iy date time datetime ra dec ROOT RAo DECo ZPT MJD FILTER
    209 
    210  $ROOT   = $1
    211  $RAo    = $2
    212  $DECo   = $3
    213  $ZPT    = $4
    214  $MJD    = $5
    215  $FILTER = $6
     237 local ix iy date time datetime ra dec ROOT RAo DECo ZPT MJD FILTER EXPTIME AIRMASS CATDIR
     238
     239 $ROOT    = $1
     240 $RAo     = $2
     241 $DECo    = $3
     242 $ZPT     = $4
     243 $EXPTIME = $5
     244 $AIRMASS = $6
     245 $MJD     = $7
     246 $FILTER  = $8
     247 $CATDIR  = $9
    216248
    217249 $TIMEFORMAT = mjd
     
    229261    $dec = $DECo - $dy / 3600.0
    230262    # echo $ra $dec $dx $dy
    231     mkinput test.in.txt $ra $dec $ZPT
     263    mkinput test.in.txt $ra $dec $ZPT 
    232264   
    233265    # ra,dec is the center of this chip
    234     exec mkcmf test.in.txt $ROOT.$ix.$iy.cmf -mjd $MJD -radec $ra $dec -type PS1_V2 -coords -photcode GPC1.$FILTER.XY$ix\$iy -no-noise -size $CHIP_DX $CHIP_DY -crpix {0.5*$CHIP_DX} {0.5*$CHIP_DY}
    235     # echo mkcmf test.in.txt $ROOT.$ix.$iy.cmf -mjd $MJD -radec $ra $dec -type PS1_V2 -coords -photcode GPC1.$FILTER.XY$ix\$iy -no-noise -size $CHIP_DX $CHIP_DY -crpix {0.5*$CHIP_DX} {0.5*$CHIP_DY}
     266    local options
     267    $options = -mjd $MJD
     268    $options = $options -radec $ra $dec
     269    $options = $options -type PS1_V2
     270    $options = $options -coords
     271    $options = $options -photcode GPC1.$FILTER.XY$ix\$iy
     272    $options = $options -no-noise
     273    $options = $options -size $CHIP_DX $CHIP_DY
     274    $options = $options -crpix {0.5*$CHIP_DX} {0.5*$CHIP_DY}
     275    $options = $options -airmass $AIRMASS
     276    $options = $options -exptime $EXPTIME
     277    exec mkcmf test.in.txt $ROOT.$ix.$iy.cmf $options
     278    echo mkcmf test.in.txt $ROOT.$ix.$iy.cmf $options
    236279   
    237280    # the fake images have inconsistent ra,dec and airmass,sidtime values
    238     exec addstar $ROOT.$ix.$iy.cmf -D CATDIR $7 -D CAMERA gpc1 -quick-airmass
     281    echo addstar $ROOT.$ix.$iy.cmf -D CATDIR $CATDIR -D CAMERA gpc1 -quick-airmass
     282    exec addstar $ROOT.$ix.$iy.cmf -D CATDIR $CATDIR -D CAMERA gpc1 -quick-airmass
     283    echo done
    239284  end
    240285 end
     
    262307  subset stars_mr = stars_mag if valid
    263308
    264   set cell_x = stars_x / ($CHIP_DX / $NCELL_X)
    265   set cell_y = stars_y / ($CHIP_DX / $NCELL_X)
     309  set cell_x = int(stars_x / ($CHIP_DX / $NCELL_X))
     310  set cell_y = int(stars_y / ($CHIP_DX / $NCELL_X))
    266311
    267312  # eddie says: M_real = M_inst + zpt + ff_offset
  • branches/eam_branches/ipp-20111122/Ohana/src/relphot/include/relphot.h

    r33157 r33169  
    5454double MIN_ERROR;
    5555double IMFIT_SYS_SIGMA_LIM;
     56double CLOUD_TOLERANCE;
    5657
    5758int    VERBOSE;
  • branches/eam_branches/ipp-20111122/Ohana/src/relphot/src/ImageOps.c

    r33157 r33169  
    312312
    313313  if (FREEZE_IMAGES) return;
     314
     315  fprintf (stderr, "limiting negative clouds to %f\n", CLOUD_TOLERANCE);
    314316
    315317  int Nsecfilt = GetPhotcodeNsecfilt ();
     
    383385      }
    384386
    385       fprintf (stderr, "%1d, %3d : %3d, %3d : %10.6f %10.6f : %6.3f  %6.3f  %6.3f  %6.3f  %6.3f\n", (int) i, (int) j, (int) c, (int) m, catalog[c].averageT[n].R, catalog[c].averageT[n].D, Msys, Mrel, Mmos, Mgrid, Mflat);
     387      // fprintf (stderr, "%1d, %3d : %3d, %3d : %10.6f %10.6f : %6.3f  %6.3f  %6.3f  %6.3f  %6.3f\n", (int) i, (int) j, (int) c, (int) m, catalog[c].averageT[n].R, catalog[c].averageT[n].D, Msys, Mrel, Mmos, Mgrid, Mflat);
    386388
    387389      list[N] = Msys - Mrel - Mmos - Mgrid + Mflat;
     
    409411    liststats (list, dlist, N, &stats);
    410412
    411     float CLOUD_TOLERANCE = 0.01;
    412413    image[i].Mcal  = stats.mean;
    413414    image[i].dMcal = stats.error;
     
    440441  int mark, Nmark;
    441442  off_t i, N;
    442   double *ImageToMeasure, *slist, *dlist;
     443  double *mlist, *slist, *dlist;
    443444  double MaxOffset, MaxScatter, MedOffset;
    444445  StatType stats;
     
    448449  if (VERBOSE) fprintf (stderr, "marking poor images\n");
    449450
    450   ALLOCATE (ImageToMeasure, double, Nimage);
     451  ALLOCATE (mlist, double, Nimage);
    451452  ALLOCATE (slist, double, Nimage);
    452453  ALLOCATE (dlist, double, Nimage);
    453454
     455  // measure stats for Mcal and dMcal
    454456  for (i = N = 0; i < Nimage; i++) {
    455457    if (image[i].flags & IMAGE_BAD) continue;
    456     ImageToMeasure[N] = fabs (image[i].Mcal);
     458    mlist[N] = image[i].Mcal;
    457459    slist[N] = image[i].dMcal;
    458460    dlist[N] = 1;
     
    460462  }
    461463  initstats ("MEAN");
    462   liststats (ImageToMeasure, dlist, N, &stats);
     464  liststats (mlist, dlist, N, &stats);
    463465  MaxOffset = MAX (IMAGE_OFFSET, 3*stats.sigma);
    464466  MedOffset = stats.median;
     
    469471  Nmark = 0;
    470472  for (i = 0; i < Nimage; i++) {
     473    // if we are keeping ubercal sacrosanct, then we should not be allowed to break them...
     474    if (KEEP_UBERCAL && (image[i].flags & ID_IMAGE_PHOTOM_UBERCAL)) continue;
     475
    471476    mark = FALSE;
    472477    image[i].flags &= ~ID_IMAGE_PHOTOM_POOR;
     
    482487  fprintf (stderr, "%d images marked poor\n", Nmark);
    483488  initstats (STATMODE);
    484   free (ImageToMeasure);
     489  free (mlist);
    485490  free (slist);
    486491  free (dlist);
  • branches/eam_branches/ipp-20111122/Ohana/src/relphot/src/StarOps.c

    r33157 r33169  
    8484            continue;
    8585          }
    86 
    87           // XXX include the ubercal measurements?  yes, of course!  these are the well-determined values, supposedly!
    8886
    8987          // XXX allow REF stars (no Image Entry) to be included in the calculation this
     
    118116          dlist[N] = MAX (catalog[i].measureT[m].dM, MIN_ERROR);
    119117
     118          // up-weight the ubercal values (or convergence can take a long time...)
     119          if (catalog[i].measureT[m].dbFlags & ID_MEAS_PHOTOM_UBERCAL) {
     120            dlist[N] = MAX (0.1*catalog[i].measureT[m].dM, MIN_ERROR);
     121          }
     122
    120123          // tie down reference photometry if the -refcode (code) option is selected
     124          // XXX re-think this...
    121125          if (refPhotcode) {
    122126            if (GetPhotcodeEquivCodebyCode(catalog[i].measureT[m].photcode) == refPhotcode[0].equiv) {
     
    140144        catalog[i].secfilt[Nsecfilt*j+Nsec].M  = stats.mean;
    141145        catalog[i].secfilt[Nsecfilt*j+Nsec].dM = stats.sigma;
    142         catalog[i].secfilt[Nsecfilt*j+Nsec].Xm = (stats.Nmeas > 1) ? 100.0*log10(stats.chisq) : NAN_S_SHORT;
     146        catalog[i].secfilt[Nsecfilt*j+Nsec].Xm = (stats.Nmeas > 1) ? 100.0*log10(stats.chisq + 1e-4) : NAN_S_SHORT;
    143147      }
    144148    }
  • branches/eam_branches/ipp-20111122/Ohana/src/relphot/src/args.c

    r33157 r33169  
    142142  }
    143143
     144  CLOUD_TOLERANCE = 0.02;
     145  if ((N = get_argument (argc, argv, "-cloud-limit"))) {
     146    remove_argument (N, &argc, argv);
     147    CLOUD_TOLERANCE = atof(argv[N]);
     148    remove_argument (N, &argc, argv);
     149  }
     150
    144151  SHOW_PARAMS = FALSE;
    145152  if ((N = get_argument (argc, argv, "-params"))) {
  • branches/eam_branches/ipp-20111122/Ohana/src/uniphot/src/match_zpts_to_images.c

    r33133 r33169  
    6969
    7070    // we have a match: set zpt and record the match
    71     image[Ni].Mcal = SCALE*code[0].C - zpts[Nz].zpt;
     71    // is the zero point supplied nominally corrected or not?
     72    // UBERCAL includes 2.5log(exptime) + K*airmass in the zero point
     73    if (UBERCAL) {
     74      image[Ni].Mcal = SCALE*code[0].C - zpts[Nz].zpt + 2.5*log10(image[Ni].exptime) + code[0].K*(image[Ni].secz - 1.000);
     75    } else {
     76      image[Ni].Mcal = SCALE*code[0].C - zpts[Nz].zpt;
     77    }
    7278    image[Ni].dMcal = zpts[Nz].zpt_err;
    7379    image[Ni].flags &= ~ID_IMAGE_PHOTOM_NOCAL; // clear the NOCAL flag
  • branches/eam_branches/ipp-20111122/Ohana/src/uniphot/src/update_catalog_setphot.c

    r33156 r33169  
    9797      catalog[0].measure[m].dMcal = dMcal;
    9898
    99       if (UBERCAL) {
     99      // if we are setting the zero points from an UBERCAL database, and this detection is from one of those images,
     100      // then tag the measurement as well
     101      if (UBERCAL && (image[id].flags & ID_IMAGE_PHOTOM_UBERCAL)) {
    100102        catalog[0].measure[m].dbFlags |=  ID_MEAS_PHOTOM_UBERCAL;
    101103      }
Note: See TracChangeset for help on using the changeset viewer.