IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Aug 2, 2012, 4:54:24 PM (14 years ago)
Author:
Serge CHASTEL
Message:

MOPS trails fitting integration (all bugs fixed now)

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/ppTranslate/src/ppMopsRead.c

    r34279 r34284  
    88#include "ppMops.h"
    99
    10 static psMetadata* sort_by_increasing_idet(psMetadata *fittedTrail, long length, long fillLength);
     10static void addDummyValues(psMetadata* md, long size);
     11static void replaceDummyValuesF32(const char* colName, psMetadata* source, psMetadata* target, psVector* indexes);
    1112
    1213/*
     
    108109      return NULL;
    109110    }
    110     //Sort table by idet if necessary
    111     psMetadata* sortedTable = sort_by_increasing_idet(table, det->num, det->num);
    112     if (sortedTable != table) {
    113       psFree(table);
    114       table = sortedTable;
    115     } // else det->table is sorted
    116111    det->table = table;
    117     psTrace("ppMops.read", 9, "PSF table:\n%s\n", psMetadataConfigFormat(det->table));
    118     psTrace("ppMops.read", 9, "End of PSF table\n");
     112    psTrace("ppMops.read", 19, "PSF table:\n%s\n", psMetadataConfigFormat(det->table));
     113    psTrace("ppMops.read", 19, "End of PSF table\n");
    119114
    120115    //fittedTrail extension (Supported for releases 2 and above)
    121116    if (skyChipPsfVersion >= 2) {
     117      //First: append all the new columns that we want to the existing table
     118      addDummyValues(det->table, size);
    122119      if (!psFitsMoveExtName(fits, "SkyChip.xfit")) {
    123120        psTrace("ppMops.read", 3, "No fitted trails extension");
    124121      } else {
    125122        psTrace("ppMops.read", 3, "Fitted trails extension found\n");
    126         det->fittedTrailsSize = psFitsTableSize(fits);
    127         if (det->fittedTrailsSize <= 0) {
     123        int fittedTrailsSize = psFitsTableSize(fits);
     124        if (fittedTrailsSize <= 0) {
    128125          psErrorStackPrint(stderr, "Unable to determine size of fitted trails extension table %d", i);
    129126          psTrace("ppMops.read", 3, "No entry in fitted trails extension!!!!\n");
    130127          psErrorClear();
    131           det->fittedTrails = NULL;
    132128        } else {
    133           det->fittedTrails = psFitsReadTableAllColumns(fits); // Table of interest
    134           if (!det->fittedTrails) {
    135             psError(PS_ERR_IO, false, "Unable to read fittedTrails table %d", i);
     129          psMetadata* fittedTrails = psFitsReadTableAllColumns(fits); // Table of interest
     130          if (!fittedTrails) {
     131            psError(PS_ERR_IO, false, "Unable to read fittedTrails table in file %d", i);
    136132            return NULL;
    137133          }
    138           psMetadata* sortedFittedTrails = sort_by_increasing_idet(det->fittedTrails, det->fittedTrailsSize, det->num);
    139           if (sortedFittedTrails != det->fittedTrails) {
    140             psFree(det->fittedTrails);
    141             det->fittedTrails = sortedFittedTrails;
    142           } // else det->fittedTrails is sorted
    143           psTrace("ppMops.read", 9, "Fitted trail:\n%s\n", psMetadataConfigFormat(det->fittedTrails));
    144           psTrace("ppMops.read", 9, "End of Fitted trail\n");
     134          //Iterate on the different names and types expected in the fittedTrails parameters
     135          psVector* idet = psMetadataLookupVector(NULL, fittedTrails, "IPP_IDET");
     136          replaceDummyValuesF32("X_EXT", fittedTrails, det->table, idet);
     137          replaceDummyValuesF32("Y_EXT", fittedTrails, det->table, idet);
     138          replaceDummyValuesF32("X_EXT_SIG", fittedTrails, det->table, idet);
     139          replaceDummyValuesF32("Y_EXT_SIG", fittedTrails, det->table, idet);
     140          replaceDummyValuesF32("EXT_INST_MAG",  fittedTrails, det->table, idet);
     141          replaceDummyValuesF32("EXT_INST_MAG_SIG",  fittedTrails, det->table, idet);
     142          replaceDummyValuesF32("NPARAMS",  fittedTrails, det->table, idet);
     143          replaceDummyValuesF32("EXT_WIDTH_MAJ",  fittedTrails, det->table, idet);
     144          replaceDummyValuesF32("EXT_WIDTH_MIN",  fittedTrails, det->table, idet);
     145          replaceDummyValuesF32("EXT_THETA",  fittedTrails, det->table, idet);
     146          replaceDummyValuesF32("EXT_WIDTH_MAJ_ERR",  fittedTrails, det->table, idet);
     147          replaceDummyValuesF32("EXT_WIDTH_MIN_ERR",  fittedTrails, det->table, idet);
     148          replaceDummyValuesF32("EXT_THETA_ERR",  fittedTrails, det->table, idet);
    145149        }
    146150      }
     
    198202
    199203      // Calculate error in RA, Dec
    200 
    201204      double xErr = xErrV->data.F32[row];
    202205      double yErr = yErrV->data.F32[row];
     
    248251}
    249252
    250 static long* sorted_idets_indexes(psVector *idets);
    251 static int compare (const void *a, const void *b);
    252 static psVector* fill_F32_vector(psVector* in, int size, long* indexes, int indexes_length);
    253 
    254 // We want the entries in table/fittedTrails to be sorted by increasing idet
    255 static psMetadata* sort_by_increasing_idet(psMetadata *aTable, long length, long fillLength) {
    256   if (aTable == NULL) {
    257     return aTable;
     253static psVector* createDummyF32(long size) {
     254  psVector* dummy = psVectorAlloc(size, PS_TYPE_F32);
     255  psVectorInit(dummy, NAN);
     256  return dummy;
     257}
     258
     259static void addDummyValues(psMetadata* md, long size) {
     260  psMetadataAdd(md, PS_LIST_TAIL, "X_EXT", PS_DATA_VECTOR, "", createDummyF32(size));
     261  psMetadataAdd(md, PS_LIST_TAIL, "Y_EXT", PS_DATA_VECTOR, "", createDummyF32(size));
     262  psMetadataAdd(md, PS_LIST_TAIL, "X_EXT_SIG", PS_DATA_VECTOR, "", createDummyF32(size));
     263  psMetadataAdd(md, PS_LIST_TAIL, "Y_EXT_SIG", PS_DATA_VECTOR, "", createDummyF32(size));
     264  psMetadataAdd(md, PS_LIST_TAIL, "EXT_INST_MAG",  PS_DATA_VECTOR, "", createDummyF32(size));
     265  psMetadataAdd(md, PS_LIST_TAIL, "EXT_INST_MAG_SIG",  PS_DATA_VECTOR, "", createDummyF32(size));
     266  psMetadataAdd(md, PS_LIST_TAIL, "NPARAMS",  PS_DATA_VECTOR, "", createDummyF32(size));
     267  psMetadataAdd(md, PS_LIST_TAIL, "EXT_WIDTH_MAJ",  PS_DATA_VECTOR, "", createDummyF32(size));
     268  psMetadataAdd(md, PS_LIST_TAIL, "EXT_WIDTH_MIN",  PS_DATA_VECTOR, "", createDummyF32(size));
     269  psMetadataAdd(md, PS_LIST_TAIL, "EXT_THETA",  PS_DATA_VECTOR, "", createDummyF32(size));
     270  psMetadataAdd(md, PS_LIST_TAIL, "EXT_WIDTH_MAJ_ERR",  PS_DATA_VECTOR, "", createDummyF32(size));
     271  psMetadataAdd(md, PS_LIST_TAIL, "EXT_WIDTH_MIN_ERR",  PS_DATA_VECTOR, "", createDummyF32(size));
     272  psMetadataAdd(md, PS_LIST_TAIL, "EXT_THETA_ERR",  PS_DATA_VECTOR, "", createDummyF32(size));
     273}
     274
     275static void replaceDummyValuesF32(const char* colName, psMetadata* source, psMetadata* target, psVector* indexes) {
     276  psVector* source_vector = psMetadataLookupVector(NULL, source, colName);
     277  psVector* target_vector = psMetadataLookupVector(NULL, target, colName);
     278  for (long index = 0; index<indexes->n; ++index) {
     279    target_vector->data.F32[indexes->data.S64[index]] = source_vector->data.F32[index];
    258280  }
    259   psVector *idets = psMetadataLookupVector(NULL, aTable, "IPP_IDET");
    260   if (!idets) {
    261     psError(PS_ERR_PROGRAMMING, true, "Failed to find IPP_IDET column");
    262     return NULL;
    263   }
    264   bool is_sorted = true;
    265   int count = 0;
    266   long previous_value = -1;
    267   while (is_sorted && (count<length) ){
    268     is_sorted = (previous_value < idets->data.S64[count]);
    269     previous_value = idets->data.S64[count];
    270     count += 1;
    271   }
    272   if (is_sorted && (length==fillLength)) {
    273     psTrace("ppMops.read", 5, "Values are already sorted by idet and the table is complete\n");
    274     return aTable;
    275   }
    276   //otherwise sort (even if it should not happen) and fill
    277   long* sortedIndexVector = sorted_idets_indexes(idets);
    278   psMetadata* sortedTable = psMetadataAlloc();
    279   char* column_name = "X_EXT";
    280   psMetadataAddVector( sortedTable, PS_LIST_TAIL, column_name, 0, NULL,
    281                        fill_F32_vector(psMetadataLookupVector(NULL, aTable, column_name), fillLength, sortedIndexVector, idets->n) );
    282   column_name = "Y_EXT";
    283   psMetadataAddVector( sortedTable, PS_LIST_TAIL, column_name, 0, NULL,
    284                        fill_F32_vector(psMetadataLookupVector(NULL, aTable, column_name), fillLength, sortedIndexVector, idets->n) );
    285   column_name = "X_EXT_SIG";
    286   psMetadataAddVector( sortedTable, PS_LIST_TAIL, column_name, 0, NULL,
    287                        fill_F32_vector(psMetadataLookupVector(NULL, aTable, column_name), fillLength, sortedIndexVector, idets->n) );
    288   column_name = "Y_EXT_SIG";
    289   psMetadataAddVector( sortedTable, PS_LIST_TAIL, column_name, 0, NULL,
    290                        fill_F32_vector(psMetadataLookupVector(NULL, aTable, column_name), fillLength, sortedIndexVector, idets->n) );
    291   column_name = "EXT_INST_MAG";
    292   psMetadataAddVector( sortedTable, PS_LIST_TAIL, column_name, 0, NULL,
    293                        fill_F32_vector(psMetadataLookupVector(NULL, aTable, column_name), fillLength, sortedIndexVector, idets->n) );
    294   column_name = "EXT_INST_MAG_SIG";
    295   psMetadataAddVector( sortedTable, PS_LIST_TAIL, column_name, 0, NULL,
    296                        fill_F32_vector(psMetadataLookupVector(NULL, aTable, column_name), fillLength, sortedIndexVector, idets->n) );
    297   column_name = "NPARAMS";
    298   psMetadataAddVector( sortedTable, PS_LIST_TAIL, column_name, 0, NULL,
    299                        fill_F32_vector(psMetadataLookupVector(NULL, aTable, column_name), fillLength, sortedIndexVector, idets->n) );
    300   column_name = "EXT_WIDTH_MAJ";
    301   psMetadataAddVector( sortedTable, PS_LIST_TAIL, column_name, 0, NULL,
    302                        fill_F32_vector(psMetadataLookupVector(NULL, aTable, column_name), fillLength, sortedIndexVector, idets->n) );
    303   column_name = "EXT_WIDTH_MIN";
    304   psMetadataAddVector( sortedTable, PS_LIST_TAIL, column_name, 0, NULL,
    305                        fill_F32_vector(psMetadataLookupVector(NULL, aTable, column_name), fillLength, sortedIndexVector, idets->n) );
    306   column_name = "EXT_THETA";
    307   psMetadataAddVector( sortedTable, PS_LIST_TAIL, column_name, 0, NULL,
    308                        fill_F32_vector(psMetadataLookupVector(NULL, aTable, column_name), fillLength, sortedIndexVector, idets->n) );
    309   column_name = "EXT_WIDTH_MAJ_ERR";
    310   psMetadataAddVector( sortedTable, PS_LIST_TAIL, column_name, 0, NULL,
    311                        fill_F32_vector(psMetadataLookupVector(NULL, aTable, column_name), fillLength, sortedIndexVector, idets->n) );
    312   column_name = "EXT_WIDTH_MIN_ERR";
    313   psMetadataAddVector( sortedTable, PS_LIST_TAIL, column_name, 0, NULL,
    314                        fill_F32_vector(psMetadataLookupVector(NULL, aTable, column_name), fillLength, sortedIndexVector, idets->n) );
    315   column_name = "EXT_THETA_ERR";
    316   psMetadataAddVector( sortedTable, PS_LIST_TAIL, column_name, 0, NULL,
    317                        fill_F32_vector(psMetadataLookupVector(NULL, aTable, column_name), fillLength, sortedIndexVector, idets->n) );
    318   return sortedTable;
    319 }
    320 
    321 static psVector* fill_F32_vector(psVector* in, int length, long* indexes, int indexes_length) {
    322   if (in == NULL) {
    323     return NULL;
    324   }
    325   psVector* out = psVectorAlloc(length, in->type.type);
    326   int current_index = 0;
    327   for (int i = 0; i < indexes_length; i++) {
    328     while ( (current_index < indexes[i]) && (current_index < length) ) {
    329       out->data.F32[current_index] = NAN;
    330       current_index += 1;
    331     }
    332     out->data.F32[current_index] = in->data.F32[i];
    333     current_index += 1;
    334   }
    335   while ( current_index < length ) {
    336     out->data.F32[current_index] = NAN;
    337     current_index += 1;
    338   }
    339   return out;
    340 }
    341 
    342 static int compare (const void *a, const void *b) {
    343   return ( *((long*) a)-*((long*) b));
    344 }
    345 
    346 static long* sorted_idets_indexes(psVector *idets) {
    347   long* sortedIndexes = (long*) malloc(idets->n * sizeof(long));
    348   for (int i=0;i<idets->n;i++) {
    349     sortedIndexes[i] = idets->data.S64[i];
    350   }
    351   qsort(sortedIndexes, idets->n, sizeof(long), (void *) compare);
    352   return sortedIndexes;
    353 }
     281}
Note: See TracChangeset for help on using the changeset viewer.