IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 42695


Ignore:
Timestamp:
Jun 24, 2024, 8:59:50 AM (2 years ago)
Author:
eugene
Message:

handle both 3pi-style convolved radial apertures and single (unconvolved) radial apertures; set initial average mags correctly to stack and warp fields, depending on the photcodes

Location:
branches/eam_branches/ipp-20230313/Ohana/src/addstar
Files:
1 added
5 edited

Legend:

Unmodified
Added
Removed

  • branches/eam_branches/ipp-20230313/Ohana/src/addstar/Makefile

    r42389 r42695  
    9696$(SRC)/find_matches_refstars.$(ARCH).o \
    9797$(SRC)/find_matches_closest_refstars.$(ARCH).o \
     98$(SRC)/photcodes.$(ARCH).o \
    9899$(SRC)/getgsc.$(ARCH).o \
    99100$(SRC)/getusno.$(ARCH).o \
    … …  
    151152$(SRC)/find_matches_refstars.$(ARCH).o \
    152153$(SRC)/find_subset.$(ARCH).o \
     154$(SRC)/photcodes.$(ARCH).o \
    153155$(SRC)/getgsc.$(ARCH).o \
    154156$(SRC)/getusno.$(ARCH).o \
    … …  
    181183$(SRC)/find_matches_refstars.$(ARCH).o \
    182184$(SRC)/find_subset.$(ARCH).o \
     185$(SRC)/photcodes.$(ARCH).o \
    183186$(SRC)/getgsc.$(ARCH).o \
    184187$(SRC)/getusno.$(ARCH).o \

  • branches/eam_branches/ipp-20230313/Ohana/src/addstar/include/addstar.h

    r41475 r42695  
    141141int     OLD_RESORT;
    142142int     READ_XRAD_DATA;
     143int     CONVOLVED_APERTURES;
    143144int     DIFF_WITH_INV;
    144145
    … …  
    343344void GetConfig (char *config, char *field, char *format, int N, void *ptr);
    344345
     346int isGPC1warp (int photcode);
     347int isGPC1stack (int photcode);
     348
    345349/**
    346350    there is an inconsistency to be resolved: fixed structures (like Image)

  • branches/eam_branches/ipp-20230313/Ohana/src/addstar/src/ReadXradFITS.c

    r42077 r42695  
    9797  // first entry and compare to the rest
    9898
    99   float fwhmValues[3];
    100   fwhmValues[0] = PSFfwhm[0];
    101   fwhmValues[1] = PSFfwhm[1];
    102   fwhmValues[2] = PSFfwhm[2];
    103 
    104   int i;
    10599  int Nap = 0;
    106   for (i = 0; i < catalog->Nmeasure; i++) {
    107     dvo_lensing_init (&catalog->lensing[i]);
    108 
    109     if (catalog->measure[i].detID < RadID[Nap]) {
    110       continue;
     100
     101  if (CONVOLVED_APERTURES) {
     102    // XXX code to support the 3 convolved apertures
     103
     104    float fwhmValues[3];
     105    fwhmValues[0] = PSFfwhm[0];
     106    fwhmValues[1] = PSFfwhm[1];
     107    fwhmValues[2] = PSFfwhm[2];
     108
     109    for (int i = 0; i < catalog->Nmeasure; i++) {
     110      dvo_lensing_init (&catalog->lensing[i]);
     111
     112      if (catalog->measure[i].detID < RadID[Nap]) {
     113        continue;
     114        // this is a psf measurement which does not have a radial aperture
     115      }
     116      if (catalog->measure[i].detID > RadID[Nap]) {
     117        myAbort("radial apertures for source not in psf list? sources out of order?  seems like a bug\n");
     118        // this could be a radial aperture which does not have a PSF source, but that is not possible
     119      }
     120
     121      // confirm the 3 FWHM values:
     122      myAssert (fwhmValues[0] == PSFfwhm[Nap+0], "FWHM mismatch %f vs %f", fwhmValues[0], PSFfwhm[Nap+0]);
     123      myAssert (fwhmValues[1] == PSFfwhm[Nap+1], "FWHM mismatch %f vs %f", fwhmValues[1], PSFfwhm[Nap+1]);
     124      myAssert (fwhmValues[2] == PSFfwhm[Nap+2], "FWHM mismatch %f vs %f", fwhmValues[2], PSFfwhm[Nap+2]);
     125
     126      // EAM 2022.02.17 : here is the comment for the PV3 load:
     127      // XXX this is all hard-wired and should make use of the headers.
     128      // psphot cmfs have 5 radial apertures:
     129      // array 0, 1, 2, 3, 4
     130      // SDSS  3, 4, 5, 6, 7
     131
     132      // EAM 2022.02.17 : here is the situation for UNIONS DR3:
     133      // we have 3 convolutions (raw, 6", 8")
     134      // for each we have 6 apertures with max radii of (4, 8, 16, 32, 48, 64) pixels = (1, 2, 4, 8, 12, 16) arcsec
     135      // I am going to save (4, 16, 32) which have index of (0, 2, 3)
     136
     137# define RAD_0 0
     138# define RAD_1 2
     139# define RAD_2 3
     140      catalog->lensing[i]. F_ApR5    = AperFlux   [(Nap + 0)*Ncol + RAD_0];
     141      catalog->lensing[i].dF_ApR5    = AperFluxErr[(Nap + 0)*Ncol + RAD_0];
     142      catalog->lensing[i].sF_ApR5    = AperFluxStd[(Nap + 0)*Ncol + RAD_0];
     143      catalog->lensing[i].fF_ApR5    = AperFill   [(Nap + 0)*Ncol + RAD_0];
     144                                   
     145      catalog->lensing[i]. F_ApR6    = AperFlux   [(Nap + 0)*Ncol + RAD_1];
     146      catalog->lensing[i].dF_ApR6    = AperFluxErr[(Nap + 0)*Ncol + RAD_1];
     147      catalog->lensing[i].sF_ApR6    = AperFluxStd[(Nap + 0)*Ncol + RAD_1];
     148      catalog->lensing[i].fF_ApR6    = AperFill   [(Nap + 0)*Ncol + RAD_1];
     149                                   
     150      catalog->lensing[i]. F_ApR7    = AperFlux   [(Nap + 0)*Ncol + RAD_2];
     151      catalog->lensing[i].dF_ApR7    = AperFluxErr[(Nap + 0)*Ncol + RAD_2];
     152      catalog->lensing[i].sF_ApR7    = AperFluxStd[(Nap + 0)*Ncol + RAD_2];
     153      catalog->lensing[i].fF_ApR7    = AperFill   [(Nap + 0)*Ncol + RAD_2];
     154
     155      catalog->lensing[i]. F_ApR5_C1 = AperFlux   [(Nap + 1)*Ncol + RAD_0];
     156      catalog->lensing[i].dF_ApR5_C1 = AperFluxErr[(Nap + 1)*Ncol + RAD_0];
     157      catalog->lensing[i]. F_ApR6_C1 = AperFlux   [(Nap + 1)*Ncol + RAD_1];
     158      catalog->lensing[i].dF_ApR6_C1 = AperFluxErr[(Nap + 1)*Ncol + RAD_1];
     159      catalog->lensing[i]. F_ApR7_C1 = AperFlux   [(Nap + 1)*Ncol + RAD_2];
     160      catalog->lensing[i].dF_ApR7_C1 = AperFluxErr[(Nap + 1)*Ncol + RAD_2];
     161
     162      catalog->lensing[i]. F_ApR5_C2 = AperFlux   [(Nap + 2)*Ncol + RAD_0];
     163      catalog->lensing[i].dF_ApR5_C2 = AperFluxErr[(Nap + 2)*Ncol + RAD_0];
     164      catalog->lensing[i]. F_ApR6_C2 = AperFlux   [(Nap + 2)*Ncol + RAD_1];
     165      catalog->lensing[i].dF_ApR6_C2 = AperFluxErr[(Nap + 2)*Ncol + RAD_1];
     166      catalog->lensing[i]. F_ApR7_C2 = AperFlux   [(Nap + 2)*Ncol + RAD_2];
     167      catalog->lensing[i].dF_ApR7_C2 = AperFluxErr[(Nap + 2)*Ncol + RAD_2];
     168
     169      catalog->lensing[i].detID = catalog->measure[i].detID;
     170
     171      // XXX set the measure, object, etc ID values here
     172      // catalog->lensing[i].objID : set in find_matches_closest.c
     173      // catalog->lensing[i].catID : set in find_matches_closest.c
     174      // catalog->lensing[i].averef : set in find_matches_closest.c
     175      // catalog->lensing[i].imageID : set in FilterStars.c and UpdateImageIDs.c
     176
     177      Nap += 3;
     178    }
     179  } else {
     180    // PSF values for unconvolve data are undefined ?
     181
     182    for (int i = 0; i < catalog->Nmeasure; i++) {
     183      dvo_lensing_init (&catalog->lensing[i]);
     184
    111185      // this is a psf measurement which does not have a radial aperture
     186      if (catalog->measure[i].detID < RadID[Nap]) { continue; }
     187
     188      // this could be a radial aperture which does not have a PSF source, but that is not possible
     189      if (catalog->measure[i].detID > RadID[Nap]) {
     190        myAbort("radial apertures for source not in psf list? sources out of order?  seems like a bug\n");
     191      }
     192
     193      // EAM 2022.02.17 : here is the comment for the PV3 load:
     194      // XXX this is all hard-wired and should make use of the headers.
     195      // psphot cmfs have 5 radial apertures:
     196      // array 0, 1, 2, 3, 4
     197      // SDSS  3, 4, 5, 6, 7
     198
     199      // EAM 2022.02.17 : here is the situation for UNIONS DR4:
     200      // we have 1 convolution (raw)
     201      // for each we have 6 apertures with max radii of (4, 8, 16, 32, 48, 64) pixels = (1, 2, 4, 8, 12, 16) arcsec
     202      // I am going to save (4, 16, 32) which have index of (0, 2, 3)
     203
     204# define RAD_0 0
     205# define RAD_1 2
     206# define RAD_2 3
     207      catalog->lensing[i]. F_ApR5    = AperFlux   [(Nap + 0)*Ncol + RAD_0];
     208      catalog->lensing[i].dF_ApR5    = AperFluxErr[(Nap + 0)*Ncol + RAD_0];
     209      catalog->lensing[i].sF_ApR5    = AperFluxStd[(Nap + 0)*Ncol + RAD_0];
     210      catalog->lensing[i].fF_ApR5    = AperFill   [(Nap + 0)*Ncol + RAD_0];
     211                                   
     212      catalog->lensing[i]. F_ApR6    = AperFlux   [(Nap + 0)*Ncol + RAD_1];
     213      catalog->lensing[i].dF_ApR6    = AperFluxErr[(Nap + 0)*Ncol + RAD_1];
     214      catalog->lensing[i].sF_ApR6    = AperFluxStd[(Nap + 0)*Ncol + RAD_1];
     215      catalog->lensing[i].fF_ApR6    = AperFill   [(Nap + 0)*Ncol + RAD_1];
     216                                   
     217      catalog->lensing[i]. F_ApR7    = AperFlux   [(Nap + 0)*Ncol + RAD_2];
     218      catalog->lensing[i].dF_ApR7    = AperFluxErr[(Nap + 0)*Ncol + RAD_2];
     219      catalog->lensing[i].sF_ApR7    = AperFluxStd[(Nap + 0)*Ncol + RAD_2];
     220      catalog->lensing[i].fF_ApR7    = AperFill   [(Nap + 0)*Ncol + RAD_2];
     221
     222      catalog->lensing[i].detID = catalog->measure[i].detID;
     223
     224      // XXX set the measure, object, etc ID values here
     225      // catalog->lensing[i].objID : set in find_matches_closest.c
     226      // catalog->lensing[i].catID : set in find_matches_closest.c
     227      // catalog->lensing[i].averef : set in find_matches_closest.c
     228      // catalog->lensing[i].imageID : set in FilterStars.c and UpdateImageIDs.c
     229
     230      Nap += 1;
    112231    }
    113     if (catalog->measure[i].detID > RadID[Nap]) {
    114       myAbort("radial apertures for source not in psf list? sources out of order?  seems like a bug\n");
    115       // this could be a radial aperture which does not have a PSF source, but that is not possible
    116     }
    117 
    118     // confirm the 3 FWHM values:
    119     myAssert (fwhmValues[0] == PSFfwhm[Nap+0], "FWHM mismatch %f vs %f", fwhmValues[0], PSFfwhm[Nap+0]);
    120     myAssert (fwhmValues[1] == PSFfwhm[Nap+1], "FWHM mismatch %f vs %f", fwhmValues[1], PSFfwhm[Nap+1]);
    121     myAssert (fwhmValues[2] == PSFfwhm[Nap+2], "FWHM mismatch %f vs %f", fwhmValues[2], PSFfwhm[Nap+2]);
    122 
    123     // EAM 2022.02.17 : here is the comment for the PV3 load:
    124     // XXX this is all hard-wired and should make use of the headers.
    125     // psphot cmfs have 5 radial apertures:
    126     // array 0, 1, 2, 3, 4
    127     // SDSS  3, 4, 5, 6, 7
    128 
    129     // EAM 2022.02.17 : here is the situation for UNIONS DR3:
    130     // we have 3 convolutions (raw, 6", 8")
    131     // for each we have 6 apertures with max radii of (4, 8, 16, 32, 48, 64) pixels = (1, 2, 4, 8, 12, 16) arcsec
    132     // I am going to save (4, 16, 32) which have index of (0, 2, 3)
    133 
    134     # define RAD_0 0
    135     # define RAD_1 2
    136     # define RAD_2 3
    137     catalog->lensing[i]. F_ApR5    = AperFlux   [(Nap + 0)*Ncol + RAD_0];
    138     catalog->lensing[i].dF_ApR5    = AperFluxErr[(Nap + 0)*Ncol + RAD_0];
    139     catalog->lensing[i].sF_ApR5    = AperFluxStd[(Nap + 0)*Ncol + RAD_0];
    140     catalog->lensing[i].fF_ApR5    = AperFill   [(Nap + 0)*Ncol + RAD_0];
    141                                    
    142     catalog->lensing[i]. F_ApR6    = AperFlux   [(Nap + 0)*Ncol + RAD_1];
    143     catalog->lensing[i].dF_ApR6    = AperFluxErr[(Nap + 0)*Ncol + RAD_1];
    144     catalog->lensing[i].sF_ApR6    = AperFluxStd[(Nap + 0)*Ncol + RAD_1];
    145     catalog->lensing[i].fF_ApR6    = AperFill   [(Nap + 0)*Ncol + RAD_1];
    146                                    
    147     catalog->lensing[i]. F_ApR7    = AperFlux   [(Nap + 0)*Ncol + RAD_2];
    148     catalog->lensing[i].dF_ApR7    = AperFluxErr[(Nap + 0)*Ncol + RAD_2];
    149     catalog->lensing[i].sF_ApR7    = AperFluxStd[(Nap + 0)*Ncol + RAD_2];
    150     catalog->lensing[i].fF_ApR7    = AperFill   [(Nap + 0)*Ncol + RAD_2];
    151 
    152     catalog->lensing[i]. F_ApR5_C1 = AperFlux   [(Nap + 1)*Ncol + RAD_0];
    153     catalog->lensing[i].dF_ApR5_C1 = AperFluxErr[(Nap + 1)*Ncol + RAD_0];
    154     catalog->lensing[i]. F_ApR6_C1 = AperFlux   [(Nap + 1)*Ncol + RAD_1];
    155     catalog->lensing[i].dF_ApR6_C1 = AperFluxErr[(Nap + 1)*Ncol + RAD_1];
    156     catalog->lensing[i]. F_ApR7_C1 = AperFlux   [(Nap + 1)*Ncol + RAD_2];
    157     catalog->lensing[i].dF_ApR7_C1 = AperFluxErr[(Nap + 1)*Ncol + RAD_2];
    158 
    159     catalog->lensing[i]. F_ApR5_C2 = AperFlux   [(Nap + 2)*Ncol + RAD_0];
    160     catalog->lensing[i].dF_ApR5_C2 = AperFluxErr[(Nap + 2)*Ncol + RAD_0];
    161     catalog->lensing[i]. F_ApR6_C2 = AperFlux   [(Nap + 2)*Ncol + RAD_1];
    162     catalog->lensing[i].dF_ApR6_C2 = AperFluxErr[(Nap + 2)*Ncol + RAD_1];
    163     catalog->lensing[i]. F_ApR7_C2 = AperFlux   [(Nap + 2)*Ncol + RAD_2];
    164     catalog->lensing[i].dF_ApR7_C2 = AperFluxErr[(Nap + 2)*Ncol + RAD_2];
    165 
    166     catalog->lensing[i].detID = catalog->measure[i].detID;
    167 
    168     // XXX set the measure, object, etc ID values here
    169     // catalog->lensing[i].objID : set in find_matches_closest.c
    170     // catalog->lensing[i].catID : set in find_matches_closest.c
    171     // catalog->lensing[i].averef : set in find_matches_closest.c
    172     // catalog->lensing[i].imageID : set in FilterStars.c and UpdateImageIDs.c
    173 
    174     Nap += 3;
    175232  }
    176233  myAssert (Nap == Nrow, "did we go too far???");

  • branches/eam_branches/ipp-20230313/Ohana/src/addstar/src/args.c

    r41475 r42695  
    9191    remove_argument (N, &argc, argv);
    9292    READ_XRAD_DATA = TRUE;
     93  }
     94  CONVOLVED_APERTURES = FALSE;
     95  if ((N = get_argument (argc, argv, "-convolved-ap"))) {
     96    remove_argument (N, &argc, argv);
     97    CONVOLVED_APERTURES = TRUE;
    9398  }
    9499  DIFF_WITH_INV = FALSE;

  • branches/eam_branches/ipp-20230313/Ohana/src/addstar/src/find_matches_closest.c

    r40291 r42695  
    248248      /* in UPDATE mode, this value is not saved; use relphot to recalculate */
    249249      if (Nsec > -1) {
     250        if (isGPC1warp(tgtcat[0].measure[Nmeas].photcode)) {
     251          if (FALSE && isnan(tgtcat[0].secfilt[n*Nsecfilt+Nsec].MpsfWrp)) {
     252            tgtcat[0].secfilt[n*Nsecfilt+Nsec].MpsfWrp  = PhotCat (&tgtcat[0].measure[Nmeas], MAG_CLASS_PSF);
     253            tgtcat[0].secfilt[n*Nsecfilt+Nsec].MkronWrp = PhotCat (&tgtcat[0].measure[Nmeas], MAG_CLASS_KRON);
     254            tgtcat[0].secfilt[n*Nsecfilt+Nsec].MapWrp   = PhotCat (&tgtcat[0].measure[Nmeas], MAG_CLASS_APER);
     255          }
     256          goto set_average;
     257        }
     258        if (isGPC1stack(tgtcat[0].measure[Nmeas].photcode)) {
     259          if (FALSE && isnan(tgtcat[0].secfilt[n*Nsecfilt+Nsec].MpsfStk)) {
     260            tgtcat[0].secfilt[n*Nsecfilt+Nsec].MpsfStk  = PhotCat (&tgtcat[0].measure[Nmeas], MAG_CLASS_PSF);
     261            tgtcat[0].secfilt[n*Nsecfilt+Nsec].MkronStk = PhotCat (&tgtcat[0].measure[Nmeas], MAG_CLASS_KRON);
     262            tgtcat[0].secfilt[n*Nsecfilt+Nsec].MapStk   = PhotCat (&tgtcat[0].measure[Nmeas], MAG_CLASS_APER);
     263          }
     264          goto set_average;
     265        }
    250266        if (isnan(tgtcat[0].secfilt[n*Nsecfilt+Nsec].MpsfChp)) {
    251           tgtcat[0].secfilt[n*Nsecfilt+Nsec].MpsfChp = PhotCat (&tgtcat[0].measure[Nmeas], MAG_CLASS_PSF);
    252         }
    253       }
     267          tgtcat[0].secfilt[n*Nsecfilt+Nsec].MpsfChp  = PhotCat (&tgtcat[0].measure[Nmeas], MAG_CLASS_PSF);
     268          tgtcat[0].secfilt[n*Nsecfilt+Nsec].MkronChp = PhotCat (&tgtcat[0].measure[Nmeas], MAG_CLASS_KRON);
     269          tgtcat[0].secfilt[n*Nsecfilt+Nsec].MapChp   = PhotCat (&tgtcat[0].measure[Nmeas], MAG_CLASS_APER);
     270        }
     271      }
     272    set_average:
    254273
    255274      /* Nm is updated, but not written out in -update mode (for existing entries)
    … …  
    352371      /* in UPDATE mode, this value is not saved; use relphot to recalculate */
    353372      if (Nsec > -1) {
    354         tgtcat[0].secfilt[Nave*Nsecfilt+Nsec].MpsfChp = PhotCat (&tgtcat[0].measure[Nmeas], MAG_CLASS_PSF);
    355       }
     373        if (isGPC1warp(tgtcat[0].measure[Nmeas].photcode)) {
     374          // tgtcat[0].secfilt[Nave*Nsecfilt+Nsec].MpsfWrp  = PhotCat (&tgtcat[0].measure[Nmeas], MAG_CLASS_PSF);
     375          // tgtcat[0].secfilt[Nave*Nsecfilt+Nsec].MkronWrp = PhotCat (&tgtcat[0].measure[Nmeas], MAG_CLASS_KRON);
     376          // tgtcat[0].secfilt[Nave*Nsecfilt+Nsec].MapWrp   = PhotCat (&tgtcat[0].measure[Nmeas], MAG_CLASS_APER);
     377          goto set_average_2;
     378        }
     379        if (isGPC1stack(tgtcat[0].measure[Nmeas].photcode)) {
     380          // tgtcat[0].secfilt[Nave*Nsecfilt+Nsec].MpsfStk  = PhotCat (&tgtcat[0].measure[Nmeas], MAG_CLASS_PSF);
     381          // tgtcat[0].secfilt[Nave*Nsecfilt+Nsec].MkronStk = PhotCat (&tgtcat[0].measure[Nmeas], MAG_CLASS_KRON);
     382          // tgtcat[0].secfilt[Nave*Nsecfilt+Nsec].MapStk   = PhotCat (&tgtcat[0].measure[Nmeas], MAG_CLASS_APER);
     383          goto set_average_2;
     384        }
     385        if (isnan(tgtcat[0].secfilt[Nave*Nsecfilt+Nsec].MpsfChp)) {
     386          tgtcat[0].secfilt[Nave*Nsecfilt+Nsec].MpsfChp  = PhotCat (&tgtcat[0].measure[Nmeas], MAG_CLASS_PSF);
     387          tgtcat[0].secfilt[Nave*Nsecfilt+Nsec].MkronChp = PhotCat (&tgtcat[0].measure[Nmeas], MAG_CLASS_KRON);
     388          tgtcat[0].secfilt[Nave*Nsecfilt+Nsec].MapChp   = PhotCat (&tgtcat[0].measure[Nmeas], MAG_CLASS_APER);
     389        }
     390      }
     391    set_average_2:
    356392
    357393      Nmeas ++;
Note: See TracChangeset for help on using the changeset viewer.