IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 31273


Ignore:
Timestamp:
Apr 12, 2011, 3:10:26 PM (15 years ago)
Author:
eugene
Message:

use MeasureTiny and AverageTiny structures for the photometry analysis to save memory

Location:
branches/eam_branches/ipp-20110404/Ohana/src/relphot
Files:
11 edited

Legend:

Unmodified
Added
Removed

  • branches/eam_branches/ipp-20110404/Ohana/src/relphot/include/relphot.h

    r31272 r31273  
    234234int rationalize_mosaics ();
    235235int LimitDensityCatalog (Catalog *subcatalog, Catalog *catalog);
     236
     237int populate_tiny_values (Catalog *catalog);
     238int free_tiny_values (Catalog *catalog);

  • branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/GridOps.c

    r31272 r31273  
    366366      mx = mlist[i][j];
    367367      c  = clist[i][j];
    368       n  = catalog[c].measure[mx].averef;
     368      n  = catalog[c].measureT[mx].averef;
    369369     
    370370      // if we have already visited this star, skip the stuff below
    … …  
    393393
    394394        // skip measurements marked as BAD
    395         if (catalog[c].measure[m].dbFlags & MEAS_BAD) {
     395        if (catalog[c].measureT[m].dbFlags & MEAS_BAD) {
    396396          Nbad ++;
    397397          continue;
    … …  
    414414        // select the color- and airmass-corrected observed magnitude for this star
    415415        // XXX need to be able to turn off the color-correction until initial average mags are found
    416         Msys = PhotCat (&catalog[c].measure[m]);
     416        Msys = PhotCatTiny (&catalog[c].measureT[m]);
    417417        if (isnan(Msys)) {
    418418          Nsys++;
    … …  
    421421
    422422        // mag-error for this measurement
    423         Merr =  MAX (catalog[c].measure[m].dM, MIN_ERROR);
     423        Merr =  MAX (catalog[c].measureT[m].dM, MIN_ERROR);
    424424
    425425        // Wsys = 1.0 / SQ(Merr);
    … …  
    550550      c = clist[i][j];
    551551     
    552       if (catalog[c].measure[m].dbFlags & MEAS_BAD) {
     552      if (catalog[c].measureT[m].dbFlags & MEAS_BAD) {
    553553        Nbad ++;
    554554        continue;
    … …  
    570570      }
    571571     
    572       n = catalog[c].measure[m].averef;
    573       Msys = PhotSys (&catalog[c].measure[m], &catalog[c].average[n], &catalog[c].secfilt[n*Nsecfilt]);
     572      n = catalog[c].measureT[m].averef;
     573      Msys = PhotSysTiny (&catalog[c].measureT[m], &catalog[c].averageT[n], &catalog[c].secfilt[n*Nsecfilt]);
    574574      if (isnan(Msys)) {
    575575        Nsys++;
    … …  
    577577      }
    578578      list[N] = Msys - Mrel - Mcal - Mmos;
    579       dlist[N] = MAX (catalog[c].measure[m].dM, MIN_ERROR);
     579      dlist[N] = MAX (catalog[c].measureT[m].dM, MIN_ERROR);
    580580      N++;
    581581    }
    … …  
    637637      c = clist[i][j];
    638638     
    639       if (catalog[c].measure[m].dbFlags & MEAS_BAD) {
     639      if (catalog[c].measureT[m].dbFlags & MEAS_BAD) {
    640640        Narea ++;
    641641        continue;
    … …  
    648648      if (isnan(Mrel)) continue;
    649649
    650       n = catalog[c].measure[m].averef;
    651       Msys = PhotSys (&catalog[c].measure[m], &catalog[c].average[n], &catalog[c].secfilt[n*Nsecfilt]);
     650      n = catalog[c].measureT[m].averef;
     651      Msys = PhotSysTiny (&catalog[c].measureT[m], &catalog[c].averageT[n], &catalog[c].secfilt[n*Nsecfilt]);
    652652
    653653      xlist[N] = Xmeas[c][m];

  • branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/ImageOps.c

    r31272 r31273  
    11# include "relphot.h"
    2 
    3 # define USE_IMAGE_ID 1
    42
    53static off_t       **bin;     // link from catalog,measure to image         
    … …  
    1210static off_t        Nimage;   // number of available images
    1311
    14 // if we search by image ID, we sort (imageIDs, imageIdx) by imageIDs to get a sorted
    15 // index
    16 
    17 # if USE_IMAGE_ID
     12// to search by image ID, we sort (imageIDs, imageIdx) by imageIDs to get a sorted index
    1813static off_t        *imageIDs; // list of all image IDs
    1914static off_t        *imageIdx; // list of index for image IDs
    20 # else
    21 static unsigned int *start;
    22 static unsigned int *stop;
    23 # endif
    2415
    2516Image *getimages (off_t *N) {
    … …  
    4031  Nimage = N;
    4132
    42 # if USE_IMAGE_ID
    4333  ALLOCATE (imageIDs, off_t, Nimage);
    4434  ALLOCATE (imageIdx, off_t, Nimage);
    … …  
    4939  }
    5040  llsortpair (imageIDs, imageIdx, Nimage);
    51 # else
    52   ALLOCATE (start,   unsigned, Nimage);
    53   ALLOCATE (stop,    unsigned, Nimage);
    54 
    55   for (i = 0; i < Nimage; i++) {
    56     start[i] = image[i].tzero - MAX(0.01*image[i].trate*image[i].NY, 1);
    57     stop[i]  = image[i].tzero + MAX(1.01*image[i].trate*image[i].NY, 1);
    58   }
    59 # endif
    6041}
    6142
    … …  
    6546  // use bisection to find the specified image ID
    6647
    67 # if USE_IMAGE_ID
    6848  off_t Nlo, Nhi, N;
    6949
    … …  
    8262      return (imageIdx[N]);
    8363  }
    84 # endif
    8564
    8665  return (-1);
    … …  
    137116 for (i = 0; i < Ncatalog; i++) {
    138117    for (j = 0; j < catalog[i].Nmeasure; j++) {
    139       ecode = GetPhotcodeEquivCodebyCode (catalog[i].measure[j].photcode);
     118      ecode = GetPhotcodeEquivCodebyCode (catalog[i].measureT[j].photcode);
    140119      found = FALSE;
    141120      for (Ns = 0; !found && (Ns < Nphotcodes); Ns++) {
    … …  
    150129}
    151130
    152 int findCCD (off_t idx, off_t meas, int cat, Measure *measure) {
     131int findCCD (off_t idx, off_t meas, int cat, MeasureTiny *measure) {
    153132
    154133  int ccdnum, found, Ns, ecode;
    … …  
    186165
    187166  // old code to add this measurement to the grid cell for this chip
    188   // ave = measure[0].averef;
    189   // ra  = catalog[cat].average[ave].R - measure[0].dR / 3600.0;
    190   // dec = catalog[cat].average[ave].D - measure[0].dD / 3600.0;
     167  // ave = measureT[0].averef;
     168  // ra  = catalog[cat].averageT[ave].R - measureT[0].dR / 3600.0;
     169  // dec = catalog[cat].averageT[ave].D - measureT[0].dD / 3600.0;
    191170  // RD_to_XY (&X, &Y, ra, dec, &image[i].coords);
    192171
    … …  
    199178}
    200179
    201 # if USE_IMAGE_ID
    202180void matchImage (Catalog *catalog, off_t meas, int cat) {
    203181
    204182  off_t idx, ID;
    205183  int status;
    206   Measure *measure;
     184  MeasureTiny *measure;
    207185 
    208   measure = &catalog[cat].measure[meas];
     186  measure = &catalog[cat].measureT[meas];
    209187
    210188  ID = measure[0].imageID;
    … …  
    241219  return;
    242220}
    243 
    244 # else
    245 // this is the time-based match
    246 void matchImage (Catalog *catalog, off_t meas, int cat) {
    247 
    248   off_t i;
    249   int status;
    250   Measure *measure;
    251  
    252   measure = &catalog[cat].measure[meas];
    253   for (i = 0; i < Nimage; i++) {
    254     if (image[0].photcode == -1) continue;
    255     if (measure[0].photcode != image[i].photcode) continue;
    256     if (measure[0].t < start[i]) continue;
    257     if (measure[0].t > stop[i]) continue;
    258    
    259     if (USE_GRID) {
    260       status = findCCD(i, meas, cat);
    261       if (!status) continue;
    262     }
    263 
    264     // index for (catalog, measure) -> image
    265     bin[cat][meas] = i;
    266 
    267     // index for image, Nentry -> catalog
    268     clist[i][Nlist[i]] = cat;
    269 
    270     // index for image, Nentry -> measure
    271     mlist[i][Nlist[i]] = meas;
    272     Nlist[i] ++;
    273 
    274     if (Nlist[i] == NLIST[i]) {
    275       NLIST[i] += 100;
    276       REALLOCATE (clist[i], off_t, NLIST[i]);
    277       REALLOCATE (mlist[i], off_t, NLIST[i]);
    278     }   
    279     return;
    280   }
    281   /*  fprintf (stderr, "can't find source image for this measurement: %d (%d)\n", measure[0].t, measure[0].photcode); */
    282 }
    283 # endif
    284221
    285222off_t getImageEntry (off_t meas, int cat) {
    … …  
    357294      c = clist[i][j];
    358295     
    359       if (catalog[c].measure[m].dbFlags & MEAS_BAD) {
     296      if (catalog[c].measureT[m].dbFlags & MEAS_BAD) {
    360297          Nbad++;
    361298          continue;
    … …  
    377314      }
    378315     
    379       n = catalog[c].measure[m].averef;
    380       Msys = PhotSys (&catalog[c].measure[m], &catalog[c].average[n], &catalog[c].secfilt[n*Nsecfilt]);
     316      n = catalog[c].measureT[m].averef;
     317      Msys = PhotSysTiny (&catalog[c].measureT[m], &catalog[c].averageT[n], &catalog[c].secfilt[n*Nsecfilt]);
    381318      if (isnan(Msys)) {
    382319        Nsys++;
    … …  
    384321      }
    385322      list[N] = Msys - Mrel - Mmos - Mgrid;
    386       dlist[N] = MAX (catalog[c].measure[m].dM, MIN_ERROR);
    387       if (catalog[c].measure[m].dM < IMFIT_SYS_SIGMA_LIM) {
     323      dlist[N] = MAX (catalog[c].measureT[m].dM, MIN_ERROR);
     324      if (catalog[c].measureT[m].dM < IMFIT_SYS_SIGMA_LIM) {
    388325        McalBright += list[N];
    389326        McalBright2 += SQ(list[N]);
    … …  
    650587  return (stats);
    651588}
    652 

  • branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/MosaicOps.c

    r31272 r31273  
    272272    for (j = 0; j < catalog[i].Nmeasure; j++) {
    273273      if (TimeSelect) {
    274         if (catalog[i].measure[j].t < TSTART) continue;
    275         if (catalog[i].measure[j].t > TSTOP) continue;
    276       }
    277       ecode = GetPhotcodeEquivCodebyCode (catalog[i].measure[j].photcode);
     274        if (catalog[i].measureT[j].t < TSTART) continue;
     275        if (catalog[i].measureT[j].t > TSTOP) continue;
     276      }
     277      ecode = GetPhotcodeEquivCodebyCode (catalog[i].measureT[j].photcode);
    278278      found = FALSE;
    279279      for (Ns = 0; !found && (Ns < Nphotcodes); Ns++) {
    … …  
    289289}
    290290
    291 # if 1
    292 // # if USE_IMAGE_ID
    293 
    294291void matchMosaics (Catalog *catalog, off_t meas, int cat) {
    295292
    296293  off_t idx, ID, mosID;
    297   Measure *measure;
    298 
    299   measure = &catalog[cat].measure[meas];
     294  MeasureTiny *measure;
     295
     296  measure = &catalog[cat].measureT[meas];
    300297
    301298  ID = measure[0].imageID;
    … …  
    336333}
    337334
    338 # else // USE_IMAGE_ID
    339 
    340 void matchMosaics (Catalog *catalog, off_t meas, int cat) {
    341 
    342   int i;
    343 
    344   for (i = 0; i < Nmosaic; i++) {
    345     if (catalog[cat].measure[meas].t < mosaic[i].start) continue;
    346     if (catalog[cat].measure[meas].t > mosaic[i].stop) continue;
    347    
    348 # ifdef GRID_V1
    349     if (USE_GRID) {
    350       ave = catalog[cat].measure[meas].averef;
    351       ra  = catalog[cat].average[ave].R - catalog[cat].measure[meas].dR / 3600.0;
    352       dec = catalog[cat].average[ave].D - catalog[cat].measure[meas].dD / 3600.0;
    353 
    354       /* X,Y always positive-definite in range 0,0 - dX, dY */
    355       RD_to_XY (&X, &Y, ra, dec, &mosaic[i].coords);
    356       setGridMeasure (meas, cat, X, Y);
    357     }
    358 # endif // GRID_V1
    359 
    360     bin[cat][meas] = i;
    361 
    362     clist[i][Nlist[i]] = cat;
    363     mlist[i][Nlist[i]] = meas;
    364     Nlist[i] ++;
    365    
    366     if (Nlist[i] == NLIST[i]) {
    367       NLIST[i] += 100;
    368       REALLOCATE (clist[i], int,   NLIST[i]);
    369       REALLOCATE (mlist[i], off_t, NLIST[i]);
    370     }   
    371     return;
    372   }
    373   // if we miss a measurement, it might be a REF type of value
    374   if (VERBOSE2) fprintf (stderr, "missed measurement "OFF_T_FMT", %d\n", meas, cat);
    375   return;
    376 }
    377 # endif // USE_IMAGE_ID
    378 
    379335float getMmos (off_t meas, int cat) {
    380336
    … …  
    437393      c = clist[i][j];
    438394     
    439       if (catalog[c].measure[m].dbFlags & MEAS_BAD) {
     395      if (catalog[c].measureT[m].dbFlags & MEAS_BAD) {
    440396          Nbad ++;
    441397          continue;
    … …  
    457413      }
    458414     
    459       n = catalog[c].measure[m].averef;
    460       Msys = PhotSys (&catalog[c].measure[m], &catalog[c].average[n], &catalog[c].secfilt[n*Nsecfilt]);
     415      n = catalog[c].measureT[m].averef;
     416      Msys = PhotSysTiny (&catalog[c].measureT[m], &catalog[c].averageT[n], &catalog[c].secfilt[n*Nsecfilt]);
    461417      if (isnan(Msys)) {
    462418        Nsys++;
    … …  
    464420      }
    465421      list[N]  = Msys - Mrel - Mcal - Mgrid;
    466       dlist[N] = MAX (catalog[c].measure[m].dM, MIN_ERROR);
     422      dlist[N] = MAX (catalog[c].measureT[m].dM, MIN_ERROR);
    467423      Mlist[N] = Msys;
    468424      dMlist[N] = list[N];
    … …  
    580536      memset (starNcount[Nstars], 0, Nphotcodes*sizeof(int));
    581537     
    582       m = catalog[i].average[j].measureOffset;
     538      m = catalog[i].averageT[j].measureOffset;
    583539
    584540      // determine the mosaic for each measurement
    585       for (k = 0; k < catalog[i].average[j].Nmeasure; k++, m++) {
     541      for (k = 0; k < catalog[i].averageT[j].Nmeasure; k++, m++) {
    586542
    587543        // skip unused measurements
    588         if (catalog[i].measure[m].dbFlags & MEAS_BAD) continue;
     544        if (catalog[i].measureT[m].dbFlags & MEAS_BAD) continue;
    589545
    590546        // skip unused measurements
    591547        int Nsec;
    592         int ecode = GetPhotcodeEquivCodebyCode (catalog[i].measure[m].photcode);
     548        int ecode = GetPhotcodeEquivCodebyCode (catalog[i].measureT[m].photcode);
    593549        found = FALSE;
    594550        for (Ns = 0; !found && (Ns < Nphotcodes); Ns ++) {
    … …  
    914870      c = clist[i][j];
    915871     
    916       if (catalog[c].measure[m].dbFlags & (ID_MEAS_AREA | ID_MEAS_NOCAL)) continue;
    917 
    918       ave = catalog[c].measure[m].averef;
    919       xlist[N] = catalog[c].average[ave].R - catalog[c].measure[m].dR / 3600.0;
    920       ylist[N] = catalog[c].average[ave].D - catalog[c].measure[m].dD / 3600.0;
     872      if (catalog[c].measureT[m].dbFlags & (ID_MEAS_AREA | ID_MEAS_NOCAL)) continue;
     873
     874      ave = catalog[c].measureT[m].averef;
     875      xlist[N] = catalog[c].averageT[ave].R - catalog[c].measureT[m].dR / 3600.0;
     876      ylist[N] = catalog[c].averageT[ave].D - catalog[c].measureT[m].dD / 3600.0;
    921877      N++;
    922878    }

  • branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/StarOps.c

    r31272 r31273  
    1818  for (i = 0; i < Ncatalog; i++) {
    1919    for (j = 0; j < catalog[i].Naverage; j++) {
    20       Nmax = MAX (Nmax, catalog[i].average[j].Nmeasure);
     20      Nmax = MAX (Nmax, catalog[i].averageT[j].Nmeasure);
    2121    }
    2222  }
    … …  
    3333  float value;
    3434
    35   ave = catalog[cat].measure[meas].averef;
    36   photcode = catalog[cat].measure[meas].photcode;
     35  ave = catalog[cat].measureT[meas].averef;
     36  photcode = catalog[cat].measureT[meas].photcode;
    3737
    3838  int ecode = GetPhotcodeEquivCodebyCode (photcode);
    … …  
    7070        // skip bad stars
    7171        if (catalog[i].secfilt[Nsecfilt*j+Nsec].flags & STAR_BAD) continue;
    72         m = catalog[i].average[j].measureOffset;
     72        m = catalog[i].averageT[j].measureOffset;
    7373
    7474        N = 0;
    75         for (k = 0; k < catalog[i].average[j].Nmeasure; k++, m++) {
     75        for (k = 0; k < catalog[i].averageT[j].Nmeasure; k++, m++) {
    7676          // skip measurements that do not match the current photcode
    77           int ecode = GetPhotcodeEquivCodebyCode (catalog[i].measure[m].photcode);
     77          int ecode = GetPhotcodeEquivCodebyCode (catalog[i].measureT[m].photcode);
    7878          if (ecode != thisCode) { continue; }
    7979
    80           if (catalog[i].measure[m].dbFlags & MEAS_BAD) {
     80          if (catalog[i].measureT[m].dbFlags & MEAS_BAD) {
    8181            Nbad ++;
    8282            continue;
    … …  
    105105          }
    106106
    107           Msys = PhotSys (&catalog[i].measure[m], &catalog[i].average[j], &catalog[i].secfilt[j*Nsecfilt]);
     107          Msys = PhotSysTiny (&catalog[i].measureT[m], &catalog[i].averageT[j], &catalog[i].secfilt[j*Nsecfilt]);
    108108          if (isnan(Msys)) {
    109109            Nsys++;
    … …  
    111111          }
    112112          list[N] = Msys - Mcal - Mmos - Mgrid;
    113           dlist[N] = MAX (catalog[i].measure[m].dM, MIN_ERROR);
     113          dlist[N] = MAX (catalog[i].measureT[m].dM, MIN_ERROR);
    114114
    115115          // tie down reference photometry if the -refcode (code) option is selected
    116116          if (refPhotcode) {
    117             if (GetPhotcodeEquivCodebyCode(catalog[i].measure[m].photcode) == refPhotcode[0].equiv) {
     117            if (GetPhotcodeEquivCodebyCode(catalog[i].measureT[m].photcode) == refPhotcode[0].equiv) {
    118118              // increase the weight by a factor of 100:
    119               dlist[N] = 0.01*catalog[i].measure[m].dM;
     119              dlist[N] = 0.01*catalog[i].measureT[m].dM;
    120120            }
    121121          }
    … …  
    159159  for (i = 0; i < Ncatalog; i++) {
    160160    for (j = 0; j < catalog[i].Naverage; j++) {
    161       Nmax = MAX (Nmax, catalog[i].average[j].Nmeasure);
     161      Nmax = MAX (Nmax, catalog[i].averageT[j].Nmeasure);
    162162    }
    163163  }
    … …  
    176176
    177177        N = 0;
    178         m = catalog[i].average[j].measureOffset;
    179         for (k = 0; k < catalog[i].average[j].Nmeasure; k++, m++) {
     178        m = catalog[i].averageT[j].measureOffset;
     179        for (k = 0; k < catalog[i].averageT[j].Nmeasure; k++, m++) {
    180180          // skip measurements that do not match the current photcode
    181           ecode = GetPhotcodeEquivCodebyCode (catalog[i].measure[m].photcode);
     181          ecode = GetPhotcodeEquivCodebyCode (catalog[i].measureT[m].photcode);
    182182          if (ecode != thisCode) { continue; }
    183183
    184           if (catalog[i].measure[m].dbFlags & MEAS_BAD) continue;
     184          if (catalog[i].measureT[m].dbFlags & MEAS_BAD) continue;
    185185
    186186          // XXX allow REF stars (no Image Entry) to be included in the calculation this
    … …  
    198198          }
    199199
    200           Msys = PhotSys (&catalog[i].measure[m], &catalog[i].average[j], &catalog[i].secfilt[j*Nsecfilt]);
     200          Msys = PhotSysTiny (&catalog[i].measureT[m], &catalog[i].averageT[j], &catalog[i].secfilt[j*Nsecfilt]);
    201201          list[N] = Msys - Mcal - Mmos - Mgrid;
    202           dlist[N] = MAX (catalog[i].measure[m].dM, MIN_ERROR);
     202          dlist[N] = MAX (catalog[i].measureT[m].dM, MIN_ERROR);
    203203          N++;
    204204        }
    … …  
    221221}
    222222
    223 // for each average object, set the average mags based on existing equiv photometry
     223// For each average object, set the average mags based on existing equiv photometry.
     224// NOTE: this function operates on the real Measure & Average structures, not the
     225// MeasureTiny & AverageTiny structures
    224226int setMave (Catalog *catalog, int Ncatalog) {
    225227
    … …  
    363365    for (j = 0; j < catalog[i].Naverage; j++) {
    364366
    365       m = catalog[i].average[j].measureOffset;
    366       for (k = 0; k < catalog[i].average[j].Nmeasure; k++, m++) {
    367         if (catalog[i].measure[m].dbFlags & MEAS_BAD) continue;
     367      m = catalog[i].averageT[j].measureOffset;
     368      for (k = 0; k < catalog[i].averageT[j].Nmeasure; k++, m++) {
     369        if (catalog[i].measureT[m].dbFlags & MEAS_BAD) continue;
    368370        Mcal  = getMcal  (m, i);
    369371        if (isnan(Mcal)) continue;
    … …  
    372374        Mgrid = getMgrid (m, i);
    373375        if (isnan(Mgrid)) continue;
     376
     377        // set the output calibration
    374378        catalog[i].measure[m].Mcal = Mcal + Mmos + Mgrid;
    375379      }
    … …  
    469473  for (i = 0; i < Ncatalog; i++) {
    470474    for (j = 0; j < catalog[i].Naverage; j++) {
    471       Nmax = MAX (Nmax, catalog[i].average[j].Nmeasure);
     475      Nmax = MAX (Nmax, catalog[i].averageT[j].Nmeasure);
    472476    }
    473477  }
    … …  
    498502
    499503        /* accumulate list of valid measurements */
    500         m = catalog[i].average[j].measureOffset;
     504        m = catalog[i].averageT[j].measureOffset;
    501505        N = 0;
    502         for (k = 0; k < catalog[i].average[j].Nmeasure; k++, m++) {
     506        for (k = 0; k < catalog[i].averageT[j].Nmeasure; k++, m++) {
    503507          // skip measurements that do not match the current photcode
    504           int ecode = GetPhotcodeEquivCodebyCode (catalog[i].measure[m].photcode);
     508          int ecode = GetPhotcodeEquivCodebyCode (catalog[i].measureT[m].photcode);
    505509          if (ecode != thisCode) { continue; }
    506510
    507           /* if (catalog[i].measure[m].dbFlags & MEAS_BAD) continue; */
     511          /* if (catalog[i].measureT[m].dbFlags & MEAS_BAD) continue; */
    508512          Mcal  = getMcal  (m, i);
    509513          if (isnan(Mcal)) { Ncal ++; continue; }
    … …  
    513517          if (isnan(Mgrid)) { Ngrid ++; continue; }
    514518
    515           Msys = PhotSys (&catalog[i].measure[m], &catalog[i].average[j], &catalog[i].secfilt[j*Nsecfilt]);
     519          Msys = PhotSysTiny (&catalog[i].measureT[m], &catalog[i].averageT[j], &catalog[i].secfilt[j*Nsecfilt]);
    516520          list[N] = Msys - Mcal - Mmos - Mgrid;
    517           dlist[N] = MAX (catalog[i].measure[m].dM, MIN_ERROR);
     521          dlist[N] = MAX (catalog[i].measureT[m].dM, MIN_ERROR);
    518522          N++;
    519523        }
    … …  
    542546        image_bad = IMAGE_BAD;
    543547        IMAGE_BAD = ID_IMAGE_PHOTOM_NOCAL;
    544         m = catalog[i].average[j].measureOffset;
     548        m = catalog[i].averageT[j].measureOffset;
    545549        N = 0;
    546         for (k = 0; k < catalog[i].average[j].Nmeasure; k++, m++) {
     550        for (k = 0; k < catalog[i].averageT[j].Nmeasure; k++, m++) {
    547551          // skip measurements that do not match the current photcode
    548           int ecode = GetPhotcodeEquivCodebyCode (catalog[i].measure[m].photcode);
     552          int ecode = GetPhotcodeEquivCodebyCode (catalog[i].measureT[m].photcode);
    549553          if (ecode != thisCode) { continue; }
    550554
    551           /* if (catalog[i].measure[m].dbFlags & MEAS_BAD) continue; */
     555          /* if (catalog[i].measureT[m].dbFlags & MEAS_BAD) continue; */
    552556          Mcal  = getMcal  (m, i);
    553557          if (isnan(Mcal)) continue;
    … …  
    557561          if (isnan(Mgrid)) continue;
    558562
    559           Msys = PhotSys (&catalog[i].measure[m], &catalog[i].average[j], &catalog[i].secfilt[j*Nsecfilt]);
     563          Msys = PhotSysTiny (&catalog[i].measureT[m], &catalog[i].averageT[j], &catalog[i].secfilt[j*Nsecfilt]);
    560564          list[N] = Msys - Mcal - Mmos - Mgrid;
    561           dlist[N] = MAX (catalog[i].measure[m].dM, MIN_ERROR);
     565          dlist[N] = MAX (catalog[i].measureT[m].dM, MIN_ERROR);
    562566          ilist[N] = m;
    563567          N++;
    … …  
    569573        for (k = 0; k < N; k++) {
    570574          if (fabs (list[k] - stats.median) > NSIGMA_REJECT*stats.sigma) {
    571             catalog[i].measure[ilist[k]].dbFlags |= ID_MEAS_POOR_PHOTOM;
     575            catalog[i].measureT[ilist[k]].dbFlags |= ID_MEAS_POOR_PHOTOM;
     576            if (final) {
     577              // for the final pass, we have a duplicate set of values in measure and measureT
     578              catalog[i].measure[ilist[k]].dbFlags |= ID_MEAS_POOR_PHOTOM;
     579            }
    572580            Ndel ++;
    573581          }
    … …  
    609617      /* calculate the average value for a single star */
    610618      if (catalog[i].secfilt[Nsecfilt*j+Nsec].flags & STAR_BAD) { N1++; continue;  }
    611       m = catalog[i].average[j].measureOffset;
     619      m = catalog[i].averageT[j].measureOffset;
    612620
    613621      N = 0;
    614       for (k = 0; k < catalog[i].average[j].Nmeasure; k++, m++) {
    615         int ecode = GetPhotcodeEquivCodebyCode (catalog[i].measure[m].photcode);
     622      for (k = 0; k < catalog[i].averageT[j].Nmeasure; k++, m++) {
     623        int ecode = GetPhotcodeEquivCodebyCode (catalog[i].measureT[m].photcode);
    616624        if (ecode != seccode) { N0++; continue;}
    617625        Mcal = getMcal  (m, i);
    … …  
    816824  for (i = 0; i < Ncatalog; i++) {
    817825    for (j = 0; j < catalog[i].Naverage; j++) {
    818       xlist[N] = catalog[i].average[j].R;
    819       ylist[N] = catalog[i].average[j].D;
     826      xlist[N] = catalog[i].averageT[j].R;
     827      ylist[N] = catalog[i].averageT[j].D;
    820828      N++;
    821829    }

  • branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/bcatalog.c

    r31272 r31273  
    22
    33extern double drand48();
     4
     5int CopyAverageTiny (AverageTiny *averageT, Average *average) {
     6
     7  averageT[0].R     = average[0].R;
     8  averageT[0].D     = average[0].D;
     9  averageT[0].flags = average[0].flags;
     10  averageT[0].Nmeasure      = average[0].Nmeasure;
     11  averageT[0].measureOffset = average[0].measureOffset;
     12
     13  // make Nmeasure & measureOffset optional?
     14
     15  return (TRUE);
     16}
     17
     18int CopyMeasureTiny (MeasureTiny *measureT, Measure *measure) {
     19
     20  measureT[0].dR       = measure[0].dR;
     21  measureT[0].dD       = measure[0].dD;
     22  measureT[0].M        = measure[0].M;
     23  measureT[0].Mcal     = measure[0].Mcal;
     24  measureT[0].dM       = measure[0].dM;
     25  measureT[0].airmass  = measure[0].airmass;
     26  measureT[0].Xccd     = measure[0].Xccd;
     27  measureT[0].Yccd     = measure[0].Yccd;
     28  measureT[0].t        = measure[0].t;
     29  measureT[0].dt       = measure[0].dt;
     30  measureT[0].averef   = measure[0].averef;
     31  measureT[0].imageID  = measure[0].imageID;
     32  measureT[0].dbFlags  = measure[0].dbFlags;
     33  measureT[0].photcode = measure[0].photcode;
     34
     35  return (TRUE);
     36}
    437
    538int bcatalog (Catalog *subcatalog, Catalog *catalog) {
    … …  
    1851  NAVERAGE = 50;
    1952  NMEASURE = 1000;
    20   ALLOCATE (subcatalog[0].average, Average, NAVERAGE);
     53  ALLOCATE (subcatalog[0].measureT, MeasureTiny, NMEASURE);
     54  ALLOCATE (subcatalog[0].averageT, AverageTiny, NAVERAGE);
    2155  ALLOCATE (subcatalog[0].secfilt, SecFilt, NAVERAGE*Nsecfilt);
    22   ALLOCATE (subcatalog[0].measure, Measure, NMEASURE);
    2356  Nmeasure = Naverage = 0;
    2457
    … …  
    3063
    3164    /* start with all stars good */
    32     subcatalog[0].average[Naverage] = catalog[0].average[i];
    33     subcatalog[0].average[Naverage].measureOffset = Nmeasure;
     65    CopyAverageTiny (&subcatalog[0].averageT[Naverage], &catalog[0].average[i]);
     66    subcatalog[0].averageT[Naverage].measureOffset = Nmeasure;
     67
    3468    for (j = 0; j < Nsecfilt; j++) {
    3569      subcatalog[0].secfilt[Nsecfilt*Naverage+j] = catalog[0].secfilt[Nsecfilt*i+j];
    … …  
    103137      }
    104138
    105       subcatalog[0].measure[Nmeasure].dbFlags &= ~ID_MEAS_SKIP_PHOTOM;
    106       subcatalog[0].measure[Nmeasure]        = catalog[0].measure[offset];
    107       subcatalog[0].measure[Nmeasure].averef = Naverage;
     139      CopyMeasureTiny (&subcatalog[0].measureT[Nmeasure], &catalog[0].measure[offset]);
     140      subcatalog[0].measureT[Nmeasure].dbFlags &= ~ID_MEAS_SKIP_PHOTOM;
     141      subcatalog[0].measureT[Nmeasure].averef = Naverage;
    108142      if (RESET) {
    109         subcatalog[0].measure[Nmeasure].Mcal = 0;
    110         subcatalog[0].measure[Nmeasure].dbFlags &= 0xff00;
    111         subcatalog[0].measure[Nmeasure].dbFlags &= ~ID_MEAS_POOR_PHOTOM;
    112         subcatalog[0].measure[Nmeasure].dbFlags &= ~ID_MEAS_AREA;
    113         subcatalog[0].measure[Nmeasure].dbFlags &= ~ID_MEAS_NOCAL;
     143        subcatalog[0].measureT[Nmeasure].Mcal = 0;
     144        subcatalog[0].measureT[Nmeasure].dbFlags &= 0xff00;
     145        subcatalog[0].measureT[Nmeasure].dbFlags &= ~ID_MEAS_POOR_PHOTOM;
     146        subcatalog[0].measureT[Nmeasure].dbFlags &= ~ID_MEAS_AREA;
     147        subcatalog[0].measureT[Nmeasure].dbFlags &= ~ID_MEAS_NOCAL;
    114148      }
    115149      Nmeasure ++;
    … …  
    117151      if (Nmeasure == NMEASURE) {
    118152        NMEASURE += 1000;
    119         REALLOCATE (subcatalog[0].measure, Measure, NMEASURE);
     153        REALLOCATE (subcatalog[0].measureT, MeasureTiny, NMEASURE);
    120154      }
    121155    }
    … …  
    134168      continue;
    135169    }
    136     subcatalog[0].average[Naverage].Nmeasure = Nm;
     170    subcatalog[0].averageT[Naverage].Nmeasure = Nm;
    137171    Naverage ++;
    138172    if (Naverage == NAVERAGE) {
    139173      NAVERAGE += 50;
    140       REALLOCATE (subcatalog[0].average, Average, NAVERAGE);
    141       REALLOCATE (subcatalog[0].secfilt, SecFilt, NAVERAGE*Nsecfilt);
    142     }
    143   }
    144   REALLOCATE (subcatalog[0].average, Average, MAX (Naverage, 1));
    145   REALLOCATE (subcatalog[0].measure, Measure, MAX (Nmeasure, 1));
     174      REALLOCATE (subcatalog[0].averageT, AverageTiny, NAVERAGE);
     175      REALLOCATE (subcatalog[0].secfilt,  SecFilt,    NAVERAGE*Nsecfilt);
     176    }
     177  }
     178  REALLOCATE (subcatalog[0].averageT, AverageTiny, MAX (Naverage, 1));
     179  REALLOCATE (subcatalog[0].measureT, MeasureTiny, MAX (Nmeasure, 1));
    146180  REALLOCATE (subcatalog[0].secfilt, SecFilt, Nsecfilt*MAX (Naverage, 1));
    147181  subcatalog[0].Naverage = Naverage;
    … …  
    211245  for (i = 0; i < Nmax; i++) {
    212246    ave = index[i];
    213     NMEASURE += subcatalog[0].average[ave].Nmeasure;
     247    NMEASURE += subcatalog[0].averageT[ave].Nmeasure;
    214248  }
    215249
    216250  // allocate the output data
    217   ALLOCATE (tmpcatalog.average, Average, Nmax);
    218   ALLOCATE (tmpcatalog.secfilt, SecFilt, Nmax * Nsecfilt);
    219   ALLOCATE (tmpcatalog.measure, Measure, NMEASURE);
     251  ALLOCATE (tmpcatalog.averageT, AverageTiny, Nmax);
     252  ALLOCATE (tmpcatalog.measureT, MeasureTiny, NMEASURE);
     253  ALLOCATE (tmpcatalog.secfilt,  SecFilt, Nmax * Nsecfilt);
    220254
    221255  off_t Nmeasure = 0;
    … …  
    224258  for (i = 0; i < Nmax; i++) {
    225259    ave = index[i];
    226     tmpcatalog.average[i] = subcatalog[0].average[ave];
    227     tmpcatalog.average[i].measureOffset = Nmeasure;
    228     for (j = 0; j < tmpcatalog.average[i].Nmeasure; j++) {
    229       off_t offset = subcatalog[0].average[ave].measureOffset + j;
    230       tmpcatalog.measure[Nmeasure] = subcatalog[0].measure[offset];
    231       tmpcatalog.measure[Nmeasure].averef = i;
     260    tmpcatalog.averageT[i] = subcatalog[0].averageT[ave];
     261    tmpcatalog.averageT[i].measureOffset = Nmeasure;
     262    for (j = 0; j < tmpcatalog.averageT[i].Nmeasure; j++) {
     263      off_t offset = subcatalog[0].averageT[ave].measureOffset + j;
     264      tmpcatalog.measureT[Nmeasure] = subcatalog[0].measureT[offset];
     265      tmpcatalog.measureT[Nmeasure].averef = i;
    232266      Nmeasure ++;
    233267    }
    … …  
    239273  }
    240274
    241   free (subcatalog[0].average);
    242   free (subcatalog[0].measure);
     275  free (subcatalog[0].averageT);
     276  free (subcatalog[0].measureT);
    243277  free (subcatalog[0].secfilt);
    244278
    245   subcatalog[0].average = tmpcatalog.average;
    246   subcatalog[0].measure = tmpcatalog.measure;
     279  subcatalog[0].averageT = tmpcatalog.averageT;
     280  subcatalog[0].measureT = tmpcatalog.measureT;
    247281  subcatalog[0].secfilt = tmpcatalog.secfilt;
    248282  subcatalog[0].Naverage = Nmax;
    … …  
    253287  return (TRUE);
    254288}
     289
     290// for the cases where we are not using a subset of the data, we still need to have a copy of these fields
     291int populate_tiny_values (Catalog *catalog) {
     292
     293  off_t i;
     294
     295  ALLOCATE (catalog[0].measureT, MeasureTiny, catalog[0].Nmeasure);
     296  ALLOCATE (catalog[0].averageT, AverageTiny, catalog[0].Naverage);
     297
     298  for (i = 0; i < catalog[0].Naverage; i++) {
     299    CopyAverageTiny (&catalog[0].averageT[i], &catalog[0].average[i]);
     300  }
     301
     302  for (i = 0; i < catalog[0].Nmeasure; i++) {
     303    CopyMeasureTiny (&catalog[0].measureT[i], &catalog[0].measure[i]);
     304  }
     305
     306  return (TRUE);
     307}
     308
     309int free_tiny_values (Catalog *catalog) {
     310
     311  free (catalog[0].averageT);
     312  free (catalog[0].measureT);
     313  return (TRUE);
     314}
     315

  • branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/plot_scatter.c

    r31272 r31273  
    3636            for (k = 0; k < catalog[i].average[j].Nmeasure; k++, m++) {
    3737                // skip measurements that do not match the current photcode
    38                 int ecode = GetPhotcodeEquivCodebyCode (catalog[i].measure[m].photcode);
     38                int ecode = GetPhotcodeEquivCodebyCode (catalog[i].measureT[m].photcode);
    3939                if (ecode != thisCode) { continue; }
    4040
    41                 if (catalog[i].measure[m].dbFlags & MEAS_BAD) continue;
     41                if (catalog[i].measureT[m].dbFlags & MEAS_BAD) continue;
    4242                Mcal = getMcal  (m, i);
    4343                if (isnan(Mcal)) continue;
    … …  
    4949                Mrel = catalog[i].secfilt[Nsecfilt*j+Nsec].M;
    5050                xlist[N] = Mrel;
    51                 ylist[N] = PhotSys  (&catalog[i].measure[m], &catalog[i].average[j], &catalog[i].secfilt[j*Nsecfilt]) - Mcal - Mmos - Mgrid - Mrel;
    52                 ilist[N] = PhotInst (&catalog[i].measure[m]);
     51                ylist[N] = PhotSysTiny (&catalog[i].measureT[m], &catalog[i].averageT[j], &catalog[i].secfilt[j*Nsecfilt]) - Mcal - Mmos - Mgrid - Mrel;
     52                ilist[N] = PhotInstTiny (&catalog[i].measureT[m]);
    5353                N++;
    5454            }

  • branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/reload_catalogs.c

    r31260 r31273  
    4949    TIMESTAMP(time2);
    5050
     51    populate_tiny_values(&catalog);
     52
    5153    initImageBins  (&catalog, 1);
    5254    initMosaicBins (&catalog, 1);
    53     initGridBins   (&catalog, 1);
     55    initGridBins   (&catalog, 1); 
    5456    TIMESTAMP(time3);
    5557
    56     findImages (&catalog, 1);
    57     findMosaics (&catalog, 1);
     58    findImages (&catalog, 1);  // FX
     59    findMosaics (&catalog, 1); //
    5860    TIMESTAMP(time4);
    5961
    … …  
    6365    dvo_catalog_save (&catalog, VERBOSE);
    6466    dvo_catalog_unlock (&catalog);
     67
     68    free_tiny_values(&catalog);
    6569    dvo_catalog_free (&catalog);
    6670    TIMESTAMP(time6);

  • branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/relphot.c

    r31272 r31273  
    9494      int star_toofew;
    9595
    96 # if (USE_DIRECT)
    97       // until we finish the grid analysis, do not reject stars out-of-hand based on ID_STAR_FEW
    98       // XXX this is kind of poor: need to have a better distinctions about STAR_BAD in setMrel vs getMrel
    99       star_toofew = STAR_TOOFEW;
    100       STAR_TOOFEW = 0;
    101       STAR_BAD  = ID_STAR_POOR;
    102 
    103       showGridCount ();
    104       setMgridDirect (catalog, Ncatalog);
    105 
    106       STAR_BAD  = ID_STAR_POOR | ID_STAR_FEW;
    107       STAR_TOOFEW = star_toofew;
    108 
    109       dump_grid ();
    110       exit (0);
    111 
    112 # else
    113 
    11496      // until we finish the grid analysis, do not reject stars out-of-hand based on ID_STAR_FEW
    11597      // XXX this is kind of poor: need to have a better distinctions about STAR_BAD in setMrel vs getMrel
    … …  
    124106      STAR_BAD  = ID_STAR_POOR | ID_STAR_FEW;
    125107      STAR_TOOFEW = star_toofew;
    126 # endif
    127108  }
    128109
    … …  
    175156
    176157  /* at this point, we have correct cal coeffs in the image/mosaic structures */
    177   for (i = 0; i < Ncatalog; i++) dvo_catalog_free (&catalog[i]);
     158  for (i = 0; i < Ncatalog; i++) {
     159    free_tiny_values (&catalog[i]);
     160    dvo_catalog_free (&catalog[i]);
     161  }
    178162  freeImageBins (Ncatalog);
    179163  freeMosaicBins (Ncatalog);

  • branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/setExclusions.c

    r31261 r31273  
    11# include "relphot.h"
     2
     3// this function sets the NOCAL and AREA dbFlags bits for the MeasureTiny elements these
     4// are used elsewhere (StarOps.c, ImageOps.c, MosaicOps.c, GridOps.c, etc) to skip bad
     5// measurements.  The only exception is 'setMave' which is called by 'relphot_objects',
     6// and uses the bits read from disk as the test
    27
    38int setExclusions (Catalog *catalog, int Ncatalog) {
    … …  
    1116  for (i = 0; i < Ncatalog; i++) {
    1217    for (j = 0; j < catalog[i].Naverage; j++) {
    13       m = catalog[i].average[j].measureOffset;
    14       for (k = 0; k < catalog[i].average[j].Nmeasure; k++, m++) {
     18      m = catalog[i].averageT[j].measureOffset;
     19      for (k = 0; k < catalog[i].averageT[j].Nmeasure; k++, m++) {
    1520
    1621        /* select measurements by photcode */
    17         ecode = GetPhotcodeEquivCodebyCode (catalog[i].measure[m].photcode);
     22        ecode = GetPhotcodeEquivCodebyCode (catalog[i].measureT[m].photcode);
    1823        found = FALSE;
    1924        for (Ns = 0; !found && (Ns < Nphotcodes); Ns++) {
    … …  
    2429        /* select measurements by time */
    2530        if (TimeSelect) {
    26           if (catalog[i].measure[m].t < TSTART) goto mark_nocal;
    27           if (catalog[i].measure[m].t > TSTOP) goto mark_nocal;
     31          if (catalog[i].measureT[m].t < TSTART) goto mark_nocal;
     32          if (catalog[i].measureT[m].t > TSTOP) goto mark_nocal;
    2833        }
    2934
    3035        /* select measurements by mag limit */
    3136        if (AreaSelect) {
    32           r = catalog[i].average[j].R + catalog[i].measure[m].dR / 3600.0;
    33           d = catalog[i].average[j].D + catalog[i].measure[m].dD / 3600.0;
     37          r = catalog[i].averageT[j].R + catalog[i].measureT[m].dR / 3600.0;
     38          d = catalog[i].averageT[j].D + catalog[i].measureT[m].dD / 3600.0;
    3439          if ((coords = getCoords (m, i)) == NULL) goto markbad;
    3540          RD_to_XY (&x, &y, r, d, coords);
    … …  
    4348
    4449      markbad:
    45         catalog[i].measure[m].dbFlags |= ID_MEAS_AREA;
     50        catalog[i].measureT[m].dbFlags |= ID_MEAS_AREA;
    4651        Narea ++;
    4752        continue;
    4853       
    4954      mark_nocal:
    50         catalog[i].measure[m].dbFlags |= ID_MEAS_NOCAL;
     55        catalog[i].measureT[m].dbFlags |= ID_MEAS_NOCAL;
    5156        Nnocal ++;
    5257        continue;

  • branches/eam_branches/ipp-20110404/Ohana/src/relphot/src/setMrelFinal.c

    r31272 r31273  
    11# include "relphot.h"
     2
     3// we've just reloaded the data from disk; we now need to apply the Image/Mosaic/Grid
     4// calibrations determined by the rest of the program.  We also need to set the final
     5// output dbFlags values
    26
    37void setMrelFinal (Catalog *catalog) {
    … …  
    4650  }
    4751
     52  // this sets flags in the measureT element, not the measure element
    4853  setExclusions (catalog, 1);  /* mark by area */
    4954
    … …  
    6469  int Nsecfilt = GetPhotcodeNsecfilt ();
    6570
    66   /* clear ID_STAR_POOR, ID_STAR_FEW, ID_MEAS_NOCAL values before writing ??? */
     71  /* clear ID_STAR_POOR, ID_STAR_FEW values before writing ??? */
     72  /* ID_MEAS_NOCAL is an internal bit, so it should be cleared */
    6773  for (i = 0; i < catalog[0].Naverage; i++) {
    6874    catalog[0].average[i].flags &= ~ID_STAR_FEW;
    … …  
    105111
    106112      /* clear SKIP for all measures at first */
    107       catalog[0].measure[m].dbFlags &= ~ID_MEAS_SKIP_PHOTOM;
     113      catalog[0].measureT[m].dbFlags &= ~ID_MEAS_SKIP_PHOTOM;
    108114
    109115      /** never use these measurements (wrong photcode, bad time range) */
    … …  
    153159
    154160    skip:
    155       catalog[0].measure[m].dbFlags |= ID_MEAS_SKIP_PHOTOM;
     161      catalog[0].measure [m].dbFlags |= ID_MEAS_SKIP_PHOTOM;
     162      catalog[0].measureT[m].dbFlags |= ID_MEAS_SKIP_PHOTOM;
    156163      Nskip ++;
    157164    }
Note: See TracChangeset for help on using the changeset viewer.