IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 40213


Ignore:
Timestamp:
Nov 22, 2017, 4:44:48 PM (9 years ago)
Author:
eugene
Message:

add selection for PSF vs APER zero points

Location:
branches/eam_branches/ohana.20170822/src/relphot
Files:
8 edited

Legend:

Unmodified
Added
Removed

  • branches/eam_branches/ohana.20170822/src/relphot/include/relphot.h

    r40212 r40213  
    425425off_t         getImageEntry       PROTO((off_t meas, int cat));
    426426
    427 float         getMcal_alt         PROTO((off_t meas, int cat, FlatCorrectionTable *flatcorr, float Xccd, float Yccd));
    428 float         getMcal             PROTO((off_t meas, int cat, FlatCorrectionTable *flatcorr, Catalog *catalog));
     427float         getMcal             PROTO((off_t meas, int cat, dvoMagClassType class));
    429428float         getMflat            PROTO((off_t meas, int cat, FlatCorrectionTable *flatcorr, Catalog *catalog));
    430429float         getMgrid            PROTO((off_t meas, int cat));

  • branches/eam_branches/ohana.20170822/src/relphot/src/GridOps.c

    r40212 r40213  
    401401
    402402        // skip images marked as BAD
    403         Mcal = getMcal  (m, c, flatcorr, catalog);
     403        Mcal = getMcal  (m, c, MAG_CLASS_PSF);
    404404        if (isnan(Mcal)) {
    405405          Ncal ++;
    … …  
    558558        continue;
    559559      }
    560       Mcal = getMcal  (m, c, flatcorr, catalog);
     560      Mcal = getMcal  (m, c, MAG_CLASS_PSF);
    561561      if (isnan(Mcal)) {
    562562        Ncal ++;
    … …  
    645645        continue;
    646646      }
    647       Mcal  = getMcal  (m, c, flatcorr, catalog);
     647      Mcal  = getMcal  (m, c, MAG_CLASS_PSF);
    648648      if (isnan(Mcal)) continue;
    649649      Mmos  = getMmos  (m, c);

  • branches/eam_branches/ohana.20170822/src/relphot/src/ImageOps.c

    r40212 r40213  
    432432}
    433433
    434 // returns image.Mcal
     434// returns image.McalPSF or image.McalAPER
    435435// NOTE: static flat-field component is included in measure.Mflat
    436 float getMcal (off_t meas, int cat, FlatCorrectionTable *flatcorr, Catalog *catalog) {
    437   OHANA_UNUSED_PARAM(flatcorr);
    438   OHANA_UNUSED_PARAM(catalog);
    439 
    440   off_t i;
    441   float value;
     436float getMcal (off_t meas, int cat, dvoMagClassType class) {
     437
     438  off_t i;
    442439
    443440  i = MeasureToImage[cat][meas];
    444   if (i == -1) return (NAN);
    445 
    446   if (image[i].flags & IMAGE_BAD) return (NAN); 
    447   value = image[i].McalPSF;
    448 
    449   // to do this, I need to pass in the catalog and flatcorr pointers
    450   // int flat_id = image[i].photom_map_id;
    451   // float offset = 0.0;
    452   // if (flat_id) {
    453   //   offset = FlatCorrectionOffset (flatcorr, flat_id, catalog[cat].measureT[meas].Xccd, catalog[cat].measureT[meas].Yccd);
    454   // }
    455   // value -= offset;
    456 
    457   return (value);
    458 }
    459 
    460 // returns image.Mcal - ff(x,y)
    461 // NOTE: static flat-field component is included in measure.Mflat
    462 float getMcal_alt (off_t meas, int cat, FlatCorrectionTable *flatcorr, float Xccd, float Yccd) {
    463   OHANA_UNUSED_PARAM(flatcorr);
    464   OHANA_UNUSED_PARAM(Xccd);
    465   OHANA_UNUSED_PARAM(Yccd);
    466 
    467   off_t i;
    468   float value;
    469 
    470   i = MeasureToImage[cat][meas];
    471   if (i == -1) return (NAN);
    472 
    473   if (image[i].flags & IMAGE_BAD) return (NAN); 
    474   value = image[i].McalPSF;
    475 
    476   // to do this, I need to pass in the catalog and flatcorr pointers
    477   // int flat_id = image[i].photom_map_id;
    478   // float offset = 0.0;
    479   // if (flat_id) {
    480   //   offset = FlatCorrectionOffset (flatcorr, flat_id, Xccd, Yccd);
    481   // }
    482   // value -= offset;
    483 
    484   return (value);
     441  if (i == -1) return NAN;
     442
     443  if (image[i].flags & IMAGE_BAD) return NAN; 
     444
     445  switch (class) {
     446    case MAG_CLASS_PSF:
     447      return image[i].McalPSF;
     448    case MAG_CLASS_APER:
     449    case MAG_CLASS_KRON:
     450      return image[i].McalAPER;
     451    default:
     452      return NAN;
     453  }
     454  return NAN; // should not be able to reach here
    485455}
    486456
    … …  
    615585void setMcal (Catalog *catalog, int PoorImages, FlatCorrectionTable *flatcorr) {
    616586
    617   off_t i, j, m, c, n, Nmax;
     587  off_t i, j, m, c, n;
    618588  int mark, bad, Nfew, Nbad, Nmos, Nrel, Ngrid, Nsys;
    619   float Msys, Mrel, Mmos, Mgrid, Mflat;
    620   // double *list, *dlist, *Mlist, *dMlist;
    621589
    622590  StatType stats;
    … …  
    635603  }
    636604
    637   Nmax = 0;
     605  off_t Nmax = 0;
    638606  for (i = 0; i < Nimage; i++) {
    639607    Nmax = MAX (Nmax, N_onImage[i]);
    640608  }
    641   StatDataSet *refStars = StatDataSetAlloc (1, Nmax);
     609
     610  StatDataSet *kronStars   = StatDataSetAlloc (1, Nmax);
     611  StatDataSet *psfStars    = StatDataSetAlloc (1, Nmax);
    642612  StatDataSet *brightStars = StatDataSetAlloc (1, Nmax);
    643613
    … …  
    690660          continue;
    691661      }
    692       Mmos  = getMmos  (m, c);
     662      float Mmos  = getMmos  (m, c);
    693663      if (isnan(Mmos)) {
    694664          Nmos ++;
    695665          continue;
    696666      }
    697       Mgrid = getMgrid (m, c);
     667      float Mgrid = getMgrid (m, c);
    698668      if (isnan(Mgrid)) {
    699669          Ngrid++;
    … …  
    701671      }
    702672
    703       // Mrel is the average magnitude for this star.  For PS1 stacks, we have too much
     673      // Mrel* is the average magnitude for this star.  For PS1 stacks, we have too much
    704674      // PSF variability.  We need to calibrate the PSF magnitudes separately from the
    705675      // Aperture-like magnitues.  (We have an option to use the kron magnitudes or the
    … …  
    707677      // magnitude type
    708678
    709       Mrel  = getMrel  (catalog, m, c, MAG_CLASS_PSF, MAG_SRC_CHP);
    710       if (isnan(Mrel)) {
     679      float MrelPSF  = getMrel  (catalog, m, c, MAG_CLASS_PSF, MAG_SRC_CHP);
     680      if (isnan(MrelPSF)) {
    711681          Nrel ++;
    712682          continue;
    713683      }
     684      float MrelKron  = getMrel  (catalog, m, c, MAG_CLASS_KRON, MAG_SRC_CHP);
    714685     
    715686      // image.Mcal is not supposed to include the flat-field correction, so we need to
    … …  
    718689      // the flat-correction.  NOTE the sign of Mflat (Image.Mcal = Measure.Mcal - Mflat)
    719690
    720       Mflat = getMflat (m, c, flatcorr, catalog);
     691      float Mflat = getMflat (m, c, flatcorr, catalog);
    721692
    722693      n = catalog[c].measureT[m].averef;
    723       Msys = PhotSysTiny (&catalog[c].measureT[m], &catalog[c].averageT[n], &catalog[c].secfilt[n*Nsecfilt], MAG_CLASS_PSF);
    724       if (isnan(Msys)) {
     694      float MsysPSF = PhotSysTiny (&catalog[c].measureT[m], &catalog[c].averageT[n], &catalog[c].secfilt[n*Nsecfilt], MAG_CLASS_PSF);
     695      if (isnan(MsysPSF)) {
    725696        Nsys++;
    726697        continue;
    727698      }
     699      float MsysKron = PhotSysTiny (&catalog[c].measureT[m], &catalog[c].averageT[n], &catalog[c].secfilt[n*Nsecfilt], MAG_CLASS_KRON);
    728700
    729701      PhotCode *code = GetPhotcodebyCode (catalog[c].measureT[m].photcode);
    … …  
    735707
    736708    skip:
    737       refStars->flxlist[Nref] = Msys - Mrel - Mmos - Mgrid + Mflat;
    738       refStars->errlist[Nref] = MAX (catalog[c].measureT[m].dM, MIN_ERROR);
    739       refStars->wgtlist[Nref] = 1;
    740       refStars->msklist[Nref] = 0;
    741       if (fabs(refStars->flxlist[Nref]) > 0.03) {
     709      psfStars->flxlist[Nref] = MsysPSF - MrelPSF - Mmos - Mgrid + Mflat;
     710      psfStars->errlist[Nref] = MAX (catalog[c].measureT[m].dM, MIN_ERROR);
     711      psfStars->wgtlist[Nref] = 1;
     712      psfStars->msklist[Nref] = 0;
     713      if (fabs(psfStars->flxlist[Nref]) > 0.03) {
    742714        // fprintf (stderr, "deviant\n");
    743715      }
    744716
     717      kronStars->flxlist[Nref] = MsysKron - MrelKron - Mmos - Mgrid + Mflat;
     718      kronStars->errlist[Nref] = MAX (catalog[c].measureT[m].dM, MIN_ERROR);
     719      kronStars->wgtlist[Nref] = 1;
     720      kronStars->msklist[Nref] = 0;
     721
    745722      if ((image[i].imageID == TEST_IMAGE1) || (image[i].imageID == TEST_IMAGE2)) {
    746         fprintf (stderr, "%1d, %3d : %3d, %3d : %10.6f %10.6f : %6.3f  %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, refStars->flxlist[Nref]);
     723        fprintf (stderr, "%1d, %3d : %3d, %3d : %10.6f %10.6f : %6.3f  %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, MsysPSF, MrelPSF, Mmos, Mgrid, Mflat, psfStars->flxlist[Nref]);
    747724      }
    748725
    749726      if (catalog[c].measureT[m].dM < IMFIT_SYS_SIGMA_LIM) {
    750         brightStars->flxlist[Nbright] = refStars->flxlist[Nref];
    751         brightStars->errlist[Nbright] = refStars->errlist[Nref];
     727        brightStars->flxlist[Nbright] = psfStars->flxlist[Nref];
     728        brightStars->errlist[Nbright] = psfStars->errlist[Nref];
    752729        brightStars->wgtlist[Nbright] = 1;
    753730        brightStars->msklist[Nbright] = 0;
    … …  
    774751    // no additional weight modification (we treat all stars on an image equally -- note an image is either ubercal-tied or not)
    775752# if (BASIC_STATS)
    776     liststats (refStars->flxlist, refStars->errlist, NULL, Nref, &stats);
     753    liststats (psfStars->flxlist, psfStars->errlist, NULL, Nref, &stats);
    777754# else
    778     liststats_irls (refStars, Nref, &stats);
     755    liststats_irls (psfStars, Nref, &stats);
    779756# endif
    780757    image[i].McalPSF    = stats.mean;
    781     image[i].McalAPER   = stats.mean;
    782758    image[i].dMcal      = stats.error;
    783759    image[i].nFitPhotom = Nref;
    … …  
    785761    Ncalibrated ++;
    786762
     763    // no additional weight modification (we treat all stars on an image equally -- note an image is either ubercal-tied or not)
     764# if (BASIC_STATS)
     765    liststats (kronStars->flxlist, kronStars->errlist, NULL, Nref, &stats);
     766# else
     767    liststats_irls (kronStars, Nref, &stats);
     768# endif
     769    image[i].McalAPER   = stats.mean;
     770
    787771    if ((image[i].imageID == TEST_IMAGE1) || (image[i].imageID == TEST_IMAGE2)) {
    788772      for (j = 0; j < Nref; j++) {
    789         fprintf (stderr, "%1d, %8d : %6.3f  %6.3f  %6.3f  %d\n", (int) i, (int) image[i].imageID, refStars->flxlist[j], refStars->errlist[j], refStars->wgtlist[j], refStars->msklist[j]);
     773        fprintf (stderr, "%1d, %8d : %6.3f  %6.3f  %6.3f  %d\n", (int) i, (int) image[i].imageID, psfStars->flxlist[j], psfStars->errlist[j], psfStars->wgtlist[j], psfStars->msklist[j]);
    790774      }
    791775    }
    … …  
    801785    if (PLOTSTUFF) {
    802786      fprintf (stderr, "Mcal for : %s : %7.4f %7.4f\n", image[i].name, image[i].McalPSF, image[i].dMcal);
    803       plot_setMcal (refStars->flxlist, Nref, &stats, CLOUD_TOLERANCE);
     787      plot_setMcal (psfStars->flxlist, Nref, &stats, CLOUD_TOLERANCE);
    804788    }
    805789
    … …  
    818802
    819803  StatDataSetFree (brightStars, 1);
    820   StatDataSetFree (refStars, 1);
     804  StatDataSetFree (psfStars, 1);
    821805
    822806  fprintf (stderr, "%d images calibrated\n", Ncalibrated);
    … …  
    10561040      c = ImageToCatalog[i][j];
    10571041
    1058       Mcal  = getMcal  (m, c, flatcorr, catalog);
     1042      Mcal  = getMcal  (m, c, MAG_CLASS_PSF);
    10591043      if (isnan(Mcal)) continue;
    10601044      Mmos  = getMmos  (m, c);

  • branches/eam_branches/ohana.20170822/src/relphot/src/MosaicOps.c

    r40212 r40213  
    925925  off_t Nmax;
    926926  int PoorImages;
    927   double *list;
    928   double *dlist;
    929   double *Mlist;
    930   double *dMlist;
     927  double *psfMagList;
     928  double *psfErrList;
     929  double *brightMagList;
     930  double *brightErrList;
     931  double *kronMagList;
     932  double *kronErrList;
    931933} SetMmosInfo;
    932934
    … …  
    959961
    960962  if (allocLists) {
    961     ALLOCATE (info->list, double, Nmax);
    962     ALLOCATE (info->dlist, double, Nmax);
    963     ALLOCATE (info->Mlist, double, Nmax);
    964     ALLOCATE (info->dMlist, double, Nmax);
     963    ALLOCATE (info->psfMagList, double, Nmax);
     964    ALLOCATE (info->psfErrList, double, Nmax);
     965    ALLOCATE (info->kronMagList, double, Nmax);
     966    ALLOCATE (info->kronErrList, double, Nmax);
     967    ALLOCATE (info->brightMagList, double, Nmax);
     968    ALLOCATE (info->brightErrList, double, Nmax);
    965969  }
    966970}
    967971
    968972void SetMmosInfoFree (SetMmosInfo *info) {
    969   free (info->list);
    970   free (info->dlist);
    971   free (info->Mlist);
    972   free (info->dMlist);
     973  free (info->psfMagList   );
     974  free (info->psfErrList   );
     975  free (info->kronMagList  );
     976  free (info->kronErrList  );
     977  free (info->brightMagList);
     978  free (info->brightErrList);
    973979}
    974980
    … …  
    10711077  liststats_setmode (&stats, "INNER_WTMEAN");
    10721078
    1073   double *list   = info->list;
    1074   double *dlist  = info->dlist;
    1075   double *Mlist  = info->Mlist;
    1076   double *dMlist = info->dMlist;
     1079  double *psfMagList    = info->psfMagList;
     1080  double *psfErrList    = info->psfErrList;
     1081  double *kronMagList   = info->kronMagList;
     1082  double *kronErrList   = info->kronErrList;
     1083  double *brightMagList = info->brightMagList;
     1084  double *brightErrList = info->brightErrList;
    10771085
    10781086  assert (Nmos >= 0);
    … …  
    11391147  int N = 0;
    11401148  for (j = 0; j < N_onMosaic[Nmos]; j++) {
    1141     float Msys, Mrel, Mcal, Mgrid, Mflat;
    11421149     
    11431150    off_t m = MosaicToMeasure[Nmos][j];
    … …  
    11461153    // NOTE : we are only using Mcal == McalPSF in this function; we set McalAPER to McalPSF
    11471154    if (fout) {
    1148       Mcal  = getMcal  (m, c, flatcorr, catalog);
    1149       Mgrid = getMgrid (m, c);
    1150       Mrel  = getMrel  (catalog, m, c, MAG_CLASS_PSF, MAG_SRC_CHP);
    1151       Mflat = getMflat (m, c, flatcorr, catalog);
     1155      float Mcal     = getMcal  (m, c, MAG_CLASS_PSF);
     1156      float Mgrid    = getMgrid (m, c);
     1157      float MrelPSF  = getMrel  (catalog, m, c, MAG_CLASS_PSF, MAG_SRC_CHP);
     1158      float Mflat    = getMflat (m, c, flatcorr, catalog);
    11521159
    11531160      off_t n = catalog[c].measureT[m].averef;
    1154       Msys = PhotSysTiny (&catalog[c].measureT[m], &catalog[c].averageT[n], &catalog[c].secfilt[n*Nsecfilt], MAG_CLASS_PSF);
    1155 
    1156       float delta = Msys - Mrel - Mcal - Mgrid + Mflat;
     1161      float MsysPSF = PhotSysTiny (&catalog[c].measureT[m], &catalog[c].averageT[n], &catalog[c].secfilt[n*Nsecfilt], MAG_CLASS_PSF);
     1162
     1163      float delta = MsysPSF - MrelPSF - Mcal - Mgrid + Mflat;
    11571164
    11581165      int isBad = (catalog[c].measureT[m].dbFlags & MEAS_BAD);
    11591166
    1160       fprintf (fout, "%f %f : %f %f %f %f %f  : %f %d\n", catalog[c].averageT[n].R, catalog[c].averageT[n].D, Msys, Mrel, Mcal, Mgrid, Mflat, delta, isBad);
     1167      fprintf (fout, "%f %f : %f %f %f %f %f  : %f %d\n", catalog[c].averageT[n].R, catalog[c].averageT[n].D, MsysPSF, MrelPSF, Mcal, Mgrid, Mflat, delta, isBad);
    11611168    }
    11621169
    … …  
    11651172      continue;
    11661173    }
    1167     Mcal  = getMcal  (m, c, flatcorr, catalog);
     1174    float Mcal = getMcal  (m, c, MAG_CLASS_PSF);
    11681175    if (isnan(Mcal)) {
    11691176      info->Ncal++;
    11701177      continue;
    11711178    }
    1172     Mgrid = getMgrid (m, c);
     1179    float Mgrid = getMgrid (m, c);
    11731180    if (isnan(Mgrid)) {
    11741181      info->Ngrid ++;
    11751182      continue;
    11761183    }
    1177     Mrel  = getMrel  (catalog, m, c, MAG_CLASS_PSF, MAG_SRC_CHP);
    1178     if (isnan(Mrel)) {
     1184
     1185    // Mrel* is the average magnitude for this star.  For PS1 stacks, we have too much
     1186    // PSF variability.  We need to calibrate the PSF magnitudes separately from the
     1187    // Aperture-like magnitues.  (We have an option to use the kron magnitudes or the
     1188    // other apertures here).  I basically need to do this analysis separately for each
     1189    // magnitude type
     1190   
     1191    float MrelPSF = getMrel  (catalog, m, c, MAG_CLASS_PSF, MAG_SRC_CHP);
     1192    if (isnan(MrelPSF)) {
    11791193      info->Nrel ++;
    11801194      continue;
    11811195    }
     1196    float MrelKron = getMrel  (catalog, m, c, MAG_CLASS_KRON, MAG_SRC_CHP);
    11821197     
    11831198    // image.Mcal is not supposed to include the flat-field correction, so we need to
    … …  
    11861201    // the flat-correction.  NOTE the sign of Mflat (Image.Mcal = Measure.Mcal - Mflat)
    11871202
    1188     Mflat = getMflat (m, c, flatcorr, catalog);
     1203    float Mflat = getMflat (m, c, flatcorr, catalog);
    11891204
    11901205    off_t n = catalog[c].measureT[m].averef;
    1191     Msys = PhotSysTiny (&catalog[c].measureT[m], &catalog[c].averageT[n], &catalog[c].secfilt[n*Nsecfilt], MAG_CLASS_PSF);
    1192     if (isnan(Msys)) {
     1206    float MsysPSF = PhotSysTiny (&catalog[c].measureT[m], &catalog[c].averageT[n], &catalog[c].secfilt[n*Nsecfilt], MAG_CLASS_PSF);
     1207    if (isnan(MsysPSF)) {
    11931208      info->Nsys++;
    11941209      continue;
    11951210    }
     1211    float MsysKron = PhotSysTiny (&catalog[c].measureT[m], &catalog[c].averageT[n], &catalog[c].secfilt[n*Nsecfilt], MAG_CLASS_KRON);
    11961212
    11971213    PhotCode *code = GetPhotcodebyCode (catalog[c].measureT[m].photcode);
    … …  
    12081224    assert (Nbright >= 0);
    12091225
    1210     list[N]  = Msys - Mrel - Mcal - Mgrid + Mflat;
    1211     dlist[N] = MAX (catalog[c].measureT[m].dM, MIN_ERROR);
     1226    psfMagList[N]  = MsysPSF - MrelPSF - Mcal - Mgrid + Mflat;
     1227    psfErrList[N]  = MAX (catalog[c].measureT[m].dM, MIN_ERROR);
     1228    kronMagList[N] = MsysKron - MrelKron - Mcal - Mgrid + Mflat;
     1229    kronErrList[N] = psfErrList[N];
    12121230    if (catalog[c].measureT[m].dM < IMFIT_SYS_SIGMA_LIM) {
    1213       Mlist[Nbright] = list[N];
    1214       dMlist[Nbright] = dlist[N];
     1231      brightMagList[Nbright] = psfMagList[N];
     1232      brightErrList[Nbright] = psfErrList[N];
    12151233      Nbright ++;
    12161234    }
    … …  
    12381256  }
    12391257
    1240   liststats (list, dlist, NULL, N, &stats);
     1258  liststats (psfMagList, psfErrList, NULL, N, &stats);
    12411259  if (VERBOSE2 && info->PoorImages) fprintf (stderr, "Mmos: %f %f %d %d\n", stats.mean, stats.sigma, stats.Nmeas, N);
    12421260
    12431261  // for now, I have no reason to measure these separately for camera-level images
    12441262  myMosaic[0].McalPSF    = stats.mean;
    1245   myMosaic[0].McalAPER   = stats.mean;
    12461263  myMosaic[0].dMcal      = stats.error;
    12471264  myMosaic[0].McalChiSq  = stats.chisq;
    12481265  myMosaic[0].nFitPhotom = N;
    12491266
     1267  liststats (kronMagList, kronErrList, NULL, N, &stats);
     1268  myMosaic[0].McalAPER   = stats.mean;
     1269
    12501270  if (testImage) {
    12511271    fprintf (stderr, "test image %d (%d) %f %f %d ... ", (int) Nmos, myMosaic[0].start, stats.mean, stats.error, myMosaic[0].nFitPhotom);
    … …  
    12541274  if (PLOTSTUFF) {
    12551275    fprintf (stderr, "Mmos: %6.3f %6.3f +/- %6.3f %5d %5d | %s\n", stats.mean, stats.median, stats.sigma, stats.Nmeas, N, image[MosaicToImage[Nmos][0]].name);
    1256     plot_setMcal (list, N, &stats, CLOUD_TOLERANCE);
     1276    plot_setMcal (psfMagList, N, &stats, CLOUD_TOLERANCE);
    12571277  }
    12581278
    12591279  // bright end scatter
    1260   liststats (Mlist, dMlist, NULL, Nbright, &stats);
     1280  liststats (brightMagList, brightErrList, NULL, Nbright, &stats);
    12611281  myMosaic[0].dMsys = stats.sigma;
    12621282
    … …  
    16941714      c = MosaicToCatalog[i][j];
    16951715
    1696       Mcal = getMcal  (m, c, flatcorr, catalog);
     1716      Mcal = getMcal  (m, c, MAG_CLASS_PSF);
    16971717      if (isnan(Mcal)) continue;
    16981718      Mgrid = getMgrid (m, c);

  • branches/eam_branches/ohana.20170822/src/relphot/src/StarOps.c

    r40212 r40213  
    484484      for (k = 0; k < catalog[i].averageT[j].Nmeasure; k++, m++) {
    485485        if (catalog[i].measureT[m].dbFlags & MEAS_BAD) continue;
    486         Mcal  = getMcal  (m, i, flatcorr, catalog);
     486        Mcal  = getMcal  (m, i, MAG_CLASS_PSF);
    487487        if (isnan(Mcal)) continue;
    488488        Mmos  = getMmos  (m, i);
    … …  
    688688         
    689689          // NOTE: we do not skip MEAS_BAD because this measurement is just an internal assessment of the outliers
    690           Mcal  = getMcal  (m, i, flatcorr, catalog);
     690          Mcal  = getMcal  (m, i, MAG_CLASS_PSF);
    691691          if (isnan(Mcal)) { Ncal ++; continue; }
    692692          Mmos  = getMmos  (m, i);
    … …  
    739739         
    740740          // NOTE: we do not skip MEAS_BAD because this measurement is just an internal assessment of the outliers
    741           Mcal  = getMcal  (m, i, flatcorr, catalog);
     741          Mcal  = getMcal  (m, i, MAG_CLASS_PSF);
    742742          if (isnan(Mcal)) continue;
    743743          Mmos  = getMmos  (m, i);
    … …  
    820820        int ecode = GetPhotcodeEquivCodebyCode (catalog[i].measureT[m].photcode);
    821821        if (ecode != seccode) { continue;}
    822         Mcal = getMcal  (m, i, flatcorr, catalog);
     822        Mcal = getMcal  (m, i, MAG_CLASS_PSF);
    823823        if (isnan(Mcal)) { continue;}
    824824        Mmos = getMmos  (m, i);

  • branches/eam_branches/ohana.20170822/src/relphot/src/plot_scatter.c

    r40212 r40213  
    4040
    4141                if (catalog[i].measureT[m].dbFlags & MEAS_BAD) continue;
    42                 Mcal = getMcal  (m, i, flatcorr, catalog);
     42                Mcal = getMcal  (m, i, MAG_CLASS_PSF);
    4343                if (isnan(Mcal)) continue;
    4444                Mmos = getMmos  (m, i);

  • branches/eam_branches/ohana.20170822/src/relphot/src/setMrelCatalog.c

    r40212 r40213  
    226226      Mcal = measureT[k].Mcal; // check that this is zero for loaded REF value
    227227    } else {
    228       // getMcal_alt returns image[].Mcal modified by flatcorr(image.photom_map_id,x,y)
    229       // NOTE: getMcal_alt does not include measure.Mflat
    230       Mcal  = getMcal_alt (meas, cat, flatcorr, measureT[k].Xccd, measureT[k].Yccd);
     228      // getMcal returns image[].Mcal; note the flat-field correction is stored in measure.Mflat
     229      Mcal  = getMcal (meas, cat, MAG_CLASS_PSF);
    231230      if (isnan(Mcal))  SKIP_THIS_MEAS(Ncal);
    232231      Mmos  = getMmos  (meas, cat);
    … …  
    530529  off_t k;
    531530
    532   float Mcal = 0, Mmos = 0, Mgrid = 0, Finst = 0;
     531  float McalPSF = 0, McalAPER = 0, Mmos = 0, Mgrid = 0, Finst = 0;
    533532
    534533  // set the primary projection cell and skycell for this coordinate
    … …  
    630629        // overlaps).  Msys + measure.Mcal is our best guess of the true magnitude
    631630        Mmos = Mgrid = 0;
    632         Mcal = measure[k].Mcal; // check that this is zero for loaded REF value
     631        McalPSF = measure[k].Mcal; // check that this is zero for loaded REF value
     632        McalAPER = McalPSF; // check that this is zero for loaded REF value
    633633      } else {
    634         Mcal  = getMcal_alt  (meas, cat, flatcorr, measure[k].Xccd, measure[k].Yccd);
    635         if (isnan(Mcal))  SKIP_THIS_MEAS_STACK(Ncal);
     634        McalPSF  = getMcal  (meas, cat, MAG_CLASS_PSF);
     635        if (isnan(McalPSF))  SKIP_THIS_MEAS_STACK(Ncal); // XXX really skip?
     636
     637        McalAPER  = getMcal  (meas, cat, MAG_CLASS_PSF);
     638        if (isnan(McalAPER))  SKIP_THIS_MEAS_STACK(Ncal); // XXX really skip?
     639
    636640        Mmos  = getMmos  (meas, cat);
    637641        if (isnan(Mmos))  SKIP_THIS_MEAS_STACK(Nmos);
    … …  
    672676
    673677    // get the zero point for the selected image
    674     float zp = PhotZeroPoint (&measure[k], &average[0], &secfilt[0]) - (Mcal + Mmos + Mgrid);
     678    // Use a different zero point for the PSF-like and APERTURE-like magnitudes
     679    float zpPSF  = PhotZeroPoint (&measure[k], &average[0], &secfilt[0]) - (McalPSF  + Mmos + Mgrid);
     680    float zpAPER = PhotZeroPoint (&measure[k], &average[0], &secfilt[0]) - (McalAPER + Mmos + Mgrid);
    675681
    676682    // flux_cgs : erg sec^1 cm^-2 Hz^-1
    … …  
    704710
    705711    // zpFactor to go from instrumental flux to Janskies
    706     float zpFactor = pow(10.0, -0.4*zp + 3.56);
     712    float zpFactorPSF  = pow(10.0, -0.4*zpPSF  + 3.56);
     713    float zpFactorAPER = pow(10.0, -0.4*zpAPER + 3.56);
    707714
    708715    // need to put in AB mag factor to get to Janskies (or uJy?)
    709     secfilt[Nsec].FpsfStk   = zpFactor * measure[k].FluxPSF; 
    710     secfilt[Nsec].dFpsfStk  = zpFactor * measure[k].dFluxPSF;
    711     secfilt[Nsec].FkronStk  = zpFactor * measure[k].FluxKron;
    712     secfilt[Nsec].dFkronStk = zpFactor * measure[k].dFluxKron;
    713     secfilt[Nsec].FapStk    = zpFactor * measure[k].FluxAp;
     716    secfilt[Nsec].FpsfStk   = zpFactorPSF * measure[k].FluxPSF; 
     717    secfilt[Nsec].dFpsfStk  = zpFactorPSF * measure[k].dFluxPSF;
     718    secfilt[Nsec].FkronStk  = zpFactorAPER * measure[k].FluxKron;
     719    secfilt[Nsec].dFkronStk = zpFactorAPER * measure[k].dFluxKron;
     720    secfilt[Nsec].FapStk    = zpFactorAPER * measure[k].FluxAp;
    714721
    715722    // NOTE: for PV3, apFluxErr is broken (see pmSourcePhotometry.c:245).  we are going to
    716723    // use PSF flux error instead here:
    717     // secfilt[Nsec].dFapStk   = zpFactor * measure[k].dFluxAp;
    718     secfilt[Nsec].dFapStk   = zpFactor * measure[k].dFluxPSF;
     724    // secfilt[Nsec].dFapStk   = zpFactorAPER * measure[k].dFluxAp;
     725    secfilt[Nsec].dFapStk   = zpFactorAPER * measure[k].dFluxPSF;
    719726
    720727    // Jy to AB mags
    … …  
    859866      Mcal = measure[k].Mcal; // check that this is zero for loaded REF value
    860867    } else {
    861       // use getMcal not getMcal_alt?
    862       Mcal  = getMcal_alt (meas, cat, NULL, measure[k].Xccd, measure[k].Yccd);
    863       // Mcal  = getMcal (meas, cat);
     868      Mcal  = getMcal (meas, cat, MAG_CLASS_PSF);
    864869      if (isnan(Mcal))  SKIP_THIS_MEAS(Ncal);
    865870    }

  • branches/eam_branches/ohana.20170822/src/relphot/src/setMrelFinal.c

    r40212 r40213  
    226226        off_t Nim = getImageEntry (m, 0);
    227227        if (Nim > -1) {
    228           if (isnan(getMcal (m, 0, flatcorr, catalog))) goto skip;
     228          if (isnan(getMcal (m, 0, MAG_CLASS_PSF))) goto skip;
    229229          if (isnan(getMmos (m, 0))) goto skip;
    230230        }
Note: See TracChangeset for help on using the changeset viewer.