Index: /branches/eam_branches/ipp-20111122/Ohana/src/addstar/src/mkcmf.c =================================================================== --- /branches/eam_branches/ipp-20111122/Ohana/src/addstar/src/mkcmf.c (revision 33387) +++ /branches/eam_branches/ipp-20111122/Ohana/src/addstar/src/mkcmf.c (revision 33388) @@ -14,5 +14,5 @@ void gauss_init (int Nbin); double rnd_gauss (double mean, double sigma); -void writeStars_PS1_V2 (FTable *ftable, double *X, double *Y, double *M, int *Flag, int Nstars); +void writeStars_PS1_V2 (FTable *ftable, double *X, double *Y, double *M, unsigned int *Flag, int Nstars); void writeStars_PS1_V1 (FTable *ftable, double *X, double *Y, double *M, int Nstars); void writeStars_PS1_DEV_1 (FTable *ftable, double *X, double *Y, double *M, int Nstars); @@ -27,5 +27,6 @@ // load a text table with X,Y,Mag (instrumental?) - int N, Nstars, NSTARS, found, *Flag; + int N, Nstars, NSTARS, found; + unsigned int *Flag; double *X, *Y, *M, Xmax, Ymax; @@ -181,5 +182,5 @@ ALLOCATE (Y, double, NSTARS); ALLOCATE (M, double, NSTARS); - ALLOCATE (Flag, int, NSTARS); + ALLOCATE (Flag, unsigned int, NSTARS); Xmax = Ymax = 0; @@ -187,5 +188,5 @@ int status; double ra, dec, mag, xobs, yobs, xraw, yraw, mraw; - int flags; + unsigned int flags; if (FROM_COORDS) { status = fscanf (f, "%lf %lf %lf %lf %lf %lf", &ra, &dec, &mag, &xraw, &yraw, &mraw); @@ -196,10 +197,10 @@ if (status == EOF) break; if (READ_FLAGS) { - status = fscanf (f, "%d", &flags); + status = fscanf (f, "%x", &flags); if (status == EOF) { fprintf (stderr, "error: missing flag for last star?\n"); exit (1); } - Flag[Nstars] = mraw; + Flag[Nstars] = flags; } else { Flag[Nstars] = 0; @@ -218,5 +219,5 @@ REALLOCATE (Y, double, NSTARS); REALLOCATE (M, double, NSTARS); - REALLOCATE (Flag, int, NSTARS); + REALLOCATE (Flag, unsigned int, NSTARS); } Nstars ++; @@ -274,5 +275,5 @@ } if (!strcmp(type, "PS1_V2")) { - writeStars_PS1_V2 (&ftable, X, Y, M, Nstars); + writeStars_PS1_V2 (&ftable, X, Y, M, Flag, Nstars); found = TRUE; } @@ -517,5 +518,5 @@ } -void writeStars_PS1_V2 (FTable *ftable, double *X, double *Y, double *M, int *Flag, int Nstars) { +void writeStars_PS1_V2 (FTable *ftable, double *X, double *Y, double *M, unsigned int *Flag, int Nstars) { int i; @@ -544,5 +545,5 @@ // randomly give poor PSFQF values - if ((BAD_PSFQF_FRAC > 0.0) && (drand48() < BAD_PSFQF_FRAC*Nstars)) { + if ((BAD_PSFQF_FRAC > 0.0) && (drand48() < BAD_PSFQF_FRAC)) { stars[i].psfQual = 0.25; } else { Index: /branches/eam_branches/ipp-20111122/Ohana/src/addstar/test/relphot.reject.dvo =================================================================== --- /branches/eam_branches/ipp-20111122/Ohana/src/addstar/test/relphot.reject.dvo (revision 33387) +++ /branches/eam_branches/ipp-20111122/Ohana/src/addstar/test/relphot.reject.dvo (revision 33388) @@ -33,4 +33,5 @@ $CHIP_DX = 1000 $CHIP_DY = 1000 + $DEFECT_FRAC = 0.03 if (not($?VERBOSE)) set VERBOSE = 0 @@ -52,19 +53,20 @@ # sequence to count images (only used in this function) - create tmpseq 0 9 + # XXX TEST create tmpseq 0 9 + create tmpseq 0 3 ##### mjd and zpt values for ubercal'ed images - $filt_uc:n = 9 + $filt_uc:n = tmpseq[] create filtN_uc 0 $filt_uc:n $filt_uc:0 = g; filtN_uc[0] = 1 $filt_uc:1 = g; filtN_uc[1] = 1 $filt_uc:2 = g; filtN_uc[2] = 1 - $filt_uc:3 = r; filtN_uc[3] = 2 - $filt_uc:4 = r; filtN_uc[4] = 2 - $filt_uc:5 = r; filtN_uc[5] = 2 - $filt_uc:6 = i; filtN_uc[6] = 3 - $filt_uc:7 = i; filtN_uc[7] = 3 - $filt_uc:8 = i; filtN_uc[8] = 3 + # $filt_uc:3 = r; filtN_uc[3] = 2 + # $filt_uc:4 = r; filtN_uc[4] = 2 + # $filt_uc:5 = r; filtN_uc[5] = 2 + # $filt_uc:6 = i; filtN_uc[6] = 3 + # $filt_uc:7 = i; filtN_uc[7] = 3 + # $filt_uc:8 = i; filtN_uc[8] = 3 # airmass slopes for this sequnece @@ -85,24 +87,24 @@ mjd_uc[1] = 55000.02 mjd_uc[2] = 55000.03 - mjd_uc[3] = 55015.01 - mjd_uc[4] = 55015.02 - mjd_uc[5] = 55015.03 - mjd_uc[6] = 55025.01 - mjd_uc[7] = 55025.02 - mjd_uc[8] = 55025.03 + # mjd_uc[3] = 55015.01 + # mjd_uc[4] = 55015.02 + # mjd_uc[5] = 55015.03 + # mjd_uc[6] = 55025.01 + # mjd_uc[7] = 55025.02 + # mjd_uc[8] = 55025.03 ##### mjd and zpt values for not-ubercal'ed images - $filt_nc:n = 9 + $filt_nc:n = tmpseq[] create filtN_nc 0 $filt_nc:n $filt_nc:0 = g; filtN_nc[0] = 1 $filt_nc:1 = g; filtN_nc[1] = 1 $filt_nc:2 = g; filtN_nc[2] = 1 - $filt_nc:3 = r; filtN_nc[3] = 2 - $filt_nc:4 = r; filtN_nc[4] = 2 - $filt_nc:5 = r; filtN_nc[5] = 2 - $filt_nc:6 = i; filtN_nc[6] = 3 - $filt_nc:7 = i; filtN_nc[7] = 3 - $filt_nc:8 = i; filtN_nc[8] = 3 + # $filt_nc:3 = r; filtN_nc[3] = 2 + # $filt_nc:4 = r; filtN_nc[4] = 2 + # $filt_nc:5 = r; filtN_nc[5] = 2 + # $filt_nc:6 = i; filtN_nc[6] = 3 + # $filt_nc:7 = i; filtN_nc[7] = 3 + # $filt_nc:8 = i; filtN_nc[8] = 3 # airmass slopes for these sequnece @@ -124,10 +126,10 @@ mjd_nc[1] = 55000.12 mjd_nc[2] = 55000.13 - mjd_nc[3] = 55015.11 - mjd_nc[4] = 55015.12 - mjd_nc[5] = 55015.13 - mjd_nc[6] = 55025.11 - mjd_nc[7] = 55025.12 - mjd_nc[8] = 55025.13 + # mjd_nc[3] = 55015.11 + # mjd_nc[4] = 55015.12 + # mjd_nc[5] = 55015.13 + # mjd_nc[6] = 55025.11 + # mjd_nc[7] = 55025.12 + # mjd_nc[8] = 55025.13 end @@ -163,9 +165,4 @@ end - # parallel-ize the database (goes to $catdir.p0) - mkparallel $catdir - - ### SERIAL database tests - # run setphot on the serial db and check that the images now match the expected values exec setphot -reset -update -ubercal testzpt.fits -D CATDIR $catdir @@ -178,4 +175,5 @@ # run relphot on the serial db and check that the images now match the expected values + if ($VERBOSE) echo relphot g,r,i -v -region 9.5 10.5 19.5 20.5 -D CATDIR $catdir -D STAR_TOOFEW 1 -D SIGMA_LIM 0.07 -statmode WT_MEAN -cloud-limit 0.5 -keep-ubercal -D IMAGE_OFFSET 0.5 -update >& log.relphot.s0 exec relphot g,r,i -v -region 9.5 10.5 19.5 20.5 -D CATDIR $catdir -D STAR_TOOFEW 1 -D SIGMA_LIM 0.07 -statmode WT_MEAN -cloud-limit 0.5 -keep-ubercal -D IMAGE_OFFSET 0.5 -update >& log.relphot.s0 for i 0 mjd_uc[] @@ -186,78 +184,5 @@ end - ### PARALLEL database tests - - # run setphot on the parallel db and check that the images now match the expected values - exec setphot -reset -parallel -update -ubercal testzpt.fits -D CATDIR $catdir.p0 - - # bring the detections back to the local db - # exec dvodist -in $catdir.p0 >& log.dvodist.in - # XXX this should now not be needed - - # check the results - for i 0 mjd_uc[] - ckexposure $catdir.p0 mjd_uc[$i] zpt_uc[$i] exptime_uc[$i] airmass_uc[$i] $filt_uc:$i setphot_uc - end - for i 0 mjd_nc[] - ckexposure $catdir.p0 mjd_nc[$i] zpt_nc[$i] exptime_nc[$i] airmass_nc[$i] $filt_nc:$i setphot_nc - end - - # run relphot on the parallel db and check that the images now match the expected values - exec relphot -parallel g,r,i -v -region 9.5 10.5 19.5 20.5 -D CATDIR $catdir.p0 -D STAR_TOOFEW 1 -D SIGMA_LIM 0.07 -statmode WT_MEAN -cloud-limit 0.5 -keep-ubercal -D IMAGE_OFFSET 0.5 -update >& log.relphot.p0 - - for i 0 mjd_uc[] - ckexposure $catdir.p0 mjd_uc[$i] zpt_uc[$i] exptime_uc[$i] airmass_uc[$i] $filt_uc:$i relphot - end - for i 0 mjd_nc[] - ckexposure $catdir.p0 mjd_nc[$i] zpt_nc[$i] exptime_nc[$i] airmass_nc[$i] $filt_nc:$i relphot - end - tapDONE -end - -macro qt - - ## SER - - # run setphot on the serial db and check that the images now match the expected values - exec setphot -reset -update -ubercal testzpt.fits -D CATDIR $catdir - for i 0 mjd_uc[] - ckexposure $catdir mjd_uc[$i] zpt_uc[$i] exptime_uc[$i] airmass_uc[$i] $filt_uc:$i setphot_uc - end - for i 0 mjd_nc[] - ckexposure $catdir mjd_nc[$i] zpt_nc[$i] exptime_nc[$i] airmass_nc[$i] $filt_nc:$i setphot_nc - end - - # run relphot on the serial db and check that the images now match the expected values - exec relphot g,r,i -v -region 9.5 10.5 19.5 20.5 -D CATDIR $catdir -D STAR_TOOFEW 1 -D SIGMA_LIM 0.07 -statmode WT_MEAN -cloud-limit 0.5 -keep-ubercal -D IMAGE_OFFSET 0.5 -update >& log.relphot.s0 - for i 0 mjd_uc[] - ckexposure $catdir mjd_uc[$i] zpt_uc[$i] exptime_uc[$i] airmass_uc[$i] $filt_uc:$i relphot - end - for i 0 mjd_nc[] - ckexposure $catdir mjd_nc[$i] zpt_nc[$i] exptime_nc[$i] airmass_nc[$i] $filt_nc:$i relphot - end - - ## PAR - - # run setphot on the parallel db and check that the images now match the expected values - exec setphot -reset -parallel -update -ubercal testzpt.fits -D CATDIR $catdir.p0 - - # check the results - for i 0 mjd_uc[] - ckexposure $catdir.p0 mjd_uc[$i] zpt_uc[$i] exptime_uc[$i] airmass_uc[$i] $filt_uc:$i setphot_uc - end - for i 0 mjd_nc[] - ckexposure $catdir.p0 mjd_nc[$i] zpt_nc[$i] exptime_nc[$i] airmass_nc[$i] $filt_nc:$i setphot_nc - end - - # run relphot on the parallel db and check that the images now match the expected values - exec relphot -parallel g,r,i -v -region 9.5 10.5 19.5 20.5 -D CATDIR $catdir.p0 -D STAR_TOOFEW 1 -D SIGMA_LIM 0.07 -statmode WT_MEAN -cloud-limit 0.5 -keep-ubercal -D IMAGE_OFFSET 0.5 -update >& log.relphot.p0 - - for i 0 mjd_uc[] - ckexposure $catdir.p0 mjd_uc[$i] zpt_uc[$i] exptime_uc[$i] airmass_uc[$i] $filt_uc:$i relphot - end - for i 0 mjd_nc[] - ckexposure $catdir.p0 mjd_nc[$i] zpt_nc[$i] exptime_nc[$i] airmass_nc[$i] $filt_nc:$i relphot - end end @@ -456,4 +381,6 @@ $options = $options -photcode GPC1.$FILTER.XY$ix\$iy $options = $options -no-noise + $options = $options -flags + $options = $options -bad-psfqf-frac $DEFECT_FRAC $options = $options -size $CHIP_DX $CHIP_DY $options = $options -crpix {0.5*$CHIP_DX} {0.5*$CHIP_DY} @@ -512,11 +439,44 @@ end - # XXX REJECT : I can deviate a single object here easily (save to a file) - # how to set flags? - # need to update mkcmf to read in flags values, set them here manually + # create bad detections / detections to ignore or avoid, save to a file? + + # create 3 types of things to ignore / avoid + # * badPhot (deviate the photometry) + # * badFlag (set a bad flag bit) + # * poorFlag (set a bad flag bit) + set badPhot = rnd(stars_r) < $DEFECT_FRAC + set badFlag = rnd(stars_r) < $DEFECT_FRAC + set poorFlag = rnd(stars_r) < $DEFECT_FRAC + + # create deviant photometry (1%) + # I'm setting the flag value 0x02 (EXTMODEL) if the object is deviated. This is a simple way to pass the info + set stars_mi = stars_mi + 0.2*badPhot + set stars_fl = zero(stars_mi) + 0x00000080*badFlag + 0x00000100*poorFlag + 0x00000002*badPhot + + # bad flags that we sould respect: + # PM_SOURCE_MODE_FAIL = 0x00000008, ///< Fit (non-linear) failed (non-converge, off-edge, run to zero) + # PM_SOURCE_MODE_SATSTAR = 0x00000080, ///< Source model peak is above saturation + # PM_SOURCE_MODE_BADPSF = 0x00000400, ///< Failed to get good estimate of object's PSF + # PM_SOURCE_MODE_DEFECT = 0x00000800, ///< Source is thought to be a defect + # PM_SOURCE_MODE_SATURATED = 0x00001000, ///< Source is thought to be saturated pixels (bleed trail) + # PM_SOURCE_MODE_CR_LIMIT = 0x00002000, ///< Source has crNsigma above limit + # PM_SOURCE_MODE_MOMENTS_FAILURE = 0x00008000, ///< could not measure the moments + # PM_SOURCE_MODE_SKY_FAILURE = 0x00010000, ///< could not measure the local sky + # PM_SOURCE_MODE_SKYVAR_FAILURE = 0x00020000, ///< could not measure the local sky variance + # PM_SOURCE_MODE_SIZE_SKIPPED = 0x10000000, ///< size could not be determined + + # poor flags that we should respect: + # PM_SOURCE_MODE_POOR = 0x00000010, ///< Fit succeeds, but low-SN, high-Chisq, or large (for PSF -- drop?) + # PM_SOURCE_MODE_PAIR = 0x00000020, ///< Source fitted with a double psf + # PM_SOURCE_MODE_BLEND = 0x00000100, ///< Source is a blend with other sources + # PM_SOURCE_MODE_BELOW_MOMENTS_SN = 0x00040000, ///< moments not measured due to low S/N + # PM_SOURCE_MODE_BLEND_FIT = 0x00400000, ///< source was fitted as a blend + # PM_SOURCE_MODE_ON_SPIKE = 0x20000000, ///< peak lands on diffraction spike + # PM_SOURCE_MODE_ON_GHOST = 0x40000000, ///< peak lands on ghost or glint + # PM_SOURCE_MODE_OFF_CHIP = 0x80000000, ///< peak lands off edge of chip # delete the output file before writing exec rm -f $1 - write -f "%10.6f %10.6f %7.3f %6.1f %6.1f %7.3f" $1 stars_r stars_d stars_mr stars_x stars_y stars_mi + write -f "%10.6f %10.6f %7.3f %6.1f %6.1f %7.3f 0x%x" $1 stars_r stars_d stars_mr stars_x stars_y stars_mi stars_fl end