IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Mar 22, 2011, 4:11:23 PM (15 years ago)
Author:
rhenders
Message:

Using new Fits class and Config class

File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/ippToPsps/src/StackBatch.c

    r30147 r31011  
    11#include "StackBatch.h"
    22#include "StackBatchEnums.h"
     3
     4#include "Fits.h"
    35
    46/**
    … …  
    1315
    1416/**
    15   Does the work.
    16   */
    17 static int run(StackBatch* this) {
    18 
    19     if (this->base.exitCode != PS_EXIT_SUCCESS) return this->base.exitCode;
    20 
    21     int status = 0;
    22     fitsfile *fitsIn;
    23 
    24     if (fits_open_file(&fitsIn, this->base.inputFiles[0], READONLY, &status)) {
    25 
    26         fits_report_error(stderr, status);
    27         return PS_EXIT_SYS_ERROR;
    28     }
    29 
    30 
    31     // get primary header and pull stuff out for later
    32     int nKeys = 0;
    33     fits_get_hdrspace(fitsIn, &nKeys, NULL, &status);
    34     float exposureTime; status=0; fits_read_key(fitsIn, TFLOAT, "EXPTIME", &exposureTime, NULL, &status);
    35     char filterType[20]; status=0; fits_read_key(fitsIn, TSTRING, "FILTER", filterType, NULL, &status);
    36 
    37 
     17  Creates the StackDetection table
     18  */
     19static bool createStackDetectionTable(
     20        StackBatch* this,
     21        Fits* fitsIn,
     22        int8_t* filterIDs,
     23        int8_t* surveyIDs,
     24        long* skycellIDs,
     25        float exposureTime
     26        ) {
     27
     28    char extensionName[25];
     29    sprintf(extensionName, "SkyChip.psf");
     30    long nDet = 0;
     31    if (!fitsIn->moveToBinaryTableAndCountRows(fitsIn, extensionName, &nDet)) return false;
     32
     33
     34    // allocate stuff
    3835    char** assocDate = (char**)calloc(this->MAXDETECT, sizeof(char**));
    3936    for (uint32_t i=0; i<this->MAXDETECT;i++) assocDate[i] = (char*)calloc(20,sizeof(char));
    40 
    41 
    42 
    43     // write StackMeta
    44     ippToPspsConfig_writeTable(this->base.config, fitsIn, this->base.fitsOut, 1, "StackMeta", true);
    45     fits_write_col(this->base.fitsOut, TLONG, STACKMETA_SKYCELLID, 1, 1, 1, &this->skycellId, &status);
    46 
    47     int8_t filterID = -1;
    48     if (!ippToPspsConfig_getFilterId(this->base.config, filterType, &filterID)) {
    49 
    50 //        this->base.exitCode = PS_EXIT_DATA_ERROR;
    51     //    return this->base.exitCode; TODO
    52     }
    53 
    5437    long* removeList = (long*)calloc(this->MAXDETECT, sizeof(long));
    5538    float* instMag = (float*)calloc(this->MAXDETECT, sizeof(float));
    56     int8_t* filterIDs = (int8_t*)calloc(this->MAXDETECT, sizeof(int8_t));
    57     int8_t* surveyIDs = (int8_t*)calloc(this->MAXDETECT, sizeof(int8_t));
    58     float floatnull = -999.0;
    59     int anynull = 0;
    60 
    61     fits_write_col(this->base.fitsOut, TBYTE, STACKMETA_FILTERID, 1, 1, 1, &filterID, &status);
    62 
    63     fits_write_col(this->base.fitsOut, TBYTE, STACKMETA_SURVEYID, 1, 1, 1, &this->base.surveyID, &status);
    64 
    65 
    66     // psf detections
    67     char extensionName[15];
    68     sprintf(extensionName, "Chip.psf");
    69     if (fits_movnam_hdu(fitsIn, BINARY_TBL, extensionName, 0, &status)) {
    70         psError(PS_ERR_IO, false, "Can't move to extension: %s\n", extensionName);
    71 
    72     }
    73     else {
    74 
    75         // some stuff is the same for all detections so we can populate here
    76         for (long s = 0; s<this->MAXDETECT; s++) {
    77 
    78             filterIDs[s] = filterID;
    79             surveyIDs[s] = this->base.surveyID;
    80             strcpy(assocDate[s], this->base.todaysDate);
     39    float* peakMag = (float*)calloc(this->MAXDETECT, sizeof(float));
     40    float* peakFlux = (float*)calloc(this->MAXDETECT, sizeof(float));
     41    int* flags1 = (int*)calloc(this->MAXDETECT, sizeof(int));
     42    int* flags2 = (int*)calloc(this->MAXDETECT, sizeof(int));
     43    uint64_t* infoFlags = (uint64_t*)calloc(this->MAXDETECT, sizeof(uint64_t));
     44
     45    bool peakFluxOk;
     46
     47    // some stuff is the same for all detections so we can populate here
     48    for (long s = 0; s<this->MAXDETECT; s++) {
     49
     50        // if running in test mode, don't use today's date
     51        if (this->base.testMode) strcpy(assocDate[s], "2010-01-01");
     52        else strcpy(assocDate[s], this->base.todaysDate);
     53    }
     54
     55    fitsIn->readColumnUsingName(fitsIn, TFLOAT, "PSF_INST_MAG", 1, 1, nDet, instMag);
     56    fitsIn->readColumnUsingName(fitsIn, TFLOAT, "PEAK_FLUX_AS_MAG", 1, 1, nDet, peakMag);
     57    fitsIn->readColumnUsingName(fitsIn, TLONG, "FLAGS", 1, 1, nDet, flags1);
     58    fitsIn->readColumnUsingName(fitsIn, TLONG, "FLAGS2", 1, 1, nDet, flags2);
     59
     60    printf("Looping through %ld psf detections\n", nDet);
     61    float mag;
     62    long unmatched = 0, totalDetections = 0, numOfDuplicates = 0, numInvalidFlux = 0, numDetectionsOut = 0;
     63
     64    for (long s=0; s<nDet; s++) {
     65
     66        // TODO implement this match in DVO
     67        if (1) {
     68
     69            //infoFlags[s] = ((uint64_t)flags1[s] << 32) | (uint64_t)flags2[s]; TODO implement after schema change to make infoFlag 64-bit
     70            infoFlags[s] = (uint64_t)flags1[s];
     71
     72            peakFluxOk = getFlux(exposureTime, peakMag[s], &peakFlux[s], 0.0, NULL);
     73
     74            mag = instMag[s];
     75            if (!peakFluxOk || !isfinite(mag) || mag < -998.0) {
     76
     77                removeList[numOfDuplicates+numInvalidFlux] = s+1;
     78                numInvalidFlux++;
     79            }
     80
     81            totalDetections++;
    8182        }
    82 
    83 
    84         long nDet = 0;
    85         if (fits_get_num_rows(fitsIn, &nDet, &status)) {
    86             fits_report_error(stderr, status);
     83        else {
     84
     85            unmatched++;
     86            continue;
    8787        }
    88 
    89         int instMagNum;
    90         status=0;fits_get_colnum(fitsIn, CASESEN, "PSF_INST_MAG", &instMagNum, &status);
    91         if (status) psError(PS_ERR_IO, false, "Unable to read col num for PSF_INST_MAG");
    92         fits_read_col(fitsIn, TFLOAT, instMagNum, 1, 1, nDet, &floatnull, instMag, &anynull, &status);
    93 
    94 
    95         printf("Looping through %ld psf detections\n", nDet);
    96         float mag;
    97         long unmatched = 0, totalDetections = 0, numOfDuplicates = 0, numInvalidFlux = 0, numDetectionsOut = 0;
    98 
    99         for (long s = 0; s<nDet; s++) {
    100 
    101             // TODO implement this match in DVO
    102             if (1) {
    103 
    104                 mag = instMag[s];
    105                 if (!isfinite(mag) || mag < -998.0) {
    106 
    107                     removeList[numOfDuplicates+numInvalidFlux] = s+1;
    108                     numInvalidFlux++;
    109                 }
    110 
    111                 totalDetections++;
    112             }
    113             else {
    114 
    115                 unmatched++;
    116                 continue;
    117             }
    118         }
    119 
    120         numDetectionsOut = totalDetections - numInvalidFlux;
    121 
    122         if (numDetectionsOut > 0) {
    123 
    124             ippToPspsConfig_writeTable(this->base.config, fitsIn, this->base.fitsOut, nDet, "StackDetection", false);
    125             fits_write_col(this->base.fitsOut, TLONG, STACKDETECTION_SKYCELLID, 1, 1, 1, &this->skycellId, &status);
    126             fits_write_col(this->base.fitsOut, TBYTE, STACKDETECTION_FILTERID, 1, 1, nDet, filterIDs, &status);
    127             fits_write_col(this->base.fitsOut, TBYTE, STACKDETECTION_SURVEYID, 1, 1, nDet, surveyIDs, &status);
    128             fits_write_col(this->base.fitsOut, TSTRING, STACKDETECTION_ASSOCDATE, 1, 1, nDet, assocDate, &status);
    129 
    130             if (numInvalidFlux) fits_delete_rowlist(this->base.fitsOut, removeList, numInvalidFlux, &status);
    131 
    132         }
    133         psLogMsg("ippToPsps", PS_LOG_INFO,
    134                 "+---------------+---------+----------+------------------+---------------+--------------+\n"
    135                 "|   Extension   | Rows in | Rows out | Missing from DVO | Duplicate IDs | Invalid Flux |\n"
    136                 "|  %12s |  %5ld  |   %5ld  |      %5ld       |    %5ld      |    %5ld     |\n",
    137                 extensionName, nDet, numDetectionsOut, unmatched, numOfDuplicates, numInvalidFlux);
    138 
    139     }
    140 
    141 
    142 
    143     // extended source
    144     sprintf(extensionName, "Chip.xsrc");
    145     if (fits_movnam_hdu(fitsIn, BINARY_TBL, extensionName, 0, &status)) {
    146         psError(PS_ERR_IO, false, "Can't move to extension: %s\n", extensionName);
    147 
    148     }
    149     else {
    150 
    151         int colNum;
    152         int type;
    153         long repeat;
    154 
    155         ippToPspsConfig_getFitsColumnMeta(
    156                 "PROF_FLUX",
    157                 &colNum,
    158                 &type,
    159                 &repeat,
    160                 fitsIn);
    161         printf("PROF_FILL = %d %d %ld\n", colNum, type, repeat);
    162         float* vector = calloc(repeat, sizeof(float));
    163 
    164         ippToPspsConfig_getColumnVector(
    165                 colNum,
    166                 1,
    167                 type,
    168                 repeat,
    169                 vector,
    170                 fitsIn);
    171 
    172 
    173 
    174         free(vector);
    175 
    176         long nDet = 0;
    177         if (fits_get_num_rows(fitsIn, &nDet, &status)) {
    178             fits_report_error(stderr, status);
    179         }
    180 
    181         printf("Looping through %ld extended source detections\n", nDet);
    182         for (long s = 0; s<nDet; s++) {
    183 
    184 
    185 
    186         }
    187 
    188         ippToPspsConfig_writeTable(this->base.config, fitsIn, this->base.fitsOut, nDet, "StackApFlx", false);
    189         fits_write_col(this->base.fitsOut, TBYTE, STACKAPFLX_FILTERID, 1, 1, nDet, filterIDs, &status);
    190         fits_write_col(this->base.fitsOut, TBYTE, STACKAPFLX_SURVEYID, 1, 1, nDet, surveyIDs, &status);
    191     }
     88    }
     89
     90    numDetectionsOut = totalDetections - numInvalidFlux;
     91    if (numDetectionsOut > 0) {
     92
     93        this->base.config->createAndPopulateTable(this->base.config, fitsIn, this->base.fitsOut, nDet, "StackDetection", false);
     94        this->base.fitsOut->writeColumn(this->base.fitsOut, TLONG, STACKDETECTION_SKYCELLID, 1, 1, nDet, skycellIDs);
     95        this->base.fitsOut->writeColumn(this->base.fitsOut, TBYTE, STACKDETECTION_FILTERID, 1, 1, nDet, filterIDs);
     96        this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKDETECTION_PEAKFLUX, 1, 1, nDet, peakFlux);
     97        this->base.fitsOut->writeColumn(this->base.fitsOut, TBYTE, STACKDETECTION_SURVEYID, 1, 1, nDet, surveyIDs);
     98        this->base.fitsOut->writeColumn(this->base.fitsOut, TSTRING, STACKDETECTION_ASSOCDATE, 1, 1, nDet, assocDate);
     99        this->base.fitsOut->writeColumn(this->base.fitsOut, TLONGLONG, STACKDETECTION_INFOFLAG, 1, 1, nDet, flags1);
     100
     101        if (numInvalidFlux) this->base.fitsOut->deleteRows(this->base.fitsOut, removeList, numInvalidFlux);
     102
     103    }
     104    psLogMsg("ippToPsps", PS_LOG_INFO,
     105            "+---------------+---------+----------+------------------+---------------+--------------+\n"
     106            "|   Extension   | Rows in | Rows out | Missing from DVO | Duplicate IDs | Invalid Flux |\n"
     107            "|  %12s |  %5ld  |   %5ld  |      %5ld       |    %5ld      |    %5ld     |\n",
     108            extensionName, nDet, numDetectionsOut, unmatched, numOfDuplicates, numInvalidFlux);
    192109
    193110
    … …  
    195112    free(removeList);
    196113    free(instMag);
     114    free(flags1);
     115    free(flags2);
     116    free(infoFlags);
     117    for (uint32_t i=0; i<this->MAXDETECT;i++) free(assocDate[i]);
     118    free(assocDate);
     119
     120    return true;
     121}
     122/**
     123  Structure to encapsulate block of flux values for StackApFlx table
     124  */
     125typedef struct {
     126
     127    float* flx;
     128    float* flxErr;
     129    float* flxStd;
     130    float* flxFill;
     131
     132} StackApFlxFluxes;
     133
     134/**
     135  Creates the StackApFlx table for extended source attributes
     136  */
     137static bool createStackApFlxTable(
     138        StackBatch* this,
     139        Fits* fitsIn,
     140        int8_t* filterIDs,
     141        int8_t* surveyIDs
     142        ) {
     143
     144    long nDet = 0;
     145    if (!fitsIn->moveToBinaryTableAndCountRows(fitsIn, "SkyChip.xrad", &nDet)) return false;
     146
     147
     148    int aperFluxColNum, aperFluxErrColNum, aperFluxStDevColNum, aperFluxFillColNum;
     149    int type; // all the same type
     150    long repeat;
     151
     152    fitsIn->getColumnMeta(fitsIn, "APER_FLUX", &aperFluxColNum, &type, &repeat);
     153    fitsIn->getColumnMeta(fitsIn, "APER_FLUX_ERR", &aperFluxErrColNum, &type, &repeat);
     154    fitsIn->getColumnMeta(fitsIn, "APER_FLUX_STDEV", &aperFluxStDevColNum, &type, &repeat);
     155    fitsIn->getColumnMeta(fitsIn, "APER_FILL", &aperFluxFillColNum, &type, &repeat);
     156    float* vector = calloc(repeat, sizeof(float));
     157
     158    // we store 10 different flux values
     159    StackApFlxFluxes* stackApFlxFluxes = (StackApFlxFluxes*)calloc(11, sizeof(StackApFlxFluxes));
     160    for (long i=0; i<11; i++) {
     161
     162        stackApFlxFluxes[i].flx = (float*)calloc(nDet, sizeof(float));
     163        stackApFlxFluxes[i].flxErr = (float*)calloc(nDet, sizeof(float));
     164        stackApFlxFluxes[i].flxStd = (float*)calloc(nDet, sizeof(float));
     165        stackApFlxFluxes[i].flxFill = (float*)calloc(nDet, sizeof(float));
     166    }
     167
     168    printf("Looping through %ld extended source detections\n", nDet);
     169    for (long s=0; s<nDet; s++) {
     170
     171        fitsIn->getColumnVector(fitsIn, aperFluxColNum,s, type, repeat, vector);
     172        for (int i=0; i<repeat; i++) stackApFlxFluxes[i].flx[s] = vector[i];
     173
     174        fitsIn->getColumnVector(fitsIn, aperFluxErrColNum,s, type, repeat, vector);
     175        for (int i=0; i<repeat; i++) stackApFlxFluxes[i].flxErr[s] = vector[i];
     176
     177        fitsIn->getColumnVector(fitsIn, aperFluxStDevColNum,s, type, repeat, vector);
     178        for (int i=0; i<repeat; i++) stackApFlxFluxes[i].flxStd[s] = vector[i];
     179
     180        fitsIn->getColumnVector(fitsIn, aperFluxFillColNum,s, type, repeat, vector);
     181        for (int i=0; i<repeat; i++) stackApFlxFluxes[i].flxFill[s] = vector[i];
     182
     183    }
     184
     185    free(vector);
     186
     187    //int status = 0;
     188    this->base.config->createAndPopulateTable(this->base.config, fitsIn, this->base.fitsOut, nDet, "StackApFlx", false);
     189    this->base.fitsOut->writeColumn(this->base.fitsOut, TBYTE, STACKAPFLX_FILTERID, 1, 1, nDet, filterIDs);
     190    this->base.fitsOut->writeColumn(this->base.fitsOut, TBYTE, STACKAPFLX_SURVEYID, 1, 1, nDet, surveyIDs);
     191
     192    // R1->r10 flux values
     193    // 1
     194    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR1, 1, 1, nDet, stackApFlxFluxes[0].flx);
     195    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR1ERR, 1, 1, nDet, stackApFlxFluxes[0].flxErr);
     196    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR1STD, 1, 1, nDet, stackApFlxFluxes[0].flxStd);
     197    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR1FILL, 1, 1, nDet, stackApFlxFluxes[0].flxFill);
     198    // 2
     199    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR2, 1, 1, nDet, stackApFlxFluxes[1].flx);
     200    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR2ERR, 1, 1, nDet, stackApFlxFluxes[1].flxErr);
     201    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR2STD, 1, 1, nDet, stackApFlxFluxes[1].flxStd);
     202    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR2FILL, 1, 1, nDet, stackApFlxFluxes[1].flxFill);
     203    // 3
     204    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR3, 1, 1, nDet, stackApFlxFluxes[2].flx);
     205    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR3ERR, 1, 1, nDet, stackApFlxFluxes[2].flxErr);
     206    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR3STD, 1, 1, nDet, stackApFlxFluxes[2].flxStd);
     207    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR3FILL, 1, 1, nDet, stackApFlxFluxes[2].flxFill);
     208    // 4
     209    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR4, 1, 1, nDet, stackApFlxFluxes[3].flx);
     210    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR4ERR, 1, 1, nDet, stackApFlxFluxes[3].flxErr);
     211    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR4STD, 1, 1, nDet, stackApFlxFluxes[3].flxStd);
     212    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR4FILL, 1, 1, nDet, stackApFlxFluxes[3].flxFill);
     213    // 5
     214    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR5, 1, 1, nDet, stackApFlxFluxes[4].flx);
     215    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR5ERR, 1, 1, nDet, stackApFlxFluxes[4].flxErr);
     216    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR5STD, 1, 1, nDet, stackApFlxFluxes[4].flxStd);
     217    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR5FILL, 1, 1, nDet, stackApFlxFluxes[4].flxFill);
     218    // 6
     219    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR6, 1, 1, nDet, stackApFlxFluxes[5].flx);
     220    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR6ERR, 1, 1, nDet, stackApFlxFluxes[5].flxErr);
     221    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR6STD, 1, 1, nDet, stackApFlxFluxes[5].flxStd);
     222    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR6FILL, 1, 1, nDet, stackApFlxFluxes[5].flxFill);
     223    // 7
     224    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR7, 1, 1, nDet, stackApFlxFluxes[6].flx);
     225    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR7ERR, 1, 1, nDet, stackApFlxFluxes[6].flxErr);
     226    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR7STD, 1, 1, nDet, stackApFlxFluxes[6].flxStd);
     227    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR7FILL, 1, 1, nDet, stackApFlxFluxes[6].flxFill);
     228    // 8
     229    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR8, 1, 1, nDet, stackApFlxFluxes[7].flx);
     230    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR8ERR, 1, 1, nDet, stackApFlxFluxes[7].flxErr);
     231    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR8STD, 1, 1, nDet, stackApFlxFluxes[7].flxStd);
     232    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR8FILL, 1, 1, nDet, stackApFlxFluxes[7].flxFill);
     233    // 9
     234    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR9, 1, 1, nDet, stackApFlxFluxes[8].flx);
     235    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR9ERR, 1, 1, nDet, stackApFlxFluxes[8].flxErr);
     236    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR9STD, 1, 1, nDet, stackApFlxFluxes[8].flxStd);
     237    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR9FILL, 1, 1, nDet, stackApFlxFluxes[8].flxFill);
     238    // 10
     239    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR10, 1, 1, nDet, stackApFlxFluxes[9].flx);
     240    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR10ERR, 1, 1, nDet, stackApFlxFluxes[9].flxErr);
     241    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR10STD, 1, 1, nDet, stackApFlxFluxes[9].flxStd);
     242    this->base.fitsOut->writeColumn(this->base.fitsOut, TFLOAT, STACKAPFLX_FLXR10FILL, 1, 1, nDet, stackApFlxFluxes[9].flxFill);
     243
     244    for (long i=0; i<11; i++) {
     245
     246        free(stackApFlxFluxes[i].flx);
     247        free(stackApFlxFluxes[i].flxErr);
     248        free(stackApFlxFluxes[i].flxStd);
     249        free(stackApFlxFluxes[i].flxFill);
     250    }
     251
     252    free(stackApFlxFluxes);
     253
     254    return true;
     255}
     256
     257/**
     258  Creates the StackModelFit table
     259  */
     260static bool createStackModelFitTable(
     261        StackBatch* this,
     262        Fits *fitsIn,
     263        int8_t* filterIDs,
     264        int8_t* surveyIDs) {
     265
     266    long nDet = 0;
     267    if (!fitsIn->moveToBinaryTableAndCountRows(fitsIn, "SkyChip.xfit", &nDet)) return false;
     268
     269    long* ippIdets = (long*)calloc(this->MAXDETECT, sizeof(long));
     270    char** modelTypes = (char**)calloc(this->MAXDETECT, sizeof(char**));
     271    for (uint32_t i=0; i<this->MAXDETECT;i++) modelTypes[i] = (char*)calloc(20,sizeof(char));
     272
     273    // read whole columns
     274    fitsIn->readColumnUsingName(fitsIn, TLONG, "IPP_IDET", 1, 1, nDet, ippIdets);
     275    fitsIn->readColumnUsingName(fitsIn, TSTRING, "MODEL_TYPE", 1, 1, nDet, modelTypes);
     276
     277    printf("Looping through %ld model fit rows\n", nDet);
     278    for (long i = 0; i<nDet; i++) {
     279
     280        printf("IPP_IDET %ld model %s\n", ippIdets[i], modelTypes[i]);
     281
     282    }
     283
     284    this->base.config->createAndPopulateTable(this->base.config, fitsIn, this->base.fitsOut, nDet, "StackModelFit", false);
     285    this->base.fitsOut->writeColumn(this->base.fitsOut, TBYTE, STACKMODELFIT_FILTERID, 1, 1, nDet, filterIDs);
     286    this->base.fitsOut->writeColumn(this->base.fitsOut, TBYTE, STACKMODELFIT_SURVEYID, 1, 1, nDet, surveyIDs);
     287
     288    // free up
     289    free(ippIdets);
     290    for (uint32_t i=0; i<this->MAXDETECT;i++) free(modelTypes[i]);
     291    free(modelTypes);
     292
     293    return true;
     294}
     295
     296/**
     297  Creates the StackToImage table
     298  */
     299static bool createStackToImageTable(
     300        StackBatch* this,
     301        Fits* fitsIn
     302        ) {
     303
     304    this->base.config->createAndPopulateTable(this->base.config, fitsIn, this->base.fitsOut, 3/*TODO*/, "StackToImage", false);
     305
     306    return true;
     307}
     308
     309/**
     310  Does the work. Writes all extensions to new FITS file.
     311  */
     312static int run(StackBatch* this) {
     313
     314    if (this->base.exitCode != PS_EXIT_SUCCESS) return this->base.exitCode;
     315
     316    // open input FITS file
     317    Fits* fitsIn = existing_Fits(this->base.inputFiles[0]);
     318    if (fitsIn->getFilePtr(fitsIn) == NULL)  return PS_EXIT_SYS_ERROR;
     319
     320    // get primary header and pull stuff out for later
     321    float exposureTime; fitsIn->getHeaderKeyValue(fitsIn, TFLOAT, "EXPTIME", &exposureTime);
     322    char filterType[20]; fitsIn->getHeaderKeyValue(fitsIn, TSTRING, "FILTER", filterType);
     323
     324    int8_t filterID = -1;
     325    if (!this->base.config->getFilterId(this->base.config, filterType, &filterID)) {
     326
     327        //        this->base.exitCode = PS_EXIT_DATA_ERROR;
     328        //    return this->base.exitCode; TODO
     329    }
     330
     331    exposureTime = 60.0;// TODO
     332    filterID = 3; // TODO
     333
     334    // write StackMeta
     335    this->base.config->createAndPopulateTable(this->base.config, fitsIn, this->base.fitsOut, 1, "StackMeta", true);
     336    this->base.fitsOut->writeColumn(this->base.fitsOut, TLONG, STACKMETA_SKYCELLID, 1, 1, 1, &this->skycellId);
     337    this->base.fitsOut->writeColumn(this->base.fitsOut, TBYTE, STACKMETA_FILTERID, 1, 1, 1, &filterID);
     338    this->base.fitsOut->writeColumn(this->base.fitsOut, TBYTE, STACKMETA_SURVEYID, 1, 1, 1, &this->base.surveyID);
     339
     340    // allocate stuff for other extensions
     341    int8_t* filterIDs = (int8_t*)calloc(this->MAXDETECT, sizeof(int8_t));
     342    int8_t* surveyIDs = (int8_t*)calloc(this->MAXDETECT, sizeof(int8_t));
     343    long* skycellIDs = (long*)calloc(this->MAXDETECT, sizeof(long));
     344
     345    // some stuff is the same for all detections so we can populate here
     346    for (long s = 0; s<this->MAXDETECT; s++) {
     347
     348        skycellIDs[s] = this->skycellId;
     349        filterIDs[s] = filterID;
     350        surveyIDs[s] = this->base.surveyID;
     351    }
     352
     353    // create all extensions
     354    createStackDetectionTable(this, fitsIn, filterIDs, surveyIDs, skycellIDs, exposureTime);
     355    createStackApFlxTable(this, fitsIn, filterIDs, surveyIDs);
     356    printf("dsdsdddds\n\n");
     357    createStackModelFitTable(this, fitsIn, filterIDs, surveyIDs);
     358    createStackToImageTable(this, fitsIn);
     359
     360    // free stuff up
    197361    free(filterIDs);
    198362    free(surveyIDs);
    199 
    200     status=0;
    201     if (fits_close_file(fitsIn, &status)) fits_report_error(stderr, status);
    202 
     363    free(skycellIDs);
     364    fitsIn->destroy(fitsIn);
    203365
    204366    return this->base.exitCode;
    205367}
    206 
    207368
    208369/**
    … …  
    240401    }
    241402
    242     this->base.parseArguments(&this->base, argc, argv);
     403    char fitsOutFile[40];
     404    sprintf(fitsOutFile, "%08d.FITS", this->skycellId);
     405    this->base.parseArguments(&this->base, argc, argv, "/stack", fitsOutFile);
    243406
    244407    if (
    … …  
    275438    this->MAXDETECT = 150000;
    276439
     440    // method pointers
    277441    this->print = print;
    278442    this->destroy = destroy;
    279443    this->base.run = run;
    280444
    281     if (!parseArguments(this, *argc, argv)) { return this; }
    282 
    283     strcat(this->base.configsDir, "/stack");
    284     sprintf(this->base.fitsOutFile, "%08d.FITS", this->skycellId);
    285 
    286     this->base.init(&this->base);
     445    parseArguments(this, *argc, argv);
    287446
    288447    return this;
    … …  
    321480}
    322481
    323 
    324 
Note: See TracChangeset for help on using the changeset viewer.