IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 39220


Ignore:
Timestamp:
Dec 3, 2015, 4:12:24 PM (11 years ago)
Author:
eugene
Message:

also merge galphot tables

Location:
branches/eam_branches/ipp-20151113/Ohana/src/dvomerge
Files:
7 edited

Legend:

Unmodified
Added
Removed

  • branches/eam_branches/ipp-20151113/Ohana/src/dvomerge/include/dvomerge.h

    r39182 r39220  
    164164StarPar   *sort_starpar           PROTO((Average *average, off_t Naverage, StarPar *starpar, off_t Nstarpar, off_t *next));
    165165
     166off_t     *build_galphot_links    PROTO((Average *average, off_t Naverage, GalPhot *galphot, off_t Ngalphot));
     167off_t     *init_galphot_links     PROTO((Average *average, off_t Naverage, GalPhot *galphot, off_t Ngalphot));
     168int        add_galp_link          PROTO((Average *average, off_t *next, off_t Ngalphot, off_t NGALPHOT));
     169GalPhot   *sort_galphot           PROTO((Average *average, off_t Naverage, GalPhot *galphot, off_t Ngalphot, off_t *next));
     170
    166171off_t     *init_missing_links     PROTO((Average *average, off_t Naverage, Missing *missing, off_t Nmissing));
    167172int        add_miss_link          PROTO((Average *average, off_t *next, off_t Nmissing));

  • branches/eam_branches/ipp-20151113/Ohana/src/dvomerge/src/LoadCatalog.c

    r38986 r39220  
    11# include "dvomerge.h"
     2# define SKIP_GALPHOT 0
    23
    34int LoadCatalog (Catalog *catalog, SkyRegion *region, char *filename, char *mode, int Nsecfilt) {
    … …  
    1011  // always load all of the data (if any exists)
    1112  // XXXX TEMP HACK : skip GALPHOT
    12   // XXXX catalog[0].catflags = DVO_LOAD_AVERAGE | DVO_LOAD_MEASURE | DVO_LOAD_MISSING | DVO_LOAD_SECFILT | DVO_LOAD_LENSING | DVO_LOAD_LENSOBJ | DVO_LOAD_STARPAR | DVO_LOAD_GALPHOT;
     13
     14# if (SKIP_GALPHOT) 
    1315  catalog[0].catflags = DVO_LOAD_AVERAGE | DVO_LOAD_MEASURE | DVO_LOAD_MISSING | DVO_LOAD_SECFILT | DVO_LOAD_LENSING | DVO_LOAD_LENSOBJ | DVO_LOAD_STARPAR;
     16# else
     17  catalog[0].catflags = DVO_LOAD_AVERAGE | DVO_LOAD_MEASURE | DVO_LOAD_MISSING | DVO_LOAD_SECFILT | DVO_LOAD_LENSING | DVO_LOAD_LENSOBJ | DVO_LOAD_STARPAR | DVO_LOAD_GALPHOT;
     18# endif
    1419 
    1520  catalog[0].catformat = dvo_catalog_catformat (CATFORMAT);  // set the default catformat from config data

  • branches/eam_branches/ipp-20151113/Ohana/src/dvomerge/src/build_links.c

    r38986 r39220  
    2424*/
    2525
     26/*** Measure ****************************************************************************************/
     27
    2628/* build the initial links assuming the table is sorted,
    2729   not partial, and has a correct set of average[].measureOffset,Nmeasure values */
    … …  
    3436
    3537  ALLOCATE (next_meas, off_t, Nmeasure);
    36   for (i = 0; i < Naverage; i++, N++) {
     38  for (i = 0; i < Naverage; i++) {
     39    if (!average[i].Nmeasure) continue;
     40    off_t m = average[i].measureOffset;
     41    myAssert (measure[m].averef == i, "not sorted");
    3742    for (j = 0; j < average[i].Nmeasure - 1; j++, N++) {
     43      myAssert (measure[m+j+1].averef == i, "not sorted");
    3844      next_meas[N] = N + 1;
    3945      if (N >= Nmeasure) {
    … …  
    5056      abort ();
    5157    }
     58    N++;
    5259  }
    5360  return (next_meas);
    … …  
    122129  }
    123130
    124   if (m == -1) {
     131  // if Nmeasure is 0, m may have been mis-set; add to the end
     132  if ((average[0].Nmeasure == 0) || (m == -1)) {
    125133    average[0].measureOffset = Nmeasure;
    126134  } else {
    … …  
    136144Measure *sort_measure (Average *average, off_t Naverage, Measure *measure, off_t Nmeasure, off_t *next_meas) {
    137145
    138   off_t i, k, n, N;
     146  off_t i, k, n, np, N;
    139147  Measure *tmpmeasure;
    140148
    141149  /* fix order of Measure (memory intensive, but fast) */
     150  np = -1; // previous entry in case of errors
    142151  N = 0;
    143152  ALLOCATE (tmpmeasure, Measure, Nmeasure);
    144153  for (i = 0; i < Naverage; i++) {
     154    if (!average[i].Nmeasure) continue;
    145155    n = average[i].measureOffset;
    146156    average[i].measureOffset = N;
    147157    int myObjID = average[i].objID;
    148158    for (k = 0; k < average[i].Nmeasure; k++, N++) {
    149       if (n == -1) abort();
     159      if (n == -1) {
     160        fprintf (stderr, "entry after %d has a problem\n", (int) np);
     161        abort();
     162      }
    150163      tmpmeasure[N] = measure[n];
    151164      myAssert (measure[n].averef == i, "error in averef?");
    152165      myAssert (measure[n].objID == myObjID, "error in objID?");
    153166      tmpmeasure[N].averef = i;
     167      np = n;
    154168      n = next_meas[n];
    155169    }
    … …  
    159173}
    160174
    161 /*******************************************************************************************/
     175/*** Missing ****************************************************************************************/
    162176
    163177/* build the initial links assuming the table is sorted */
    … …  
    229243}
    230244
    231 /*******************************************************************************************/
     245/*** Lensing ****************************************************************************************/
    232246
    233247/* build the initial links assuming the table is sorted,
    … …  
    352366  Lensing *tmplensing;
    353367
    354   /*
    355   for (i = 0; i < Naverage; i++) {
    356     if (average[i].Nlensing != 4) {
    357       fprintf (stderr, "check %d %d %d\n", (int) i, (int) average[i].Nlensing, (int) average[i].lensingOffset);
    358     }
    359   }
    360   */
    361 
    362368  /* fix order of Lensing (memory intensive, but fast) */
    363369  np = -1;
    … …  
    368374    n = average[i].lensingOffset;
    369375    average[i].lensingOffset = N;
     376    int myObjID = average[i].objID;
    370377    for (k = 0; k < average[i].Nlensing; k++, N++) {
    371378      if (n == -1) {
    … …  
    374381      }
    375382      tmplensing[N] = lensing[n];
    376       if (lensing[n].averef != i) abort();
     383      myAssert (lensing[n].averef == i, "error in averef");
     384      myAssert (lensing[n].objID == myObjID, "error in objID?");
    377385      tmplensing[N].averef = i;
    378386      np = n;
    … …  
    384392}
    385393
    386 /******** StarPar ********************************************************************************/
     394/*** StarPar ********************************************************************************/
    387395
    388396/* build the initial links assuming the table is sorted,
    … …  
    398406  for (i = 0; i < Naverage; i++) {
    399407    if (!average[i].Nstarpar) continue;
     408    off_t m = average[i].starparOffset;
     409    myAssert (starpar[m].averef == i, "not sorted");
    400410    for (j = 0; j < average[i].Nstarpar - 1; j++, N++) {
     411      myAssert (starpar[m+j+1].averef == i, "not sorted");
    401412      next_star[N] = N + 1;
    402413      if (N >= Nstarpar) {
    … …  
    486497  }
    487498
    488   if (m == -1) {
     499  // if Nmeasure is 0, m may have been mis-set; add to the end
     500  if ((average[0].Nstarpar == 0) || (m == -1)) {
    489501    average[0].starparOffset = Nstarpar;
    490502  } else {
    … …  
    503515  StarPar *tmpstarpar;
    504516
    505   /*
    506   for (i = 0; i < Naverage; i++) {
    507     if (average[i].Nstarpar != 4) {
    508       fprintf (stderr, "check %d %d %d\n", (int) i, (int) average[i].Nstarpar, (int) average[i].starparOffset);
    509     }
    510   }
    511   */
    512 
    513517  /* fix order of StarPar (memory intensive, but fast) */
    514   np = -1;
     518  np = -1; // previous entry in case of errors
    515519  N = 0;
    516520  ALLOCATE (tmpstarpar, StarPar, Nstarpar);
    … …  
    519523    n = average[i].starparOffset;
    520524    average[i].starparOffset = N;
     525    int myObjID = average[i].objID;
    521526    for (k = 0; k < average[i].Nstarpar; k++, N++) {
    522527      if (n == -1) {
    … …  
    525530      }
    526531      tmpstarpar[N] = starpar[n];
    527       if (starpar[n].averef != i) abort();
     532      myAssert (starpar[n].averef == i, "error in averef?");
     533      myAssert (starpar[n].objID == myObjID, "error in objID?");
    528534      tmpstarpar[N].averef = i;
    529535      np = n;
    … …  
    535541}
    536542
    537 
     543/*** GalPhot ****************************************************************************************/
     544
     545/* build the initial links assuming the table is sorted,
     546   not partial, and has a correct set of average[].galphotOffset,Ngalphot values */
     547off_t *init_galphot_links (Average *average, off_t Naverage, GalPhot *galphot, off_t Ngalphot) {
     548
     549  off_t i, j, N;
     550  off_t *next_galp;
     551
     552  N = 0;
     553
     554  ALLOCATE (next_galp, off_t, Ngalphot);
     555  for (i = 0; i < Naverage; i++) {
     556    if (!average[i].Ngalphot) continue;
     557    off_t m = average[i].galphotOffset;
     558    myAssert (galphot[m].averef == i, "not sorted");
     559    for (j = 0; j < average[i].Ngalphot - 1; j++, N++) {
     560      myAssert (galphot[m+j+1].averef == i, "not sorted");
     561      next_galp[N] = N + 1;
     562      if (N >= Ngalphot) {
     563        fprintf (stderr, "WARNING: N out of bounds (1)\n");
     564      }
     565    }
     566    next_galp[N] = -1;
     567    if (N >= Ngalphot) {
     568      fprintf (stderr, "WARNING: N out of bounds (2)\n");
     569    }
     570
     571    if (N >= Ngalphot) {
     572      fprintf (stderr, "overflow in init_galphot_links\n");
     573      abort ();
     574    }
     575    N++;
     576  }
     577  return (next_galp);
     578}
     579
     580/* construct galphot links which are valid FOR THIS LOAD
     581 * - if we have a full load, we will get links which can
     582 *   be used by other programs (eg, relphot, etc)
     583 * - if we have a partial load, the links are only valid
     584 *   for that partial load
     585 */
     586
     587off_t *build_galphot_links (Average *average, off_t Naverage, GalPhot *galphot, off_t Ngalphot) {
     588
     589  off_t i, m, k, Nm, averef;
     590  off_t *next_galp;
     591
     592  ALLOCATE (next_galp, off_t, Ngalphot);
     593
     594  /* reset the Nm, offset values for average */
     595  for (i = 0; i < Naverage; i++) {
     596    average[i].galphotOffset = -1;
     597    average[i].Ngalphot     =  0;
     598  }
     599
     600  for (Nm = 0; Nm < Ngalphot; Nm++) {
     601    averef = galphot[Nm].averef;
     602    m = average[averef].galphotOffset; 
     603    next_galp[Nm] = -1;
     604
     605    if (m == -1) { /* no links yet for source */
     606      average[averef].galphotOffset = Nm;
     607      average[averef].Ngalphot     = 1;
     608      continue;
     609    }
     610
     611    for (k = 0; next_galp[m] != -1; k++) {
     612      m = next_galp[m];
     613      if (m >= Ngalphot) {
     614        fprintf (stderr, "WARNING: m out of bounds (1)\n");
     615      }
     616    }
     617
     618    average[averef].Ngalphot = k + 2;
     619    next_galp[m] = Nm;
     620    if (m >= Ngalphot) {
     621      fprintf (stderr, "WARNING: m out of bounds (2)\n");
     622    }
     623  }
     624  return (next_galp);
     625}
     626
     627/* average[].galphotOffset, average[].Ngalphot are valid within an addstar run */
     628int add_galp_link (Average *average, off_t *next_galp, off_t Ngalphot, off_t NGALPHOT) {
     629
     630  off_t k, m;
     631
     632  /* if we have trouble, check validity of next_galp[m] : m < Ngalphot */
     633  m = average[0].galphotOffset; 
     634
     635  for (k = 0; k < average[0].Ngalphot - 1; k++)  {
     636    m = next_galp[m];
     637    if (m >= NGALPHOT) {
     638      fprintf (stderr, "WARNING: m out of bounds (3)\n");
     639    }
     640  }
     641
     642  /* set up references */
     643  next_galp[Ngalphot] = -1;
     644  if (Ngalphot >= NGALPHOT) {
     645    fprintf (stderr, "WARNING: Ngalphot out of bounds (1)\n");
     646  }
     647
     648  // if Ngalphot is 0, m may have been mis-set; add to the end
     649  if ((average[0].Ngalphot == 0) || (m == -1)) {
     650    average[0].galphotOffset = Ngalphot;
     651  } else {
     652    next_galp[m] = Ngalphot;
     653    if (m >= NGALPHOT) {
     654      fprintf (stderr, "WARNING: m out of bounds (4)\n");
     655    }
     656  }
     657
     658  return (TRUE);
     659}
     660
     661GalPhot *sort_galphot (Average *average, off_t Naverage, GalPhot *galphot, off_t Ngalphot, off_t *next_galp) {
     662
     663  off_t i, k, n, np, N;
     664  GalPhot *tmpgalphot;
     665
     666  /* fix order of GalPhot (memory intensive, but fast) */
     667  np = -1;
     668  N = 0;
     669  ALLOCATE (tmpgalphot, GalPhot, Ngalphot);
     670  for (i = 0; i < Naverage; i++) {
     671    if (!average[i].Ngalphot) continue;
     672    n = average[i].galphotOffset;
     673    average[i].galphotOffset = N;
     674    int myObjID = average[i].objID;
     675    for (k = 0; k < average[i].Ngalphot; k++, N++) {
     676      if (n == -1) {
     677        fprintf (stderr, "entry after %d has a problem\n", (int) np);
     678        abort();
     679      }
     680      tmpgalphot[N] = galphot[n];
     681      myAssert (galphot[n].averef == i, "error in averef?");
     682      myAssert (galphot[n].objID == myObjID, "error in objID?");
     683      tmpgalphot[N].averef = i;
     684      np = n;
     685      n = next_galp[n];
     686    }
     687  }
     688  free (galphot);
     689  return (tmpgalphot);
     690}
     691

  • branches/eam_branches/ipp-20151113/Ohana/src/dvomerge/src/dvo_image_merge_dbs.c

    r38553 r39220  
    234234    catalog[0].lensing[i].imageID = newID;
    235235  }
     236
     237  // NOTE galphot also carries an imageID, but it is actually the EXTERNAL image ID reference.  Do NOT update it.
     238
    236239  return TRUE;
    237240}

  • branches/eam_branches/ipp-20151113/Ohana/src/dvomerge/src/dvoverify_utils.c

    r38986 r39220  
    350350      m = catalog.average[i].lensobjOffset;
    351351      for (j = 0; j < catalog.average[i].Nlensobj; j++) {
     352        if ((catalog.lensobj[m+j].objID == 0xffffffff) && (catalog.lensobj[m+j].catID == 0xffffffff)) continue;
    352353        objIDsOK &= (catalog.average[i].objID == catalog.lensobj[m+j].objID);
    353354        catIDsOK &= (catalog.average[i].catID == catalog.lensobj[m+j].catID);
     355        if (!(catalog.average[i].objID == catalog.lensobj[m+j].objID)) {
     356          fprintf (stderr, "!");
     357        }
     358        if (!(catalog.average[i].catID == catalog.lensobj[m+j].catID)) {
     359          fprintf (stderr, "?");
     360        }
    354361      }
    355362    }

  • branches/eam_branches/ipp-20151113/Ohana/src/dvomerge/src/merge_catalogs_new.c

    r38441 r39220  
    1111 
    1212  off_t i, j, offset;
    13   off_t NAVERAGE, NMEASURE, NLENSING, NSTARPAR, Naverage, Nmeasure, Nlensing, Nstarpar, NsecfiltIn, NsecfiltOut, Nm;
     13  off_t NAVERAGE, NMEASURE, NLENSING, NSTARPAR, NGALPHOT, Naverage, Nmeasure, Nlensing, Nstarpar, Ngalphot, NsecfiltIn, NsecfiltOut, Nm;
    1414
    1515  Naverage = output[0].Naverage;
    … …  
    1717  Nlensing = output[0].Nlensing;
    1818  Nstarpar = output[0].Nstarpar;
     19  Ngalphot = output[0].Ngalphot;
    1920
    2021  NsecfiltOut = output[0].Nsecfilt;
    … …  
    2526  NLENSING = Nlensing + 1000;
    2627  NSTARPAR = Nstarpar + 1000;
     28  NGALPHOT = Ngalphot + 1000;
    2729
    2830  REALLOCATE (output[0].average, Average, NAVERAGE);
    … …  
    3133  REALLOCATE (output[0].lensing, Lensing, NLENSING);
    3234  REALLOCATE (output[0].starpar, StarPar, NSTARPAR);
     35  REALLOCATE (output[0].galphot, GalPhot, NGALPHOT);
    3336
    3437  // add all of these entries to the existing catalog
    … …  
    100103      output[0].average[Naverage].Nstarpar = Nm;
    101104
     105      Nm = 0;
     106      for (j = 0; j < input[0].average[i].Ngalphot; j++) {
     107          offset = input[0].average[i].galphotOffset + j;
     108
     109          output[0].galphot[Ngalphot] = input[0].galphot[offset];
     110          output[0].galphot[Ngalphot].averef = Naverage;
     111
     112          Ngalphot ++;
     113          Nm ++;
     114          if (Ngalphot == NGALPHOT) {
     115              NGALPHOT += 1000;
     116              REALLOCATE (output[0].galphot, GalPhot, NGALPHOT);
     117          }
     118      }
     119      output[0].average[Naverage].Ngalphot = Nm;
     120
    102121      Naverage ++;
    103122      if (Naverage == NAVERAGE) {
    … …  
    111130  REALLOCATE (output[0].lensing, Lensing, MAX (Nlensing, 1));
    112131  REALLOCATE (output[0].starpar, StarPar, MAX (Nstarpar, 1));
     132  REALLOCATE (output[0].galphot, GalPhot, MAX (Ngalphot, 1));
    113133  REALLOCATE (output[0].secfilt, SecFilt, NsecfiltOut*MAX (Naverage, 1));
    114134  output[0].Naverage = Naverage;
    115135  output[0].Nlensing = Nlensing;
    116136  output[0].Nstarpar = Nstarpar;
     137  output[0].Ngalphot = Ngalphot;
    117138  output[0].Nsecfilt_mem = Naverage * NsecfiltOut;
    118139
    … …  
    124145 
    125146  if (VERBOSE) {
    126       fprintf (stderr, OFF_T_FMT": using "OFF_T_FMT" stars ("OFF_T_FMT" measures, "OFF_T_FMT" lensing) for catalog\n",
     147      fprintf (stderr, OFF_T_FMT": using "OFF_T_FMT" stars ("OFF_T_FMT" measures, "OFF_T_FMT" lensing, "OFF_T_FMT" starpar, "OFF_T_FMT" galphot) for catalog\n",
    127148                i,
    128149                output[0].Naverage,
    129150                output[0].Nmeasure,
    130                 output[0].Nlensing);
     151                output[0].Nlensing,
     152                output[0].Nstarpar,
     153                output[0].Ngalphot);
    131154  }
    132155  return (TRUE);

  • branches/eam_branches/ipp-20151113/Ohana/src/dvomerge/src/merge_catalogs_old.c

    r38986 r39220  
    1414  double *X1, *Y1, *X2, *Y2;
    1515  double dX, dY, dR;
    16   off_t *N1, *N2, *next_meas, *next_lens, *next_star;
    17   off_t Nave, NAVE, Nmeas, NMEAS, Nmatch, Nlens, NLENS, Nstar, NSTAR;
     16  off_t *N1, *N2, *next_meas, *next_lens, *next_star, *next_galp;
     17  off_t Nave, NAVE, Nmeas, NMEAS, Nmatch, Nlens, NLENS, Nstar, NSTAR, Ngalp, NGALP;
    1818  int NsecfiltIn;
    1919  int NsecfiltOut;
    … …  
    5454  NLENS = Nlens = output[0].Nlensing;
    5555  NSTAR = Nstar = output[0].Nstarpar;
     56  NGALP = Ngalp = output[0].Ngalphot;
    5657
    5758  // current max obj ID for this catalog
    … …  
    118119    next_lens = init_lensing_links (output[0].average, Nave, output[0].lensing, Nlens);
    119120    next_star = init_starpar_links (output[0].average, Nave, output[0].starpar, Nstar);
     121    next_galp = init_galphot_links (output[0].average, Nave, output[0].galphot, Ngalp);
    120122  } else {
    121123    next_meas = build_measure_links (output[0].average, Nave, output[0].measure, Nmeas);
    122124    next_lens = build_lensing_links (output[0].average, Nave, output[0].lensing, Nlens);
    123125    next_star = build_starpar_links (output[0].average, Nave, output[0].starpar, Nstar);
     126    next_galp = build_galphot_links (output[0].average, Nave, output[0].galphot, Ngalp);
    124127  }   
    125128
    … …  
    193196      REALLOCATE (output[0].starpar, StarPar, NSTAR);
    194197    }
     198    if (Ngalp + input[0].average[N].Ngalphot >= NGALP) {
     199      NGALP = Ngalp + input[0].average[N].Ngalphot + 1000;
     200      REALLOCATE (next_galp, off_t, NGALP);
     201      REALLOCATE (output[0].galphot, GalPhot, NGALP);
     202    }
    195203
    196204    // 4) average properties from the input and the output db need to be properly merged.
    … …  
    214222      if (REPLACE_BY_PHOTCODE) {
    215223        // index to first measure for this object
    216         // XXX this does not support lensing measurements
     224        // XXX this does not support lensing, starpar, or galphot measurements
    217225        if (replace_match (&output[0].average[n], output[0].measure, next_meas, &input[0].average[N], &input[0].measure[offset])) {
    218226          input[0].found_t[N] = output[0].average[n].measureOffset; 
    … …  
    307315    }
    308316
     317    // if galphot measurements exist, add them too
     318    if (output[0].galphot) {
     319      for (Nin = 0; Nin < input[0].average[N].Ngalphot; Nin++) {
     320        /* add to end of galphot list */
     321        add_galp_link (&output[0].average[n], next_galp, Ngalp, NGALP);
     322       
     323        // set the new galphot
     324        off_t galpoff = input[0].average[N].galphotOffset + Nin;
     325        output[0].galphot[Ngalp] = input[0].galphot[galpoff];
     326
     327        output[0].galphot[Ngalp].averef   = n;
     328        output[0].galphot[Ngalp].objID    = output[0].average[n].objID;
     329        output[0].galphot[Ngalp].catID    = output[0].catID;
     330        output[0].average[n].Ngalphot ++;
     331        Ngalp ++;
     332      }
     333    }
     334
    309335    // update the average properties to reflect the incoming entries:
    310336    // if the original value is NAN but the input value is not, accept the input:
    … …  
    357383      REALLOCATE (output[0].starpar, StarPar, NSTAR);
    358384    }
     385    if (Ngalp + input[0].average[N].Ngalphot >= NGALP) {
     386      NGALP = Ngalp + input[0].average[N].Ngalphot + 1000;
     387      REALLOCATE (next_galp, off_t, NGALP);
     388      REALLOCATE (output[0].galphot, GalPhot, NGALP);
     389    }
    359390    if (Nave >= NAVE) {
    360391      NAVE = Nave + 1000;
    … …  
    480511    }
    481512
     513    /** add galphot for this input average object **/
     514    if (output[0].galphot) {
     515      output[0].average[Nave].galphotOffset  = Ngalp;
     516      for (Nin = 0; Nin < input[0].average[N].Ngalphot; Nin ++) {
     517        // supply the galphot values from this detection
     518        off_t galpoff = input[0].average[N].galphotOffset + Nin;
     519        output[0].galphot[Ngalp]           = input[0].galphot[galpoff];
     520
     521        // the following galphot elements cannot be set until here:
     522        output[0].galphot[Ngalp].averef   = Nave;
     523        output[0].galphot[Ngalp].objID    = output[0].average[Nave].objID;
     524        output[0].galphot[Ngalp].catID    = output[0].catID;
     525
     526        // as we add galphot, update Ngalphot to match
     527        output[0].average[Nave].Ngalphot ++;
     528
     529        /* we set next[Ngalp] to -1 here, and update correctly below */
     530        next_galp[Ngalp] = -1;
     531        Ngalp ++;
     532      }
     533      int Ngroup = input[0].average[N].Ngalphot;
     534      for (j = 0; j < Ngroup - 1; j++) {
     535        next_galp[Ngalp - Ngroup + j] = Ngalp - Ngroup + j + 1;
     536      }
     537    }
     538
    482539    Nave ++;
    483540  }
    … …  
    489546  REALLOCATE (output[0].lensing, Lensing, Nlens);
    490547  REALLOCATE (output[0].starpar, StarPar, Nstar);
     548  REALLOCATE (output[0].galphot, GalPhot, Ngalp);
    491549 
    492550# define NOSORT 0
    … …  
    498556    output[0].lensing = sort_lensing (output[0].average, Nave, output[0].lensing, Nlens, next_lens);
    499557    output[0].starpar = sort_starpar (output[0].average, Nave, output[0].starpar, Nstar, next_star);
     558    output[0].galphot = sort_galphot (output[0].average, Nave, output[0].galphot, Ngalp, next_galp);
    500559  }
    501560
    … …  
    515574  output[0].Nlensing = Nlens;
    516575  output[0].Nstarpar = Nstar;
     576  output[0].Ngalphot = Ngalp;
    517577  output[0].Nsecfilt_mem = Nave*NsecfiltOut;
    518   if (VERBOSE) fprintf (stderr, "Nstars, Nave, Nmeas, Nlens: "OFF_T_FMT" "OFF_T_FMT" "OFF_T_FMT" "OFF_T_FMT", ("OFF_T_FMT" matches)\n",  Nstars,  Nave,  Nmeas,  Nlens, Nmatch);
     578  if (VERBOSE) fprintf (stderr, "Nstars, Nave, Nmeas, Nlens, Ngalp: "OFF_T_FMT" "OFF_T_FMT" "OFF_T_FMT" "OFF_T_FMT" "OFF_T_FMT", ("OFF_T_FMT" matches)\n",  Nstars,  Nave,  Nmeas,  Nlens, Ngalp, Nmatch);
    519579
    520580  free (next_meas);
    521581  free (next_lens);
    522582  free (next_star);
     583  free (next_galp);
    523584
    524585  free (X2);
Note: See TracChangeset for help on using the changeset viewer.