IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 37084


Ignore:
Timestamp:
Jul 17, 2014, 4:21:39 PM (12 years ago)
Author:
eugene
Message:

add lensing parameters

Location:
branches/eam_branches/ipp-20140717/Ohana/src/dvomerge
Files:
11 edited

Legend:

Unmodified
Added
Removed
  • branches/eam_branches/ipp-20140717/Ohana/src/dvomerge/include/dvomerge.h

    r35765 r37084  
    129129int        merge_catalogs_old     PROTO((SkyRegion *region, Catalog *output, Catalog *input, double RADIUS, int *secflitMap));
    130130
     131off_t     *build_measure_links    PROTO((Average *average, off_t Naverage, Measure *measure, off_t Nmeasure));
    131132off_t     *init_measure_links     PROTO((Average *average, off_t Naverage, Measure *measure, off_t Nmeasure));
     133int        add_meas_link          PROTO((Average *average, off_t *next, off_t Nmeasure, off_t NMEASURE));
     134Measure   *sort_measure           PROTO((Average *average, off_t Naverage, Measure *measure, off_t Nmeasure, off_t *next));
     135
     136off_t     *build_lensing_links    PROTO((Average *average, off_t Naverage, Lensing *lensing, off_t Nlensing));
     137off_t     *init_lensing_links     PROTO((Average *average, off_t Naverage, Lensing *lensing, off_t Nlensing));
     138int        add_lens_link          PROTO((Average *average, off_t *next, off_t Nlensing, off_t NLENSING));
     139Lensing   *sort_lensing           PROTO((Average *average, off_t Naverage, Lensing *lensing, off_t Nlensing, off_t *next));
     140
    132141off_t     *init_missing_links     PROTO((Average *average, off_t Naverage, Missing *missing, off_t Nmissing));
    133 off_t     *build_measure_links    PROTO((Average *average, off_t Naverage, Measure *measure, off_t Nmeasure));
    134 int        add_meas_link          PROTO((Average *average, off_t *next, off_t Nmeasure, off_t NMEASURE));
    135142int        add_miss_link          PROTO((Average *average, off_t *next, off_t Nmissing));
    136 Measure   *sort_measure           PROTO((Average *average, off_t Naverage, Measure *measure, off_t Nmeasure, off_t *next));
    137143Missing   *sort_missing           PROTO((Average *average, off_t Naverage, Missing *missing, off_t Nmissing, off_t *next_miss));
    138144
  • branches/eam_branches/ipp-20140717/Ohana/src/dvomerge/src/LoadCatalog.c

    r29938 r37084  
    88
    99    // always load all of the data (if any exists)
    10     catalog[0].catflags = LOAD_AVES | LOAD_MEAS | LOAD_MISS | LOAD_SECF;
     10    catalog[0].catflags = LOAD_AVES | LOAD_MEAS | LOAD_MISS | LOAD_SECF | LOAD_LENSING | LOAD_LENSOBJ;
    1111
    1212    catalog[0].catformat = dvo_catalog_catformat (CATFORMAT);  // set the default catformat from config data
  • branches/eam_branches/ipp-20140717/Ohana/src/dvomerge/src/build_links.c

    r27435 r37084  
    11# include "dvomerge.h"
     2
     3/*
     4
     5There are two modes for the measure table: sorted and unsorted.
     6
     7In sorted mode, all measures associated with a given average are in a single block.  The block
     8is pointed to by average->measureOffset and the range is average->Nmeasure
     9
     10In unsorted mode, it is not possible to go directly from average to measure without scanning.
     11In this case, it is necessary to use the value measure->averef to find the corresponding
     12average entry. 
     13
     14Note that average->measureOffset and measure->averef are only valid for a given load of the
     15data: they refer to the sequence number in the data blocks.
     16
     17next_meas is a list of the equivalent sequence of the measure block as if it were sorted.
     18
     19to find the sequence of measurements for a given average:
     20n_0 = average->measureOffset
     21n_1 = next_meas[n_0]
     22n_i = next_meas[n_i-1]
     23
     24*/
    225
    326/* build the initial links assuming the table is sorted,
     
    629
    730  off_t i, j, N;
    8   off_t *next;
     31  off_t *next_meas;
    932
    1033  N = 0;
    1134
    12   ALLOCATE (next, off_t, Nmeasure);
     35  ALLOCATE (next_meas, off_t, Nmeasure);
    1336  for (i = 0; i < Naverage; i++, N++) {
    1437    for (j = 0; j < average[i].Nmeasure - 1; j++, N++) {
    15       next[N] = N + 1;
     38      next_meas[N] = N + 1;
    1639      if (N >= Nmeasure) {
    1740        fprintf (stderr, "WARNING: N out of bounds (1)\n");
    1841      }
    1942    }
    20     next[N] = -1;
     43    next_meas[N] = -1;
    2144    if (N >= Nmeasure) {
    2245      fprintf (stderr, "WARNING: N out of bounds (2)\n");
     
    2851    }
    2952  }
    30   return (next);
    31 }
    32 
    33 /* average[].measureOffset, average[].Nmeasure are valid within an addstar run */
    34 int add_meas_link (Average *average, off_t *next, off_t Nmeasure, off_t NMEASURE) {
    35 
    36   off_t k, m;
    37 
    38   /* if we have trouble, check validity of next[m] : m < Nmeasure */
    39   m = average[0].measureOffset; 
    40 
    41   for (k = 0; k < average[0].Nmeasure - 1; k++)  {
    42     m = next[m];
    43     if (m >= NMEASURE) {
    44       fprintf (stderr, "WARNING: m out of bounds (3)\n");
    45     }
    46   }
    47 
    48   /* set up references */
    49   next[Nmeasure] = -1;
    50   if (Nmeasure >= NMEASURE) {
    51     fprintf (stderr, "WARNING: Nmeasure out of bounds (1)\n");
    52   }
    53 
    54   if (m == -1) {
    55     average[0].measureOffset = Nmeasure;
    56   } else {
    57     next[m] = Nmeasure;
    58     if (m >= NMEASURE) {
    59       fprintf (stderr, "WARNING: m out of bounds (4)\n");
    60     }
    61   }
    62 
    63   return (TRUE);
    64 }
    65 
    66 Measure *sort_measure (Average *average, off_t Naverage, Measure *measure, off_t Nmeasure, off_t *next) {
    67 
    68   off_t i, k, n, N;
    69   Measure *tmpmeasure;
    70 
    71   /* fix order of Measure (memory intensive, but fast) */
    72   N = 0;
    73   ALLOCATE (tmpmeasure, Measure, Nmeasure);
    74   for (i = 0; i < Naverage; i++) {
    75     n = average[i].measureOffset;
    76     average[i].measureOffset = N;
    77     for (k = 0; k < average[i].Nmeasure; k++, N++) {
    78       tmpmeasure[N] = measure[n];
    79       tmpmeasure[N].averef = i;
    80       n = next[n];
    81     }
    82   }
    83   free (measure);
    84   return (tmpmeasure);
     53  return (next_meas);
    8554}
    8655
     
    9564
    9665  off_t i, m, k, Nm, averef;
    97   off_t *next;
    98 
    99   ALLOCATE (next, off_t, Nmeasure);
     66  off_t *next_meas;
     67
     68  ALLOCATE (next_meas, off_t, Nmeasure);
    10069
    10170  /* reset the Nm, offset values for average */
     
    10877    averef = measure[Nm].averef;
    10978    m = average[averef].measureOffset; 
    110     next[Nm] = -1;
     79    next_meas[Nm] = -1;
    11180
    11281    if (m == -1) { /* no links yet for source */
     
    11685    }
    11786
    118     for (k = 0; next[m] != -1; k++) {
    119       m = next[m];
     87    for (k = 0; next_meas[m] != -1; k++) {
     88      m = next_meas[m];
    12089      if (m >= Nmeasure) {
    12190        fprintf (stderr, "WARNING: m out of bounds (1)\n");
     
    12493
    12594    average[averef].Nmeasure = k + 2;
    126     next[m] = Nm;
     95    next_meas[m] = Nm;
    12796    if (m >= Nmeasure) {
    12897      fprintf (stderr, "WARNING: m out of bounds (2)\n");
    12998    }
    13099  }
    131   return (next);
    132 }
     100  return (next_meas);
     101}
     102
     103/* average[].measureOffset, average[].Nmeasure are valid within an addstar run */
     104int add_meas_link (Average *average, off_t *next_meas, off_t Nmeasure, off_t NMEASURE) {
     105
     106  off_t k, m;
     107
     108  /* if we have trouble, check validity of next_meas[m] : m < Nmeasure */
     109  m = average[0].measureOffset; 
     110
     111  for (k = 0; k < average[0].Nmeasure - 1; k++)  {
     112    m = next_meas[m];
     113    if (m >= NMEASURE) {
     114      fprintf (stderr, "WARNING: m out of bounds (3)\n");
     115    }
     116  }
     117
     118  /* set up references */
     119  next_meas[Nmeasure] = -1;
     120  if (Nmeasure >= NMEASURE) {
     121    fprintf (stderr, "WARNING: Nmeasure out of bounds (1)\n");
     122  }
     123
     124  if (m == -1) {
     125    average[0].measureOffset = Nmeasure;
     126  } else {
     127    next_meas[m] = Nmeasure;
     128    if (m >= NMEASURE) {
     129      fprintf (stderr, "WARNING: m out of bounds (4)\n");
     130    }
     131  }
     132
     133  return (TRUE);
     134}
     135
     136Measure *sort_measure (Average *average, off_t Naverage, Measure *measure, off_t Nmeasure, off_t *next_meas) {
     137
     138  off_t i, k, n, N;
     139  Measure *tmpmeasure;
     140
     141  /* fix order of Measure (memory intensive, but fast) */
     142  N = 0;
     143  ALLOCATE (tmpmeasure, Measure, Nmeasure);
     144  for (i = 0; i < Naverage; i++) {
     145    n = average[i].measureOffset;
     146    average[i].measureOffset = N;
     147    for (k = 0; k < average[i].Nmeasure; k++, N++) {
     148      if (n == -1) abort();
     149      tmpmeasure[N] = measure[n];
     150      if (measure[n].averef != i) abort();
     151      tmpmeasure[N].averef = i;
     152      n = next_meas[n];
     153    }
     154  }
     155  free (measure);
     156  return (tmpmeasure);
     157}
     158
     159/*******************************************************************************************/
    133160
    134161/* build the initial links assuming the table is sorted */
     
    136163
    137164  off_t i, j, N;
    138   off_t *next;
     165  off_t *next_miss;
    139166
    140167  N = 0;
    141168
    142   ALLOCATE (next, off_t, Nmissing);
     169  ALLOCATE (next_miss, off_t, Nmissing);
    143170  for (i = 0; i < Naverage; i++) {
    144171    for (j = 0; j < average[i].Nmissing - 1; j++, N++) {
    145       next[N] = N + 1;
     172      next_miss[N] = N + 1;
    146173    }
    147174    if (average[i].Nmissing > 0) {
    148       next[N] = -1;
     175      next_miss[N] = -1;
    149176      if (N >= Nmissing) {
    150177        fprintf (stderr, "overflow in init_missing_links");
     
    155182
    156183  }
    157   return (next);
    158 }
    159 
    160 int add_miss_link (Average *average, off_t *next, off_t Nmissing) {
     184  return (next_miss);
     185}
     186
     187int add_miss_link (Average *average, off_t *next_miss, off_t Nmissing) {
    161188
    162189  off_t k, m;
     
    165192  if (average[0].Nmissing < 1) {
    166193    average[0].missingOffset = Nmissing;
    167     next[Nmissing] = -1;
     194    next_miss[Nmissing] = -1;
    168195    return (TRUE);
    169196  }
    170197
    171198  m = average[0].missingOffset; 
    172   for (k = 0; k < average[0].Nmissing - 1; k++) m = next[m];
     199  for (k = 0; k < average[0].Nmissing - 1; k++) m = next_miss[m];
    173200  /* set up references */
    174   next[Nmissing] = -1;
    175   next[m] = Nmissing;
     201  next_miss[Nmissing] = -1;
     202  next_miss[m] = Nmissing;
    176203  return (TRUE);
    177204}
     
    180207   we must always save the missing table, if it exists */
    181208
    182 Missing *sort_missing (Average *average, off_t Naverage, Missing *missing, off_t Nmissing, off_t *next) {
     209Missing *sort_missing (Average *average, off_t Naverage, Missing *missing, off_t Nmissing, off_t *next_miss) {
    183210
    184211  off_t i, k, n, N;
     
    193220    for (k = 0; k < average[i].Nmissing; k++, N++) {
    194221      tmpmissing[N] = missing[n];
    195       n = next[n];
     222      n = next_miss[n];
    196223    }
    197224  }
     
    199226  return (tmpmissing);
    200227}
     228
     229/*******************************************************************************************/
     230
     231/* build the initial links assuming the table is sorted,
     232   not partial, and has a correct set of average[].lensingOffset,Nlensing values */
     233off_t *init_lensing_links (Average *average, off_t Naverage, Lensing *lensing, off_t Nlensing) {
     234
     235  off_t i, j, N;
     236  off_t *next_lens;
     237
     238  N = 0;
     239
     240  ALLOCATE (next_lens, off_t, Nlensing);
     241  for (i = 0; i < Naverage; i++) {
     242    if (!average[i].Nlensing) continue;
     243    for (j = 0; j < average[i].Nlensing - 1; j++, N++) {
     244      next_lens[N] = N + 1;
     245      if (N >= Nlensing) {
     246        fprintf (stderr, "WARNING: N out of bounds (1)\n");
     247      }
     248    }
     249    next_lens[N] = -1;
     250    if (N >= Nlensing) {
     251      fprintf (stderr, "WARNING: N out of bounds (2)\n");
     252    }
     253
     254    if (N >= Nlensing) {
     255      fprintf (stderr, "overflow in init_lensing_links\n");
     256      abort ();
     257    }
     258    N++;
     259  }
     260  return (next_lens);
     261}
     262
     263/* construct lensing links which are valid FOR THIS LOAD
     264 * - if we have a full load, we will get links which can
     265 *   be used by other programs (eg, relphot, etc)
     266 * - if we have a partial load, the links are only valid
     267 *   for that partial load
     268 */
     269
     270off_t *build_lensing_links (Average *average, off_t Naverage, Lensing *lensing, off_t Nlensing) {
     271
     272  off_t i, m, k, Nm, averef;
     273  off_t *next_lens;
     274
     275  ALLOCATE (next_lens, off_t, Nlensing);
     276
     277  /* reset the Nm, offset values for average */
     278  for (i = 0; i < Naverage; i++) {
     279    average[i].lensingOffset = -1;
     280    average[i].Nlensing     =  0;
     281  }
     282
     283  for (Nm = 0; Nm < Nlensing; Nm++) {
     284    averef = lensing[Nm].averef;
     285    m = average[averef].lensingOffset; 
     286    next_lens[Nm] = -1;
     287
     288    if (m == -1) { /* no links yet for source */
     289      average[averef].lensingOffset = Nm;
     290      average[averef].Nlensing     = 1;
     291      continue;
     292    }
     293
     294    for (k = 0; next_lens[m] != -1; k++) {
     295      m = next_lens[m];
     296      if (m >= Nlensing) {
     297        fprintf (stderr, "WARNING: m out of bounds (1)\n");
     298      }
     299    }
     300
     301    average[averef].Nlensing = k + 2;
     302    next_lens[m] = Nm;
     303    if (m >= Nlensing) {
     304      fprintf (stderr, "WARNING: m out of bounds (2)\n");
     305    }
     306  }
     307  return (next_lens);
     308}
     309
     310/* average[].lensingOffset, average[].Nlensing are valid within an addstar run */
     311int add_lens_link (Average *average, off_t *next_lens, off_t Nlensing, off_t NLENSING) {
     312
     313  off_t k, m;
     314
     315  /* if we have trouble, check validity of next_lens[m] : m < Nlensing */
     316  m = average[0].lensingOffset; 
     317
     318  for (k = 0; k < average[0].Nlensing - 1; k++)  {
     319    m = next_lens[m];
     320    if (m >= NLENSING) {
     321      fprintf (stderr, "WARNING: m out of bounds (3)\n");
     322    }
     323  }
     324
     325  /* set up references */
     326  next_lens[Nlensing] = -1;
     327  if (Nlensing >= NLENSING) {
     328    fprintf (stderr, "WARNING: Nlensing out of bounds (1)\n");
     329  }
     330
     331  if (m == -1) {
     332    average[0].lensingOffset = Nlensing;
     333  } else {
     334    next_lens[m] = Nlensing;
     335    if (m >= NLENSING) {
     336      fprintf (stderr, "WARNING: m out of bounds (4)\n");
     337    }
     338  }
     339
     340  return (TRUE);
     341}
     342
     343Lensing *sort_lensing (Average *average, off_t Naverage, Lensing *lensing, off_t Nlensing, off_t *next_lens) {
     344
     345  off_t i, k, n, N;
     346  Lensing *tmplensing;
     347
     348  /* fix order of Lensing (memory intensive, but fast) */
     349  N = 0;
     350  ALLOCATE (tmplensing, Lensing, Nlensing);
     351  for (i = 0; i < Naverage; i++) {
     352    if (!average[i].Nlensing) continue;
     353    n = average[i].lensingOffset;
     354    average[i].lensingOffset = N;
     355    for (k = 0; k < average[i].Nlensing; k++, N++) {
     356      if (n == -1) abort();
     357      tmplensing[N] = lensing[n];
     358      if (lensing[n].averef != i) abort();
     359      tmplensing[N].averef = i;
     360      n = next_lens[n];
     361    }
     362  }
     363  free (lensing);
     364  return (tmplensing);
     365}
     366
     367
  • branches/eam_branches/ipp-20140717/Ohana/src/dvomerge/src/dvomergeContinue.c

    r34260 r37084  
    131131      merge_catalogs_old (&outsky[0].regions[j], &outcatalog, &incatalog, RADIUS, secfiltMap);
    132132
    133       outcatalog.catflags = LOAD_AVES | LOAD_MEAS | LOAD_MISS | LOAD_SECF;
    134 
    135133      // if we receive a signal which would cause us to exit, wait until the full catalog is written
    136134      SetProtect (TRUE);
  • branches/eam_branches/ipp-20140717/Ohana/src/dvomerge/src/dvomergeContinue_threaded.c

    r34260 r37084  
    7474      merge_catalogs_old (&threadData->outsky->regions[j], &outcatalog, &incatalog, RADIUS, threadData->secfiltMap);
    7575
    76       outcatalog.catflags = LOAD_AVES | LOAD_MEAS | LOAD_MISS | LOAD_SECF;
    77 
    7876      // if we receive a signal which would cause us to exit, wait until the full catalog is written
    7977      SetProtect (TRUE);
  • branches/eam_branches/ipp-20140717/Ohana/src/dvomerge/src/dvomergeCreate.c

    r33963 r37084  
    178178    SkyListFree (inlist);
    179179
    180     outcatalog.catflags = LOAD_AVES | LOAD_MEAS | LOAD_MISS | LOAD_SECF;
    181180    dvo_catalog_save (&outcatalog, VERBOSE);
    182181    dvo_catalog_unlock (&outcatalog);
  • branches/eam_branches/ipp-20140717/Ohana/src/dvomerge/src/dvomergeFromList.c

    r35590 r37084  
    143143    merge_catalogs_old (&skyregion, &outcatalog, &incatalog, RADIUS, secfiltMap);
    144144   
    145     outcatalog.catflags = LOAD_AVES | LOAD_MEAS | LOAD_MISS | LOAD_SECF;
    146 
    147145    dmhObjectAdd (outstat[0].history, &outcatalog.header, inStats);
    148146
  • branches/eam_branches/ipp-20140717/Ohana/src/dvomerge/src/dvomergeUpdate_catalogs.c

    r37046 r37084  
    165165      merge_catalogs_old (outlist[0].regions[j], &outcatalog, &incatalog, RADIUS, secfiltMap);
    166166
    167       outcatalog.catflags = LOAD_AVES | LOAD_MEAS | LOAD_MISS | LOAD_SECF;
    168 
    169167      if (outstat[j].missed) {
    170168        dmhObjectAdd (outstat[j].history, &outcatalog.header, inStats);
  • branches/eam_branches/ipp-20140717/Ohana/src/dvomerge/src/dvomergeUpdate_threaded.c

    r35765 r37084  
    7272      merge_catalogs_old (&threadData->outsky->regions[j], &outcatalog, &incatalog, RADIUS, threadData->secfiltMap);
    7373
    74       outcatalog.catflags = LOAD_AVES | LOAD_MEAS | LOAD_MISS | LOAD_SECF;
    75 
    7674      // if we receive a signal which would cause us to exit, wait until the full catalog is written
    7775      SetProtect (TRUE);
  • branches/eam_branches/ipp-20140717/Ohana/src/dvomerge/src/merge_catalogs_new.c

    r28855 r37084  
    88// input entries always define new objects
    99
    10 #define notyet 1
    11 
    1210int merge_catalogs_new (SkyRegion *region, Catalog *output, Catalog *input, int *secfiltMap) {
    1311 
    1412  off_t i, j, offset;
    15   off_t NAVERAGE, NMEASURE, Naverage, Nmeasure, NsecfiltIn, NsecfiltOut, Nm;
     13  off_t NAVERAGE, NMEASURE, NLENSING, Naverage, Nmeasure, Nlensing, NsecfiltIn, NsecfiltOut, Nm;
    1614
    1715  Naverage = output[0].Naverage;
     
    6361      }
    6462      output[0].average[Naverage].Nmeasure = Nm;
     63
     64      Nm = 0;
     65      for (j = 0; j < input[0].average[i].Nlensing; j++) {
     66          offset = input[0].average[i].lensingOffset + j;
     67
     68          output[0].lensing[Nlensing] = input[0].lensing[offset];
     69          output[0].lensing[Nlensing].averef = Naverage;
     70
     71          Nlensing ++;
     72          Nm ++;
     73          if (Nlensing == NLENSING) {
     74              NLENSING += 1000;
     75              REALLOCATE (output[0].lensing, Lensing, NLENSING);
     76          }
     77      }
     78      output[0].average[Naverage].Nlensing = Nm;
     79
    6580      Naverage ++;
    6681      if (Naverage == NAVERAGE) {
     
    7287  REALLOCATE (output[0].average, Average, MAX (Naverage, 1));
    7388  REALLOCATE (output[0].measure, Measure, MAX (Nmeasure, 1));
     89  REALLOCATE (output[0].lensing, Lensing, MAX (Nlensing, 1));
    7490  REALLOCATE (output[0].secfilt, SecFilt, NsecfiltOut*MAX (Naverage, 1));
    7591  output[0].Naverage = Naverage;
    76   output[0].Nmeasure = Nmeasure;
     92  output[0].Nlensing = Nlensing;
    7793  output[0].Nsecf_mem = Naverage * NsecfiltOut;
    7894
     
    84100 
    85101  if (VERBOSE) {
    86       fprintf (stderr, OFF_T_FMT": using "OFF_T_FMT" stars ("OFF_T_FMT" measures) for catalog\n",
     102      fprintf (stderr, OFF_T_FMT": using "OFF_T_FMT" stars ("OFF_T_FMT" measures, "OFF_T_FMT" lensing) for catalog\n",
    87103                i,
    88104                output[0].Naverage,
    89                 output[0].Nmeasure);
     105                output[0].Nmeasure,
     106                output[0].Nlensing);
    90107  }
    91108  return (TRUE);
  • branches/eam_branches/ipp-20140717/Ohana/src/dvomerge/src/merge_catalogs_old.c

    r37046 r37084  
    1414  double *X1, *Y1, *X2, *Y2;
    1515  double dX, dY, dR;
    16   off_t *N1, *N2,  *next_meas;
    17   off_t Nave, NAVE, Nmeas, NMEAS, Nmatch;
     16  off_t *N1, *N2, *next_meas, *next_lens;
     17  off_t Nave, NAVE, Nmeas, NMEAS, Nmatch, Nlens, NLENS;
    1818  int NsecfiltIn;
    1919  int NsecfiltOut;
     
    5252  Nmatch = 0;
    5353  NMEAS = Nmeas = output[0].Nmeasure;
     54  NLENS = Nlens = output[0].Nlensing;
    5455
    5556  // current max obj ID for this catalog
     
    113114  }
    114115
    115 
    116116  /* set up pointers for linked list of measure */
    117117  if (output[0].sorted && (output[0].Nmeasure >= output[0].Nmeas_disk)) {
     
    119119    // is sorted while processed
    120120    next_meas = init_measure_links (output[0].average, Nave, output[0].measure, Nmeas);
     121    next_lens = init_lensing_links (output[0].average, Nave, output[0].lensing, Nlens);
    121122  } else {
    122123    next_meas = build_measure_links (output[0].average, Nave, output[0].measure, Nmeas);
     124    next_lens = build_lensing_links (output[0].average, Nave, output[0].lensing, Nlens);
    123125  }   
    124126
     
    182184      REALLOCATE (output[0].measure, Measure, NMEAS);
    183185    }
     186    if (Nlens + input[0].average[N].Nmeasure >= NLENS) {
     187      NLENS = Nlens + input[0].average[N].Nmeasure + 1000;
     188      REALLOCATE (next_lens, off_t, NLENS);
     189      REALLOCATE (output[0].lensing, Lensing, NLENS);
     190    }
    184191
    185192    // 4) average properties from the input and the output db need to be properly merged.
     
    191198      if (REPLACE_BY_PHOTCODE) {
    192199        // index to first measure for this object
     200        // XXX this does not support lensing measurements
    193201        int Mout = output[0].average[n].measureOffset; 
    194202        if (replace_match(&output[0].average[n], &output[0].measure[Mout], &input[0].average[N], &input[0].measure[offset])) {
     
    248256    }
    249257
     258    // if lensing measurements exist, add them too
     259    if (output[0].lensing) {
     260      for (Nin = 0; Nin < input[0].average[N].Nlensing; Nin++) {
     261        /* add to end of lensing list */
     262        add_lens_link (&output[0].average[n], next_lens, Nlens, NLENS);
     263       
     264        // set the new lensing
     265        off_t lensoff = input[0].average[N].lensingOffset + Nin;
     266        output[0].lensing[Nlens] = input[0].lensing[lensoff];
     267
     268        output[0].lensing[Nlens].averef   = n;
     269        output[0].lensing[Nlens].objID    = output[0].average[n].objID;
     270        output[0].lensing[Nlens].catID    = output[0].catID;
     271        Nlens ++;
     272      }
     273    }
     274
    250275    // update the average properties to reflect the incoming entries:
    251276    // if the original value is NAN but the input value is not, accept the input:
     
    288313      REALLOCATE (output[0].measure, Measure, NMEAS);
    289314    }
     315    if (Nlens + input[0].average[N].Nlensing >= NLENS) {
     316      NLENS = Nlens + input[0].average[N].Nlensing + 1000;
     317      REALLOCATE (next_lens, off_t, NLENS);
     318      REALLOCATE (output[0].lensing, Lensing, NLENS);
     319    }
    290320    if (Nave >= NAVE) {
    291321      NAVE = Nave + 1000;
     
    355385      next_meas[Nmeas - Ngroup + j] = Nmeas - Ngroup + j + 1;
    356386    }
     387
     388    /** add lensing for this input average object **/
     389    if (output[0].lensing) {
     390      for (Nin = 0; Nin < input[0].average[N].Nlensing; Nin ++) {
     391        // supply the lensing values from this detection
     392        off_t lensoff = input[0].average[N].lensingOffset + Nin;
     393        output[0].lensing[Nlens]           = input[0].lensing[lensoff];
     394
     395        // the following lensing elements cannot be set until here:
     396        output[0].lensing[Nlens].averef   = Nave;
     397        output[0].lensing[Nlens].objID    = output[0].average[Nave].objID;
     398        output[0].lensing[Nlens].catID    = output[0].catID;
     399
     400        // as we add lensing, update Nlensing to match
     401        output[0].average[Nave].Nlensing ++;
     402
     403        /* we set next[Nlens] to -1 here, and update correctly below */
     404        next_lens[Nlens] = -1;
     405        Nlens ++;
     406      }
     407      int Ngroup = input[0].average[N].Nlensing;
     408      for (j = 0; j < Ngroup - 1; j++) {
     409        next_lens[Nlens - Ngroup + j] = Nlens - Ngroup + j + 1;
     410      }
     411    }
     412
    357413    Nave ++;
    358414  }
     
    362418  REALLOCATE (output[0].average, Average, Nave);
    363419  REALLOCATE (output[0].measure, Measure, Nmeas);
     420  REALLOCATE (output[0].lensing, Lensing, Nlens);
    364421 
    365422# define NOSORT 0
     
    369426    output[0].sorted = TRUE;
    370427    output[0].measure = sort_measure (output[0].average, Nave, output[0].measure, Nmeas, next_meas);
     428    output[0].lensing = sort_lensing (output[0].average, Nave, output[0].lensing, Nlens, next_lens);
    371429  }
    372430
     
    384442  output[0].Naverage = Nave;
    385443  output[0].Nmeasure = Nmeas;
     444  output[0].Nlensing = Nlens;
    386445  output[0].Nsecf_mem = Nave*NsecfiltOut;
    387   if (VERBOSE) fprintf (stderr, "Nstars, Nave, Nmeas: "OFF_T_FMT" "OFF_T_FMT" "OFF_T_FMT", ("OFF_T_FMT" matches)\n",  Nstars,  Nave,  Nmeas, Nmatch);
     446  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);
    388447
    389448  free (next_meas);
     449  free (next_lens);
    390450
    391451  free (X2);
Note: See TracChangeset for help on using the changeset viewer.