IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 33448


Ignore:
Timestamp:
Mar 9, 2012, 7:14:51 AM (14 years ago)
Author:
eugene
Message:

clean top-level function to move the former default mode into relastro_images; adapt to -parallel mode, including the use of MeasureTiny in relastro_images; add relastro_client, etc; make sure high-speed mode also gets a copy of the photcode table, skytable, and images; for high-speed, make an empty WHERE the default (and print the WHERE actually used); move star maps to a block in relastro_images, not load_catalog; drop the UserCatalog global (-catalog RA DEC now just defines a UserPatch value); dropped unused -dophot and -area options

Location:
branches/eam_branches/ipp-20111122/Ohana/src/relastro
Files:
1 added
25 edited

Legend:

Unmodified
Added
Removed

  • branches/eam_branches/ipp-20111122/Ohana/src/relastro/Makefile

    r32854 r33448  
    1 default: relastro
     1default: relastro relastro_client
    22help:
    33        @echo "make options: relastro (default)"
    … …  
    1818
    1919relastro: $(BIN)/relastro.$(ARCH)
    20 install: $(DESTBIN)/relastro $(DESTBIN)/testparallax
     20relastro_client: $(BIN)/relastro_client.$(ARCH)
     21
     22install: $(DESTBIN)/relastro $(DESTBIN)/relastro_client $(DESTBIN)/testparallax
    2123
    2224RELASTRO = \
    … …  
    5355$(SRC)/select_images.$(ARCH).o       \
    5456$(SRC)/relastro.$(ARCH).o            \
     57$(SRC)/relastro_images.$(ARCH).o    \
    5558$(SRC)/relastro_objects.$(ARCH).o    \
    5659$(SRC)/save_catalogs.$(ARCH).o       \
    57 $(SRC)/write_coords.$(ARCH).o        \
    5860$(SRC)/CoordOps.$(ARCH).o            \
    5961$(SRC)/FixProblemImages.$(ARCH).o    \
    … …  
    6466$(SRC)/relastro_merge_source.$(ARCH).o  \
    6567$(SRC)/resort_catalog.$(ARCH).o  \
     68$(SRC)/BrightCatalog.$(ARCH).o  \
    6669$(SRC)/relastroVisual.$(ARCH).o
     70
     71#$(SRC)/write_coords.$(ARCH).o        \
    6772
    6873$(RELASTRO): $(INC)/relastro.h $(KAPA_INCS)
    6974$(BIN)/relastro.$(ARCH): $(RELASTRO) $(KAPA_LIBS)
     75
     76RELASTRO_CLIENT = \
     77$(SRC)/ConfigInit.$(ARCH).o          \
     78$(SRC)/FitPM.$(ARCH).o               \
     79$(SRC)/FitPar.$(ARCH).o              \
     80$(SRC)/FitPMandPar.$(ARCH).o         \
     81$(SRC)/FitSimple.$(ARCH).o           \
     82$(SRC)/ImageOps.$(ARCH).o            \
     83$(SRC)/MosaicOps.$(ARCH).o           \
     84$(SRC)/ParFactor.$(ARCH).o           \
     85$(SRC)/SetSignals.$(ARCH).o          \
     86$(SRC)/Shutdown.$(ARCH).o            \
     87$(SRC)/UpdateObjects.$(ARCH).o       \
     88$(SRC)/UpdateObjectOffsets.$(ARCH).o \
     89$(SRC)/UpdateMeasures.$(ARCH).o      \
     90$(SRC)/GetAstromError.$(ARCH).o      \
     91$(SRC)/args.$(ARCH).o                \
     92$(SRC)/bcatalog.$(ARCH).o            \
     93$(SRC)/dvo_astrom_ops.$(ARCH).o      \
     94$(SRC)/fitpoly.$(ARCH).o             \
     95$(SRC)/initialize.$(ARCH).o          \
     96$(SRC)/liststats.$(ARCH).o           \
     97$(SRC)/load_catalogs.$(ARCH).o       \
     98$(SRC)/load_images.$(ARCH).o         \
     99$(SRC)/mkpolyterm.$(ARCH).o          \
     100$(SRC)/select_images.$(ARCH).o       \
     101$(SRC)/relastro_client.$(ARCH).o        \
     102$(SRC)/relastro_objects.$(ARCH).o    \
     103$(SRC)/save_catalogs.$(ARCH).o       \
     104$(SRC)/CoordOps.$(ARCH).o            \
     105$(SRC)/high_speed_catalogs.$(ARCH).o  \
     106$(SRC)/high_speed_objects.$(ARCH).o  \
     107$(SRC)/high_speed_utils.$(ARCH).o  \
     108$(SRC)/relastro_merge_source.$(ARCH).o  \
     109$(SRC)/resort_catalog.$(ARCH).o  \
     110$(SRC)/plot_scatter.$(ARCH).o        \
     111$(SRC)/plotstuff.$(ARCH).o           \
     112$(SRC)/relastroVisual.$(ARCH).o \
     113$(SRC)/BrightCatalog.$(ARCH).o
     114
     115#$(SRC)/FitChip.$(ARCH).o             \
     116#$(SRC)/FitMosaic.$(ARCH).o           \
     117#$(SRC)/UpdateChips.$(ARCH).o         \
     118#$(SRC)/UpdateMosaic.$(ARCH).o        \
     119#$(SRC)/UpdateSimple.$(ARCH).o        \
     120#$(SRC)/write_coords.$(ARCH).o        \
     121#$(SRC)/FixProblemImages.$(ARCH).o    \
     122#$(SRC)/StarMaps.$(ARCH).o    \
     123
     124$(RELASTRO_CLIENT): $(INC)/relastro.h $(KAPA_INCS)
     125$(BIN)/relastro_client.$(ARCH): $(RELASTRO_CLIENT) $(KAPA_LIBS)
    70126
    71127testparallax: $(BIN)/testparallax.$(ARCH)

  • branches/eam_branches/ipp-20111122/Ohana/src/relastro/include/relastro.h

    r32740 r33448  
    1515typedef enum {FIT_NONE, FIT_AVERAGE, FIT_PM_ONLY, FIT_PAR_ONLY, FIT_PM_AND_PAR} FitMode;
    1616
    17 typedef enum {TARGET_NONE, TARGET_OBJECTS, TARGET_SIMPLE, TARGET_CHIPS, TARGET_MOSAICS, TARGET_HIGH_SPEED, TARGET_MERGE_SOURCE} FitTarget;
     17typedef enum {TARGET_NONE, TARGET_SIMPLE, TARGET_CHIPS, TARGET_MOSAICS, TARGET_HIGH_SPEED, TARGET_MERGE_SOURCE, TARGET_UPDATE_OBJECTS, TARGET_UPDATE_OFFSETS, TARGET_LOAD_OBJECTS} FitTarget;
     18
     19typedef struct {
     20  Average     *average;       // array of (minimal) average data
     21  MeasureTiny *measure;       // array of (minimal) measure data
     22  SecFilt     *secfilt;       // array of secfilt data (matched to average by Nsecfilt)
     23  off_t       Naverage;
     24  off_t       Nmeasure;
     25} BrightCatalog;
     26
     27typedef struct {
     28  Catalog *catalog;           // array of catalogs generated
     29  int NCATALOG;               // number of catalogs allocated
     30  int Ncatalog;               // number of catalogs generated
     31  int Nsecfilt;               // number of catalogs generated
     32  off_t *NAVERAGE;            // allocated Averages per catalog
     33  off_t *NMEASURE;            // allocated Measures per catalog
     34  int   *index;               // lookup table catID -> catalog[i]
     35  int   *catIDs;              // lookup table catID <- catalog[i]
     36  int    maxID;               // max catID value to date
     37} CatalogSplitter;
    1838
    1939typedef struct {
    … …  
    86106
    87107/* global variables set in parameter file */
    88 char   ImageCat[256];
    89 char   GSCFILE[256];
    90 char   CATDIR[256];
     108# define MAX_PATH_LENGTH 1024
     109char   ImageCat[MAX_PATH_LENGTH];
     110char   GSCFILE[MAX_PATH_LENGTH];
     111char   CATDIR[MAX_PATH_LENGTH];
    91112char   *HIGH_SPEED_DIR;
    92113char   CATMODE[16];    /* raw, mef, split, mysql */
    93114char   CATFORMAT[16];  /* internal, elixir, loneos, panstarrs */
    94 char   SKY_TABLE[256];
     115char   SKY_TABLE[MAX_PATH_LENGTH];
    95116int    SKY_DEPTH;  /** XXX EAM : depth of catalog tables, fix usage */
     117
     118int          HOST_ID;
     119char        *HOSTDIR;
     120char        *BCATALOG;
    96121
    97122unsigned int OBJ_ID_SRC;
    … …  
    109134
    110135double RADIUS; // match radius for high-speed objects
     136
     137int    PARALLEL;
     138int    PARALLEL_MANUAL;
     139int    PARALLEL_SERIAL;
    111140
    112141int    VERBOSE;
    … …  
    129158double MaxDensityValue;
    130159
    131 int UserCatalog;
    132 double UserCatalogRA, UserCatalogDEC;
    133 
    134160char          *PHOTCODE_KEEP_LIST, *PHOTCODE_SKIP_LIST;
    135161int           NphotcodesKeep,      NphotcodesSkip;
    … …  
    142168SkyRegionSelection SELECTION;
    143169
    144 int AreaSelect;
    145 double AreaXmin, AreaXmax, AreaYmin, AreaYmax;
    146 
    147 int ImagSelect, ImagMin, ImagMax;
     170int ImagSelect;
     171double ImagMin, ImagMax;
    148172
    149173double  PlotMmin, PlotMmax, PlotdMmin, PlotdMmax;
    150174
    151 int DophotSelect, DophotValue;
    152175int PhotFlagSelect, PhotFlagPoor, PhotFlagBad;
    153176
    … …  
    181204char         *GetPhotnamebyCode   PROTO((PhotCodeData *photcodes, int code));
    182205void          InterpolateGrid     PROTO((float *buffer, int Nx, int Ny, Coords *ccd, Coords *gcoords));
    183 off_t          *SelectRefMosaic     PROTO((Mosaic **refmosaic, off_t *Nimage));
     206off_t        *SelectRefMosaic     PROTO((Mosaic **refmosaic, off_t *Nimage));
    184207int           args                PROTO((int argc, char **argv));
     208int           args_client         PROTO((int argc, char **argv));
    185209int           bcatalog            PROTO((Catalog *subcatalog, Catalog *catalog));
    186210void          clean_images        PROTO((void));
    … …  
    226250void          initMrel            PROTO((Catalog *catalog, int Ncatalog));
    227251void          initialize          PROTO((int argc, char **argv));
     252void          initialize_client   PROTO((int argc, char **argv));
    228253void          initstats           PROTO((char *mode));
    229254int           liststats           PROTO((double *value, double *dvalue, int N, StatType *stats));
    230 Catalog      *load_catalogs       PROTO((SkyList *skylist, int *Ncatalog, int subselect));
    231 SkyList      *load_images         PROTO((FITS_DB *db, SkyRegion *region));
     255Catalog      *load_catalogs       PROTO((SkyList *skylist, int *Ncatalog, int subselect, int hostID, char *hostpath));
     256int           load_images         PROTO((FITS_DB *db, SkyList *skylist));
    232257Image        *select_images       PROTO((SkyList *skylist, Image *timage, off_t Ntimage, off_t **LineNumber, off_t *Nimage));
    233258
    … …  
    314339void FlagOutliers(Catalog *catalog);
    315340int MeasFilterTest(Measure *measure, int applySigmaLim);
     341int MeasFilterTestTiny(MeasureTiny *measure, int applySigmaLim);
    316342
    317343int sun_ecliptic (double jd, double *lambda, double *beta, double *epsilon, double *Radius);
    … …  
    333359Mosaic *getMosaicForImage (off_t im);
    334360
    335 double getMeanR (Measure *measure, Average *average, SecFilt *secfilt);
    336 double getMeanD (Measure *measure, Average *average, SecFilt *secfilt);
    337 int setMeanR (double ra_fit, Measure *measure, Average *average, SecFilt *secfilt);
    338 int setMeanD (double dec_fit, Measure *measure, Average *average, SecFilt *secfilt);
     361double getMeanR (MeasureTiny *measure, Average *average, SecFilt *secfilt);
     362double getMeanD (MeasureTiny *measure, Average *average, SecFilt *secfilt);
     363int setMeanR (double ra_fit, MeasureTiny *measure, Average *average, SecFilt *secfilt);
     364int setMeanD (double dec_fit, MeasureTiny *measure, Average *average, SecFilt *secfilt);
     365int setMeanR_Big (double ra_fit, Measure *measure, Average *average, SecFilt *secfilt);
     366int setMeanD_Big (double dec_fit, Measure *measure, Average *average, SecFilt *secfilt);
    339367
    340368float GetAstromError (Measure *measure, int mode);
    341 int relastro_objects (void);
    342 int UpdateObjectOffsets (SkyList *skylist);
     369float GetAstromErrorTiny (MeasureTiny *measure, int mode);
     370int relastro_objects (SkyList *skylist, int hostID, char *hostpath);
     371int relastro_images (SkyList *skylist);
     372int UpdateObjectOffsets (SkyList *skylist, int hostID, char *hostpath);
    343373
    344374int relastroVisualPlotChipFit(StarData *raw, StarData *ref, double dRmax, int numObj);
    … …  
    361391void resetImageRaw (Catalog *catalog, int Ncatalog, off_t im);
    362392
    363 int high_speed_catalogs ();
     393int high_speed_catalogs (SkyTable *sky, SkyList *skylist, int hostID, char *hostpath);
    364394int high_speed_objects (SkyRegion *region, Catalog *catalog);
    365395int MeasMatchesPhotcode(Measure *measure, PhotCode **photcodeSet, int Nset);
    … …  
    369399int createStarMapPoints();
    370400int checkStarMap(int N);
     401int createStarMap (Catalog *catalog, int Ncatalog);
    371402
    372403int GetScatterRawRef(float *dLsig, float *dMsig, float *dRsig, int *nKeep, StarData *raw, StarData *ref, int Nstars, float SigmaLimit);
    … …  
    378409void setupAreaSelection(SkyRegion *region);
    379410
    380 int relastro_merge_source ();
     411int relastro_merge_source (SkyTable *sky);
    381412void resort_catalog (Catalog *catalog);
     413
     414BrightCatalog *BrightCatalogLoad(char *filename);
     415int BrightCatalogSave(char *filename, BrightCatalog *catalog);
     416BrightCatalog *BrightCatalogMerge (Catalog *catalog, int Ncatalog);
     417CatalogSplitter *BrightCatalogSplitInit (int Nsecfilt);
     418int BrightCatalogSplit (CatalogSplitter *catalogs, BrightCatalog *bcatalog);
     419int BrightCatalogSplitFree (CatalogSplitter *catalogs);

  • branches/eam_branches/ipp-20111122/Ohana/src/relastro/src/BrightCatalog.c

    r33396 r33448  
    11# include "relastro.h"
    22
    3 // relphot_client reads from the local catalogs and generates a set of bright, subset catalogs
     3// NOTE: this file is nearly identical to relphot/src/BrightCatalog.c, but here we use Average, not AverageTiny
     4
     5// relastro_client -load reads from the local catalogs and generates a set of bright, subset catalogs
    46// we are going to save these as a single big FITS file, with only 3 extensions:
    57
    68// measure, average, secfilt.  these only need to contain the fields relevant to the
    7 // MeasureTiny, AverageTiny, and Secfilt tables, but with cat_id appended to each row so we can
     9// MeasureTiny, Average, and Secfilt tables, but with cat_id appended to each row so we can
    810// reassign correctly on read
    911
    … …  
    131133    // read the fits table bytes
    132134    if (!gfits_fread_ftable_data (f, &ftable, FALSE)) return (NULL);
    133  
    134     // need to create and assign to flat-field correction
    135     GET_COLUMN(R,             "RA",        double);
    136     GET_COLUMN(D,             "DEC",       double);
    137     GET_COLUMN(Nmeasure,      "NMEAS",       int);
    138     GET_COLUMN(measureOffset, "MEASURE_OFF", int);
    139     GET_COLUMN(flags,         "FLAGS",       int);
    140     GET_COLUMN(catID,         "CAT_ID",      int);
     135
     136    GET_COLUMN(R,              "RA",          double);
     137    GET_COLUMN(D,              "DEC",         double);
     138    GET_COLUMN(dR,             "RA_ERR",      float);
     139    GET_COLUMN(dD,             "DEC_ERR",     float);
     140    GET_COLUMN(uR,             "U_RA",        float);
     141    GET_COLUMN(uD,             "U_DEC",       float);
     142    GET_COLUMN(duR,            "V_RA_ERR",    float);
     143    GET_COLUMN(duD,            "V_DEC_ERR",   float);
     144    GET_COLUMN(P,              "PAR",         float);
     145    GET_COLUMN(dP,             "PAR_ERR",     float);
     146    GET_COLUMN(ChiSqAve,       "CHISQ_POS",   float);
     147    GET_COLUMN(ChiSqPM,        "CHISQ_PM",    float);
     148    GET_COLUMN(ChiSqPar,       "CHISQ_PAP",   float);
     149    GET_COLUMN(Tmean,          "MEAN_EPOCH",  int);
     150    GET_COLUMN(Trange,         "TIME_RANGE",  int);
     151    GET_COLUMN(Xp,             "SIGMA_POS",   float);
     152    GET_COLUMN(Npos,           "NUMBER_POS",  short);
     153    GET_COLUMN(Nmeasure,       "NMEASURE",    short);
     154    GET_COLUMN(Nmissing,       "NMISSING",    short);
     155    GET_COLUMN(Nextend,        "NEXTEND",     short);
     156    GET_COLUMN(measureOffset,  "OFF_MEASURE", int);
     157    GET_COLUMN(missingOffset,  "OFF_MISSING", int);
     158    GET_COLUMN(extendOffset,   "OFF_EXTEND",  int);
     159    GET_COLUMN(flags,          "FLAGS",       int);
     160    GET_COLUMN(photFlagsUpper, "PHOTFLAGS_U", int);
     161    GET_COLUMN(photFlagsLower, "PHOTFLAGS_L", int);
     162    GET_COLUMN(objID,          "OBJ_ID",      int);
     163    GET_COLUMN(catID,          "CAT_ID",      int);
     164    GET_COLUMN(extID,          "EXT_ID",      int64_t);
    141165    // XXX free the fits table data here
    142166
    143     AverageTiny *average = NULL;
    144     ALLOCATE (average, AverageTiny, Nrow);
     167    Average *average = NULL;
     168    ALLOCATE (average, Average, Nrow);
    145169    for (i = 0; i < Nrow; i++) {
    146       average[i].R              = R[i];
    147       average[i].D              = D[i];
    148       average[i].Nmeasure       = Nmeasure[i];
    149       average[i].measureOffset  = measureOffset[i];
    150       average[i].flags          = flags[i];
    151       average[i].catID          = catID[i];
     170      average[i].R               = R[i]               ;
     171      average[i].D               = D[i]               ;
     172      average[i].dR              = dR[i]              ;
     173      average[i].dD              = dD[i]              ;
     174      average[i].uR              = uR[i]              ;
     175      average[i].uD              = uD[i]              ;
     176      average[i].duR             = duR[i]             ;
     177      average[i].duD             = duD[i]             ;
     178      average[i].P               = P[i]               ;
     179      average[i].dP              = dP[i]              ;
     180      average[i].ChiSqAve        = ChiSqAve[i]        ;
     181      average[i].ChiSqPM         = ChiSqPM[i]         ;
     182      average[i].ChiSqPar        = ChiSqPar[i]        ;
     183      average[i].Tmean           = Tmean[i]           ;
     184      average[i].Trange          = Trange[i]          ;
     185      average[i].Xp              = Xp[i]              ;
     186      average[i].Npos            = Npos[i]            ;
     187      average[i].Nmeasure        = Nmeasure[i]        ;
     188      average[i].Nmissing        = Nmissing[i]        ;
     189      average[i].Nextend         = Nextend[i]         ;
     190      average[i].measureOffset   = measureOffset[i]   ;
     191      average[i].missingOffset   = missingOffset[i]   ;
     192      average[i].extendOffset    = extendOffset[i]    ; 
     193      average[i].flags           = flags[i]           ;
     194      average[i].photFlagsUpper  = photFlagsUpper[i]  ;
     195      average[i].photFlagsLower  = photFlagsLower[i]  ;
     196      average[i].objID           = objID[i]           ;
     197      average[i].catID           = catID[i]           ;
     198      average[i].extID           = extID[i]           ;
    152199    }
    153200    fprintf (stderr, "loaded data for %lld averages\n", (long long) Nrow);
    154201
    155     free (R             );
    156     free (D             );
    157     free (Nmeasure      );
    158     free (measureOffset );
    159     free (flags         );
    160     free (catID         );
     202    free (R);
     203    free (D);
     204    free (dR);
     205    free (dD);
     206    free (uR);
     207    free (uD);
     208    free (duR);
     209    free (duD);
     210    free (P);
     211    free (dP);
     212    free (ChiSqAve);
     213    free (ChiSqPM);
     214    free (ChiSqPar);
     215    free (Tmean);
     216    free (Trange);
     217    free (Xp);
     218    free (Npos);
     219    free (Nmeasure);
     220    free (Nmissing);
     221    free (Nextend);
     222    free (measureOffset);
     223    free (missingOffset);
     224    free (extendOffset);
     225    free (flags);
     226    free (photFlagsUpper);
     227    free (photFlagsLower);
     228    free (objID);
     229    free (catID);
     230    free (extID);
    161231
    162232    catalog->average = average;
    … …  
    257327    gfits_define_bintable_column (&theader, "J", "IMAGE_ID", "image",                       NULL,    1.0, 1.0*0x8000);
    258328    gfits_define_bintable_column (&theader, "J", "DB_FLAGS", "flags",                       NULL,    1.0, 1.0*0x8000);
     329    gfits_define_bintable_column (&theader, "J", "PHOT_FLAGS", "photflags",                 NULL,    1.0, 1.0*0x8000);
    259330    gfits_define_bintable_column (&theader, "J", "CAT_ID",   "catalog",                     NULL,    1.0, 1.0*0x8000);
    260331    gfits_define_bintable_column (&theader, "I", "PHOTCODE", "photcode",                    NULL,    1.0, 1.0*0x80);
    … …  
    343414  }
    344415
    345   /*** AverageTiny ***/
     416  /*** Average ***/
    346417  {
    347     gfits_create_table_header (&theader, "BINTABLE", "AVERAGE_TINY");
    348 
    349     gfits_define_bintable_column (&theader, "D", "RA",          "ra (J2000)",            "degree", 1.0, 0.0);
    350     gfits_define_bintable_column (&theader, "D", "DEC",         "dec (J2000)",           "degree", 1.0, 0.0);
    351     gfits_define_bintable_column (&theader, "J", "NMEAS",       "number of measures",     NULL,    1.0, 0.0);
    352     gfits_define_bintable_column (&theader, "J", "MEASURE_OFF", "index to measurements",  NULL,    1.0, 0.0);
    353     gfits_define_bintable_column (&theader, "J", "FLAGS",       "flags",                  NULL,    1.0, 0.0);
    354     gfits_define_bintable_column (&theader, "J", "CAT_ID",      "catalog ref",            NULL,    1.0, 0.0);
     418    gfits_create_table_header (&theader, "BINTABLE", "AVERAGE");
     419
     420    gfits_define_bintable_column (&theader, "D", "RA",          "RA",                                                     "", 1.0, 0.0);
     421    gfits_define_bintable_column (&theader, "D", "DEC",         "DEC",                                                    "", 1.0, 0.0);
     422    gfits_define_bintable_column (&theader, "E", "RA_ERR",      "RA error",                                               "", 1.0, 0.0);
     423    gfits_define_bintable_column (&theader, "E", "DEC_ERR",     "DEC error",                                              "", 1.0, 0.0);
     424    gfits_define_bintable_column (&theader, "E", "U_RA",        "RA*cos(D) proper-motion",                                "", 1.0, 0.0);
     425    gfits_define_bintable_column (&theader, "E", "U_DEC",       "DEC proper-motion",                                      "", 1.0, 0.0);
     426    gfits_define_bintable_column (&theader, "E", "V_RA_ERR",    "RA*cos(D) p-m error",                                    "", 1.0, 0.0);
     427    gfits_define_bintable_column (&theader, "E", "V_DEC_ERR",   "DEC p-m error",                                          "", 1.0, 0.0);
     428    gfits_define_bintable_column (&theader, "E", "PAR",         "parallax",                                               "", 1.0, 0.0);
     429    gfits_define_bintable_column (&theader, "E", "PAR_ERR",     "parallax error",                                         "", 1.0, 0.0);
     430    gfits_define_bintable_column (&theader, "E", "CHISQ_POS",   "astrometry analysis chisq",                              "", 1.0, 0.0);
     431    gfits_define_bintable_column (&theader, "E", "CHISQ_PM",    "astrometry analysis chisq",                              "", 1.0, 0.0);
     432    gfits_define_bintable_column (&theader, "E", "CHISQ_PAP",   "astrometry analysis chisq",                              "", 1.0, 0.0);
     433    gfits_define_bintable_column (&theader, "J", "MEAN_EPOCH",  "mean epoch (PM-PAR ref)",                                "", 1.0, 1.0*0x8000);
     434    gfits_define_bintable_column (&theader, "J", "TIME_RANGE",  "mean epoch (PM-PAR ref)",                                "", 1.0, 1.0*0x8000);
     435    gfits_define_bintable_column (&theader, "E", "SIGMA_POS",   "position scatter",                                       "", 1.0, 0.0);   
     436    gfits_define_bintable_column (&theader, "I", "NUMBER_POS",  "number of detections used for astrometry",               "", 1.0, 1.0*0x80);
     437    gfits_define_bintable_column (&theader, "I", "NMEASURE",    "number of psf measurements",                             "", 1.0, 1.0*0x80);
     438    gfits_define_bintable_column (&theader, "I", "NMISSING",    "number of missings",                                     "", 1.0, 1.0*0x80);
     439    gfits_define_bintable_column (&theader, "I", "NEXTEND",     "number of extended measurements",                        "", 1.0, 1.0*0x80);
     440    gfits_define_bintable_column (&theader, "J", "OFF_MEASURE", "offset to first psf measurement",                        "", 1.0, 1.0*0x8000);
     441    gfits_define_bintable_column (&theader, "J", "OFF_MISSING", "offset to first missing obs",                            "", 1.0, 1.0*0x8000);
     442    gfits_define_bintable_column (&theader, "J", "OFF_EXTEND",  "offset to first extended measurement",                   "", 1.0, 1.0*0x8000);
     443    gfits_define_bintable_column (&theader, "J", "FLAGS",       "average object flags (star; ghost; etc)",                "", 1.0, 1.0*0x8000);
     444    gfits_define_bintable_column (&theader, "J", "PHOTFLAGS_U", "upper bit of 2 bit summary of per-measure photflags",    "", 1.0, 1.0*0x8000);
     445    gfits_define_bintable_column (&theader, "J", "PHOTFLAGS_L", "lower bit of 2 bit summary of per-measure photflags",    "", 1.0, 1.0*0x8000);
     446    gfits_define_bintable_column (&theader, "J", "OBJ_ID",      "unique ID for object in table",                          "", 1.0, 1.0*0x8000);
     447    gfits_define_bintable_column (&theader, "J", "CAT_ID",      "unique ID for table in which object was first realized", "", 1.0, 1.0*0x8000);
     448    gfits_define_bintable_column (&theader, "K", "EXT_ID",      "external ID for object (eg PSPS objID)",                 "", 1.0, 1.0*0x80000000);
    355449
    356450    // generate the output array that carries the data
    … …  
    358452
    359453    // create intermediate storage arrays
    360     double *R             ; ALLOCATE (R,             double, catalog->Naverage);
    361     double *D             ; ALLOCATE (D,             double, catalog->Naverage);
    362     int   *Nmeasure       ; ALLOCATE (Nmeasure,      int,    catalog->Naverage);
    363     int   *measureOffset  ; ALLOCATE (measureOffset, int,    catalog->Naverage);
    364     int   *flags          ; ALLOCATE (flags,         int,    catalog->Naverage);
    365     int   *catID          ; ALLOCATE (catID,         int,    catalog->Naverage);
     454    double   *R             ; ALLOCATE (R             , double  , catalog->Naverage);
     455    double   *D             ; ALLOCATE (D             , double  , catalog->Naverage);
     456    float    *dR            ; ALLOCATE (dR            , float   , catalog->Naverage);
     457    float    *dD            ; ALLOCATE (dD            , float   , catalog->Naverage);
     458    float    *uR            ; ALLOCATE (uR            , float   , catalog->Naverage);
     459    float    *uD            ; ALLOCATE (uD            , float   , catalog->Naverage);
     460    float    *duR           ; ALLOCATE (duR           , float   , catalog->Naverage);
     461    float    *duD           ; ALLOCATE (duD           , float   , catalog->Naverage);
     462    float    *P             ; ALLOCATE (P             , float   , catalog->Naverage);
     463    float    *dP            ; ALLOCATE (dP            , float   , catalog->Naverage);
     464    float    *ChiSqAve      ; ALLOCATE (ChiSqAve      , float   , catalog->Naverage);
     465    float    *ChiSqPM       ; ALLOCATE (ChiSqPM       , float   , catalog->Naverage);
     466    float    *ChiSqPar      ; ALLOCATE (ChiSqPar      , float   , catalog->Naverage);
     467    int      *Tmean         ; ALLOCATE (Tmean         , int     , catalog->Naverage);
     468    int      *Trange        ; ALLOCATE (Trange        , int     , catalog->Naverage);
     469    float    *Xp            ; ALLOCATE (Xp            , float   , catalog->Naverage);
     470    short    *Npos          ; ALLOCATE (Npos          , short   , catalog->Naverage);
     471    short    *Nmeasure      ; ALLOCATE (Nmeasure      , short   , catalog->Naverage);
     472    short    *Nmissing      ; ALLOCATE (Nmissing      , short   , catalog->Naverage);
     473    short    *Nextend       ; ALLOCATE (Nextend       , short   , catalog->Naverage);
     474    int      *measureOffset ; ALLOCATE (measureOffset , int     , catalog->Naverage);
     475    int      *missingOffset ; ALLOCATE (missingOffset , int     , catalog->Naverage);
     476    int      *extendOffset  ; ALLOCATE (extendOffset  , int     , catalog->Naverage);
     477    int      *flags         ; ALLOCATE (flags         , int     , catalog->Naverage);
     478    int      *photFlagsUpper; ALLOCATE (photFlagsUpper, int     , catalog->Naverage);
     479    int      *photFlagsLower; ALLOCATE (photFlagsLower, int     , catalog->Naverage);
     480    int      *objID         ; ALLOCATE (objID         , int     , catalog->Naverage);
     481    int      *catID         ; ALLOCATE (catID         , int     , catalog->Naverage);
     482    uint64_t *extID         ; ALLOCATE (extID         , uint64_t, catalog->Naverage);
    366483
    367484    // assign the storage arrays
    368     AverageTiny *average = catalog->average;
     485    Average *average = catalog->average;
    369486    for (i = 0; i < catalog->Naverage; i++) {
    370       R[i]              = average[i].R       ;
    371       D[i]              = average[i].D       ;
    372       Nmeasure[i]       = average[i].Nmeasure;
    373       measureOffset[i]  = average[i].measureOffset;
    374       flags[i]          = average[i].flags;
    375       catID[i]          = average[i].catID;
     487      R[i]               = average[i].R               ;
     488      D[i]               = average[i].D               ;
     489      dR[i]              = average[i].dR              ;
     490      dD[i]              = average[i].dD              ;
     491      uR[i]              = average[i].uR              ;
     492      uD[i]              = average[i].uD              ;
     493      duR[i]             = average[i].duR             ;
     494      duD[i]             = average[i].duD             ;
     495      P[i]               = average[i].P               ;
     496      dP[i]              = average[i].dP              ;
     497      ChiSqAve[i]        = average[i].ChiSqAve        ;
     498      ChiSqPM[i]         = average[i].ChiSqPM         ;
     499      ChiSqPar[i]        = average[i].ChiSqPar        ;
     500      Tmean[i]           = average[i].Tmean           ;
     501      Trange[i]          = average[i].Trange          ;
     502      Xp[i]              = average[i].Xp              ;
     503      Npos[i]            = average[i].Npos            ;
     504      Nmeasure[i]        = average[i].Nmeasure        ;
     505      Nmissing[i]        = average[i].Nmissing        ;
     506      Nextend[i]         = average[i].Nextend         ;
     507      measureOffset[i]   = average[i].measureOffset   ;
     508      missingOffset[i]   = average[i].missingOffset   ;
     509      extendOffset[i]    = average[i].extendOffset    ; 
     510      flags[i]           = average[i].flags           ;
     511      photFlagsUpper[i]  = average[i].photFlagsUpper  ;
     512      photFlagsLower[i]  = average[i].photFlagsLower  ;
     513      objID[i]           = average[i].objID           ;
     514      catID[i]           = average[i].catID           ;
     515      extID[i]           = average[i].extID           ;
    376516    }
    377517
    378518    // add the columns to the output array
    379     gfits_set_bintable_column (&theader, &ftable, "RA",          R,             catalog->Naverage);
    380     gfits_set_bintable_column (&theader, &ftable, "DEC",         D,             catalog->Naverage);
    381     gfits_set_bintable_column (&theader, &ftable, "NMEAS",       Nmeasure,      catalog->Naverage);
    382     gfits_set_bintable_column (&theader, &ftable, "MEASURE_OFF", measureOffset, catalog->Naverage);
    383     gfits_set_bintable_column (&theader, &ftable, "FLAGS",       flags,         catalog->Naverage);
    384     gfits_set_bintable_column (&theader, &ftable, "CAT_ID",      catID,         catalog->Naverage);
    385 
    386     free (R             );
    387     free (D             );
    388     free (Nmeasure      );
    389     free (measureOffset );
    390     free (flags         );
    391     free (catID         );
     519    gfits_set_bintable_column (&theader, &ftable, "RA",          R,               catalog->Naverage);
     520    gfits_set_bintable_column (&theader, &ftable, "DEC",         D,               catalog->Naverage);
     521    gfits_set_bintable_column (&theader, &ftable, "RA_ERR",      dR,              catalog->Naverage);
     522    gfits_set_bintable_column (&theader, &ftable, "DEC_ERR",     dD,              catalog->Naverage);
     523    gfits_set_bintable_column (&theader, &ftable, "U_RA",        uR,              catalog->Naverage);
     524    gfits_set_bintable_column (&theader, &ftable, "U_DEC",       uD,              catalog->Naverage);
     525    gfits_set_bintable_column (&theader, &ftable, "V_RA_ERR",    duR,             catalog->Naverage);
     526    gfits_set_bintable_column (&theader, &ftable, "V_DEC_ERR",   duD,             catalog->Naverage);
     527    gfits_set_bintable_column (&theader, &ftable, "PAR",         P,               catalog->Naverage);
     528    gfits_set_bintable_column (&theader, &ftable, "PAR_ERR",     dP,              catalog->Naverage);
     529    gfits_set_bintable_column (&theader, &ftable, "CHISQ_POS",   ChiSqAve,        catalog->Naverage);
     530    gfits_set_bintable_column (&theader, &ftable, "CHISQ_PM",    ChiSqPM,         catalog->Naverage);
     531    gfits_set_bintable_column (&theader, &ftable, "CHISQ_PAP",   ChiSqPar,        catalog->Naverage);
     532    gfits_set_bintable_column (&theader, &ftable, "MEAN_EPOCH",  Tmean,           catalog->Naverage);
     533    gfits_set_bintable_column (&theader, &ftable, "TIME_RANGE",  Trange,          catalog->Naverage);
     534    gfits_set_bintable_column (&theader, &ftable, "SIGMA_POS",   Xp,              catalog->Naverage);
     535    gfits_set_bintable_column (&theader, &ftable, "NUMBER_POS",  Npos,            catalog->Naverage);
     536    gfits_set_bintable_column (&theader, &ftable, "NMEASURE",    Nmeasure,        catalog->Naverage);
     537    gfits_set_bintable_column (&theader, &ftable, "NMISSING",    Nmissing,        catalog->Naverage);
     538    gfits_set_bintable_column (&theader, &ftable, "NEXTEND",     Nextend,         catalog->Naverage);
     539    gfits_set_bintable_column (&theader, &ftable, "OFF_MEASURE", measureOffset,   catalog->Naverage);
     540    gfits_set_bintable_column (&theader, &ftable, "OFF_MISSING", missingOffset,   catalog->Naverage);
     541    gfits_set_bintable_column (&theader, &ftable, "OFF_EXTEND",  extendOffset,    catalog->Naverage);
     542    gfits_set_bintable_column (&theader, &ftable, "FLAGS",       flags,           catalog->Naverage);
     543    gfits_set_bintable_column (&theader, &ftable, "PHOTFLAGS_U", photFlagsUpper,  catalog->Naverage);
     544    gfits_set_bintable_column (&theader, &ftable, "PHOTFLAGS_L", photFlagsLower,  catalog->Naverage);
     545    gfits_set_bintable_column (&theader, &ftable, "OBJ_ID",      objID,           catalog->Naverage);
     546    gfits_set_bintable_column (&theader, &ftable, "CAT_ID",      catID,           catalog->Naverage);
     547    gfits_set_bintable_column (&theader, &ftable, "EXT_ID",      extID,           catalog->Naverage);
     548
     549    free (R);
     550    free (D);
     551    free (dR);
     552    free (dD);
     553    free (uR);
     554    free (uD);
     555    free (duR);
     556    free (duD);
     557    free (P);
     558    free (dP);
     559    free (ChiSqAve);
     560    free (ChiSqPM);
     561    free (ChiSqPar);
     562    free (Tmean);
     563    free (Trange);
     564    free (Xp);
     565    free (Npos);
     566    free (Nmeasure);
     567    free (Nmissing);
     568    free (Nextend);
     569    free (measureOffset);
     570    free (missingOffset);
     571    free (extendOffset);
     572    free (flags);
     573    free (photFlagsUpper);
     574    free (photFlagsLower);
     575    free (objID);
     576    free (catID);
     577    free (extID);
    392578
    393579    gfits_fwrite_Theader (f, &theader);
    … …  
    486672
    487673  ALLOCATE (bcatalog[0].measure, MeasureTiny, Nmeas);
    488   ALLOCATE (bcatalog[0].average, AverageTiny, Naves);
    489   ALLOCATE (bcatalog[0].secfilt, SecFilt, Naves*Nsecfilt);
    490 
     674  ALLOCATE (bcatalog[0].average, Average,    Naves);
     675  ALLOCATE (bcatalog[0].secfilt, SecFilt,     Naves*Nsecfilt);
     676 
    491677  int Nm = 0;
    492678  int Na = 0;
    … …  
    494680    if (!catalog[i].Naverage) continue;
    495681    for (j = 0; j < catalog[i].Naverage; j++) {
    496       // CopyAverageTiny (&bcatalog[0].average[Na], &catalog[i].average[j]);
    497       bcatalog[0].average[Na] = catalog[i].averageT[j];
    498       bcatalog[0].average[Na].catID = catalog[i].catID;
     682      bcatalog[0].average[Na] = catalog[i].average[j];
    499683      for (k = 0; k < Nsecfilt; k++) {
    500684        bcatalog[0].secfilt[Nsecfilt*Na + k] = catalog[i].secfilt[Nsecfilt*j + k];
    … …  
    504688    }
    505689    for (j = 0; j < catalog[i].Nmeasure; j++) {
    506       // CopyMeasureTiny (&bcatalog[0].measure[Nm], &catalog[i].measure[j]);
    507690      bcatalog[0].measure[Nm] = catalog[i].measureT[j];
    508       bcatalog[0].measure[Nm].catID = catalog[i].catID;
    509691      Nm++;
    510692      assert (Nm <= Nmeas);
    … …  
    548730    catalogs->catalog[i].Naverage = 0;
    549731    catalogs->catalog[i].Nmeasure = 0;
    550     ALLOCATE (catalogs->catalog[i].averageT, AverageTiny, catalogs->NAVERAGE[i]);
     732    ALLOCATE (catalogs->catalog[i].average,  Average,    catalogs->NAVERAGE[i]);
    551733    ALLOCATE (catalogs->catalog[i].measureT, MeasureTiny, catalogs->NMEASURE[i]);
    552734    ALLOCATE (catalogs->catalog[i].secfilt,  SecFilt,     catalogs->NAVERAGE[i]*Nsecfilt);
    … …  
    557739// distribute a bright catalog across separate catalogs
    558740int BrightCatalogSplitFree (CatalogSplitter *catalogs) {
    559 
    560   // XXX don't free the catalogs :
    561   // for (i = 0; i < NCATALOG; i++) {
    562   //   free (catalogs->catalog[i].averageT);
    563   //   free (catalogs->catalog[i].measureT);
    564   //   free (catalogs->catalog[i].secfilt);
    565   // }
    566   // free (catalogs->catalog);
    567741
    568742  free (catalogs->catIDs);
    … …  
    609783    assert (Ncat < catalogs->NCATALOG);
    610784
     785    catalogs->catalog[Ncat].Nsecfilt = Nsecfilt;
     786    catalogs->catalog[Ncat].catID = catID;
     787
    611788    catalogs->index[catID] = Ncat;
    612789    assert (catID > 0);
    … …  
    639816        catalogs->catalog[j].Naverage = 0;
    640817        catalogs->catalog[j].Nmeasure = 0;
    641         ALLOCATE (catalogs->catalog[j].averageT, AverageTiny, catalogs->NAVERAGE[j]);
     818        ALLOCATE (catalogs->catalog[j].average,  Average,    catalogs->NAVERAGE[j]);
    642819        ALLOCATE (catalogs->catalog[j].measureT, MeasureTiny, catalogs->NMEASURE[j]);
    643820        ALLOCATE (catalogs->catalog[j].secfilt,  SecFilt,     catalogs->NAVERAGE[j]*Nsecfilt);
    … …  
    658835
    659836    int Na = catalogs->catalog[Nc].Naverage;
    660     catalogs->catalog[Nc].averageT[Na] = bcatalog->average[i];
     837    catalogs->catalog[Nc].average[Na] = bcatalog->average[i];
    661838
    662839    // secfilt entries are grouped in blocks for each average
    … …  
    669846    if (catalogs->catalog[Nc].Naverage >= catalogs->NAVERAGE[Nc]) {
    670847      catalogs->NAVERAGE[Nc] += 100;
    671       REALLOCATE (catalogs->catalog[Nc].averageT, AverageTiny, catalogs->NAVERAGE[Nc]);
    672       REALLOCATE (catalogs->catalog[Nc].secfilt,  SecFilt,     catalogs->NAVERAGE[Nc]*Nsecfilt);
     848      REALLOCATE (catalogs->catalog[Nc].average,  Average, catalogs->NAVERAGE[Nc]);
     849      REALLOCATE (catalogs->catalog[Nc].secfilt,  SecFilt, catalogs->NAVERAGE[Nc]*Nsecfilt);
    673850    }         
    674851  }

  • branches/eam_branches/ipp-20111122/Ohana/src/relastro/src/ConfigInit.c

    r31635 r33448  
    2626  if (!ScanConfig (config, "RELASTRO_IMFIT_SYS_SIGMA_LIM", "%lf", 0, &IMFIT_SYS_SIGMA_LIM)) IMFIT_SYS_SIGMA_LIM = 0.01;
    2727
    28   // XXX these are used in relphot to identify poor stars / images -- define
    29   // an equivalent concept here?
    30   // GetConfig (config, "STAR_SCATTER",           "%lf", 0, &STAR_SCATTER);
    31   // GetConfig (config, "IMAGE_SCATTER",          "%lf", 0, &IMAGE_SCATTER);
    32   // GetConfig (config, "IMAGE_OFFSET",           "%lf", 0, &IMAGE_OFFSET);
    33   // GetConfig (config, "STAR_CHISQ",             "%lf", 0, &STAR_CHISQ);
    34 
    3528  GetConfig (config, "PM_DT_MIN",              "%lf", 0, &PM_DT_MIN);
    3629  GetConfig (config, "PM_TOOFEW",              "%d",  0, &PM_TOOFEW);
    … …  
    4235  GetConfig (config, "ADDSTAR_RADIUS",         "%lf", 0, &ADDSTAR_RADIUS);
    4336
     37  // force CATDIR to be absolute (so parallel mode will work)
     38  GetConfig (config, "CATDIR",                 "%s",  0, CATDIR);
     39  char *tmpcatdir = abspath (CATDIR, MAX_PATH_LENGTH);
     40  strcpy (CATDIR, tmpcatdir);
     41  free (tmpcatdir);
     42
    4443  GetConfig (config, "GSCFILE",                "%s",  0, GSCFILE);
    45   GetConfig (config, "CATDIR",                 "%s",  0, CATDIR);
    4644  ScanConfig(config, "CATMODE",                "%s",  0, CATMODE);
    4745  ScanConfig(config, "CATFORMAT",              "%s",  0, CATFORMAT);

  • branches/eam_branches/ipp-20111122/Ohana/src/relastro/src/FixProblemImages.c

    r33422 r33448  
    6767    }
    6868
    69     catalog = load_catalogs (&sublist, &Ncat, FALSE);
     69    // XXX use a different function here
     70    // catalog = load_catalogs (&sublist, &Ncat, FALSE);
    7071    assert (Ncat == Ncatlist);
    7172

  • branches/eam_branches/ipp-20111122/Ohana/src/relastro/src/GetAstromError.c

    r30616 r33448  
    11# include "relastro.h"
    22# define WEIGHTED_ERRORS 1
     3
     4float GetAstromErrorTiny (MeasureTiny *measure, int mode) {
     5
     6  PhotCode *code;
     7  float dPobs, dPsys, dPtotal, dM, AS, MS, dX, dY;
     8
     9  if (!WEIGHTED_ERRORS) {
     10    // if we don't understand the errors at all, this at least lets us get things roughly
     11    // right:
     12    return 0.1;
     13  }
     14
     15  switch (mode) {
     16    case ERROR_MODE_RA:
     17      dPobs = FromShortPixels(measure[0].dXccd);  // dXccd is a value in pixels
     18      break;
     19    case ERROR_MODE_DEC:
     20      dPobs = FromShortPixels(measure[0].dYccd);  // dYccd is a value in pixels
     21      break;
     22    case ERROR_MODE_POS:
     23      dX = FromShortPixels(measure[0].dXccd);  // dXccd is a value in pixels
     24      dY = FromShortPixels(measure[0].dYccd);  // dYccd is a value in pixels
     25      dPobs = hypot (dX, dY);
     26      break;
     27    default:
     28      abort();
     29  }
     30
     31  code  = GetPhotcodebyCode (measure[0].photcode);
     32  AS    = code[0].astromErrScale;
     33  MS    = code[0].astromErrMagScale;
     34  dPsys = code[0].astromErrSys;
     35  dM    = measure[0].dM;
     36  dPtotal = sqrt(SQ(dPsys) + SQ(AS*dPobs) + SQ(MS*dM));
     37
     38  dPtotal = MAX (dPtotal, MIN_ERROR);
     39  return (dPtotal);
     40}
    341
    442float GetAstromError (Measure *measure, int mode) {

  • branches/eam_branches/ipp-20111122/Ohana/src/relastro/src/ImageOps.c

    r33422 r33448  
    333333}
    334334
    335 void dump_measures(Average *average, MeasureTiny *measure) {
     335void dump_measures (Average *average, Measure *measure) {
    336336
    337337  off_t j, off;
    … …  
    406406    c = clist[im][i];
    407407
    408     MeasureTiny *measure = &catalog[c].measureT[m];
     408    Measure *measure = &catalog[c].measure[m];
    409409
    410410    X = measure[0].Xccd;
    … …  
    431431    if (fabs(dR) > 3.0*ADDSTAR_RADIUS) {
    432432      fprintf (stderr, "measurement is far from average location (R): %f %f (%f %f)\n", catalog[c].average[n].R, catalog[c].average[n].D, dR, dD);
    433       dump_measures (&catalog[c].average[n], catalog[c].measureT);
     433      dump_measures (&catalog[c].average[n], catalog[c].measure);
    434434      // abort ();
    435435    }
    436436    if (fabs(dD) > 3.0*ADDSTAR_RADIUS) {
    437437      fprintf (stderr, "measurement is far from average location (D): %f %f (%f %f)\n", catalog[c].average[n].R, catalog[c].average[n].D, dR, dD);
    438       dump_measures (&catalog[c].average[n], catalog[c].measureT);
     438      dump_measures (&catalog[c].average[n], catalog[c].measure);
    439439      // abort ();
    440440    }
    … …  
    443443    if (fabs(catalog[c].measureT[m].dR - dR) > DPOS_MAX_ASEC) {
    444444      fprintf (stderr, "measurement is far from original location (R): %f %f (%f %f)\n", catalog[c].average[n].R, catalog[c].average[n].D, dR, dD);
    445       dump_measures (&catalog[c].average[n], catalog[c].measureT);
     445      dump_measures (&catalog[c].average[n], catalog[c].measure);
    446446      // abort();
    447447    }
    … …  
    581581    c = clist[im][i];
    582582
    583     MeasureTiny *measure = &catalog[c].measure[m];
     583    MeasureTiny *measure = &catalog[c].measureT[m];
    584584
    585585    /* apply the current image transformation or use the current value of R+dR, D+dD? */
    … …  
    589589    raw[i].Mag  = measure[0].M;
    590590    raw[i].dMag = measure[0].dM;
    591     raw[i].dPos = GetAstromError (&measure[0], ERROR_MODE_POS);
     591    raw[i].dPos = GetAstromErrorTiny (&measure[0], ERROR_MODE_POS);
    592592
    593593    n = measure[0].averef;
    … …  
    654654    m = mlist[im][i];
    655655    c = clist[im][i];
    656     n = catalog[c].measure[m].averef;
    657 
    658     MeasureTiny *measure = &catalog[c].measure[m];
     656
     657    MeasureTiny *measure = &catalog[c].measureT[m];
     658    n = measure[0].averef;
    659659
    660660    /* apply the current image transformation or use the current value of R+dR, D+dD? */
    … …  
    664664    ref[i].Mag  = measure[0].M;
    665665    ref[i].dMag = measure[0].dM;
    666     ref[i].dPos = GetAstromError (&measure[0], ERROR_MODE_POS);
     666    ref[i].dPos = GetAstromErrorTiny (&measure[0], ERROR_MODE_POS);
    667667
    668668    ref[i].mask = FALSE;
    … …  
    694694void FlagOutliers2D(Catalog *catalog);
    695695
     696// operates on Full values (not tiny)
    696697void FlagOutliers (Catalog *catalog) {
    697698
    … …  
    731732    // pointer to this set of measurements
    732733    m = catalog[0].average[j].measureOffset;
    733     MeasureTiny *measure = &catalog[0].measureT[m];
     734    Measure *measure = &catalog[0].measure[m];
    734735
    735736    /* accumulate list of valid measurements */
    … …  
    825826    // pointer to this set of measurements
    826827    m = catalog[0].average[j].measureOffset;
    827     MeasureTiny *measure = &catalog[0].measureT[m];
     828    Measure *measure = &catalog[0].measure[m];
    828829
    829830    /* accumulate list of valid measurements */
    … …  
    915916}
    916917
    917 
    918 
    919918/** Determine whether a measurement should be included in the analysis, based on supplied filter criteria */
    920919// we only optionally apply the sigma limit: for object averages, this should not be used (should it?)
    921 int MeasFilterTest(MeasureTiny *measure, int applySigmaLim) {
     920int MeasFilterTestTiny(MeasureTiny *measure, int applySigmaLim) {
    922921  int found, k;
    923922  long mask;
    … …  
    978977  return TRUE;
    979978}
     979
     980/** Determine whether a measurement should be included in the analysis, based on supplied filter criteria */
     981// we only optionally apply the sigma limit: for object averages, this should not be used (should it?)
     982int MeasFilterTest(Measure *measure, int applySigmaLim) {
     983  int found, k;
     984  long mask;
     985  PhotCode *code;
     986  float mag;
     987
     988  if (!finite(measure[0].dR) || !finite(measure[0].dD)) return FALSE;
     989  if (!finite(measure[0].M)) return FALSE; //XXX is this necessary for all relastro tasks?
     990  if (!finite(measure[0].dM)) return FALSE; //XXX is this necessary for all relastro tasks?
     991 
     992  /* select measurements by photcode, or equiv photcode, if specified */
     993  if (NphotcodesKeep > 0) {
     994    found = FALSE;
     995    for (k = 0; (k < NphotcodesKeep) && !found; k++) {
     996      if (photcodesKeep[k][0].code == measure[0].photcode) found = TRUE;
     997      if (photcodesKeep[k][0].code == GetPhotcodeEquivCodebyCode(measure[0].photcode)) found = TRUE;
     998    }
     999    if (!found) return FALSE;
     1000  }
     1001 
     1002  if (NphotcodesSkip > 0) {
     1003    found = FALSE;
     1004    for (k = 0; (k < NphotcodesSkip) && !found; k++) {
     1005      if (photcodesSkip[k][0].code == measure[0].photcode) found = TRUE;
     1006      if (photcodesSkip[k][0].code == GetPhotcodeEquivCodebyCode(measure[0].photcode)) found = TRUE;
     1007    }
     1008    if (found) return FALSE;
     1009  } 
     1010 
     1011  /* select measurements by time */
     1012  if (TimeSelect) {
     1013    if (measure[0].t < TSTART) return FALSE;
     1014    if (measure[0].t > TSTOP) return FALSE;
     1015  }
     1016 
     1017  /* select measurements by quality */
     1018  if (PhotFlagSelect) {
     1019    if (PhotFlagBad) {
     1020      mask = PhotFlagBad;
     1021    } else {
     1022      code = GetPhotcodebyCode (measure[0].photcode);
     1023      mask = code[0].astromBadMask;
     1024    }
     1025    if (mask & measure[0].photFlags) return FALSE;
     1026  }
     1027
     1028  /* select measurements by measurement error */
     1029  if (applySigmaLim && (SIGMA_LIM > 0) && (measure[0].dM > SIGMA_LIM)) {
     1030    return FALSE;
     1031  }
     1032 
     1033  /* select measurements by mag limit */
     1034  if (ImagSelect) {
     1035    mag = PhotInst (measure);
     1036    if (mag < ImagMin || mag > ImagMax) return FALSE;
     1037  }
     1038 
     1039  return TRUE;
     1040}

  • branches/eam_branches/ipp-20111122/Ohana/src/relastro/src/StarMaps.c

    r33422 r33448  
    6868  for (i = 0; i < catalog[0].Nmeasure; i++) {
    6969   
    70     MeasureTiny *measure = &catalog[0].measureT[i]
     70    MeasureTiny *measure = &catalog[0].measureT[i];
    7171
    7272    N = getImageByID(measure[0].imageID);
    … …  
    8484
    8585  return (TRUE);
     86}
     87
     88int createStarMap (Catalog *catalog, int Ncatalog) {
     89
     90  int i;
     91
     92  initStarMaps();
     93
     94  for (i = 0; i < Ncatalog; i++) {
     95    // check coverage per chip (operates on (Average, MeasureTiny, Secfilt)
     96    updateStarMaps (&catalog[i]);
     97  }
     98
     99  createStarMapPoints();
     100  return TRUE;
    86101}
    87102

  • branches/eam_branches/ipp-20111122/Ohana/src/relastro/src/UpdateChips.c

    r32346 r33448  
    143143  }
    144144
    145   double Rmid;
    146   if (UserCatalog) {
    147       Rmid = UserCatalogRA;
    148   } else {
    149       Rmid = 0.5*(UserPatch.Rmin + UserPatch.Rmax);
    150   }
     145  double Rmid = 0.5*(UserPatch.Rmin + UserPatch.Rmax);
    151146
    152147  R = ohana_normalize_angle_to_midpoint (R, Rmid);

  • branches/eam_branches/ipp-20111122/Ohana/src/relastro/src/UpdateMeasures.c

    r32346 r33448  
    11# include "relastro.h"
    22
     3// this function operates on Measure, not MeasureTiny
    34int UpdateMeasures (Catalog *catalog, int Ncatalog) {
    45

  • branches/eam_branches/ipp-20111122/Ohana/src/relastro/src/UpdateObjectOffsets.c

    r33422 r33448  
    77//   original coordinates
    88
    9 int UpdateObjectOffsets (SkyList *skylist) {
     9int UpdateObjectOffsets_parallel (SkyList *sky);
     10
     11int UpdateObjectOffsets (SkyList *skylist, int hostID, char *hostpath) {
    1012
    1113  int i;
    1214  Catalog catalog;
    1315
     16  // XXX need to decide how to determine PARALLEL mode...
     17  if (PARALLEL & !hostID) {
     18    UpdateObjectOffsets_parallel (skylist);
     19    return TRUE;
     20  }
     21
    1422  // load data from each region file, only use bright stars
    1523  for (i = 0; i < skylist[0].Nregions; i++) {
    1624
     25    if (hostID && (skylist[0].regions[i]->hostID != hostID)) continue;
     26
     27    // define the catalog file name
     28    char hostfile[1024];
     29    snprintf (hostfile, 1024, "%s/%s.cpt", hostpath, skylist[0].regions[i]->name);
     30    catalog.filename = hostID ? hostfile : skylist[0].filename[i];
     31
    1732    // set up the basic catalog info
    18     catalog.filename  = skylist[0].filename[i];
    1933    catalog.catformat = dvo_catalog_catformat (CATFORMAT);    // set the default catformat from config data
    2034    catalog.catmode   = dvo_catalog_catmode (CATMODE);        // set the default catmode from config data
    … …  
    3347    }
    3448
     49    // set the values in MeasureTiny needed by UpdateObjects
     50    populate_tiny_values(&catalog, DVO_TV_MEASURE);
     51
    3552    // match measurements with images
    3653    initImageBins (&catalog, 1, FALSE);
    … …  
    4057    UpdateMeasures (&catalog, 1);
    4158
    42   XXX: note that this version of UpdateObjects is / should be expecting catalog.measure, not catalog.measureT
     59    // uses MeasureTiny values, but will also update Measure values when populated
    4360    UpdateObjects (&catalog, 1);
    4461
    … …  
    5168  return (TRUE);
    5269}
     70
     71// CATDIR is supplied globally
     72# define DEBUG 1
     73int UpdateObjectOffsets_parallel (SkyList *sky) {
     74
     75  // launch the setphot_client jobs to the parallel hosts
     76
     77  // load the list of hosts
     78  HostTable *table = HostTableLoad (CATDIR, sky->hosts);
     79  if (!table) {
     80    fprintf (stderr, "ERROR: failure reading Host Table %s for database %s\n", sky->hosts, CATDIR);
     81    exit (1);
     82  }   
     83
     84  int i;
     85  for (i = 0; i < table->Nhosts; i++) {
     86
     87    // ensure that the paths are absolute path names
     88    char *tmppath = abspath (table->hosts[i].pathname, MAX_PATH_LENGTH);
     89    free (table->hosts[i].pathname);
     90    table->hosts[i].pathname = tmppath;
     91
     92    // options / arguments that can affect relastro_client -load:
     93    // VERBOSE, VERBOSE2
     94    // RESET (-reset)
     95    // TimeSelect -time
     96    // (note that psfQual is applied rigidly at 0.85, as is the galaxy test)
     97    // ImagSelect, ImagMin, ImagMax
     98    // MaxDensityUse, MaxDensityValue
     99
     100    char command[1024];
     101    snprintf (command, 1024, "relastro_client -update-offsets -hostID %d -D CATDIR %s -hostdir %s -region %f %f %f %f -statmode %s -minerror %f",
     102              table->hosts[i].hostID, CATDIR, table->hosts[i].pathname, UserPatch.Rmin, UserPatch.Rmax, UserPatch.Dmin, UserPatch.Dmax, STATMODE, MIN_ERROR);
     103
     104    char tmpline[1024];
     105    if (FIT_MODE == FIT_PM_ONLY)         { snprintf (tmpline, 1024, "%s -pm",    command);           strcpy (command, tmpline); }
     106    if (FIT_MODE == FIT_PAR_ONLY)        { snprintf (tmpline, 1024, "%s -par",   command);           strcpy (command, tmpline); }
     107    if (FIT_MODE == FIT_PM_AND_PAR)      { snprintf (tmpline, 1024, "%s -pmpar", command);           strcpy (command, tmpline); }
     108
     109    if (VERBOSE)       { snprintf (tmpline, 1024, "%s -v",              command);                    strcpy (command, tmpline); }
     110    if (VERBOSE2)      { snprintf (tmpline, 1024, "%s -vv",             command);                    strcpy (command, tmpline); }
     111    if (RESET)         { snprintf (tmpline, 1024, "%s -reset",          command);                    strcpy (command, tmpline); }
     112    if (ImagSelect)    { snprintf (tmpline, 1024, "%s -instmag %f %f",  command, ImagMin, ImagMax);  strcpy (command, tmpline); }
     113    if (MaxDensityUse) { snprintf (tmpline, 1024, "%s -max-density %f", command, MaxDensityValue);   strcpy (command, tmpline); }
     114   
     115    if (FlagOutlier)   { snprintf (tmpline, 1024, "%s -clip %d",        command, CLIP_THRESH);       strcpy (command, tmpline); }
     116   
     117    if (PHOTCODE_KEEP_LIST) { snprintf (tmpline, 1024, "%s +photcode %s", command, PHOTCODE_KEEP_LIST); strcpy (command, tmpline); }
     118    if (PHOTCODE_SKIP_LIST) { snprintf (tmpline, 1024, "%s -photcode %s", command, PHOTCODE_SKIP_LIST); strcpy (command, tmpline); }
     119    if (PhotFlagSelect)    { snprintf (tmpline, 1024, "%s +photflags",   command);                     strcpy (command, tmpline); }
     120    if (PhotFlagBad)       { snprintf (tmpline, 1024, "%s +photflagbad %d", command, PhotFlagBad);     strcpy (command, tmpline); }
     121    if (PhotFlagPoor)      { snprintf (tmpline, 1024, "%s +photflagpoor %d", command, PhotFlagPoor);   strcpy (command, tmpline); }
     122    // XXX note that the above pass in the flag as decimal -- also note that args.c cannot handle 0xHEX values
     123
     124    if (TimeSelect) {
     125      char *tstart = ohana_sec_to_date (TSTART);
     126      char *tstop  = ohana_sec_to_date (TSTOP);
     127      snprintf (tmpline, 1024, "%s -time %s %s", command, tstart, tstop);
     128      free (tstart);
     129      free (tstop);
     130      strcpy (command, tmpline);
     131    }
     132    fprintf (stderr, "command: %s\n", command);
     133
     134    if (PARALLEL_MANUAL) continue;
     135
     136    if (PARALLEL_SERIAL) {
     137      int status = system (command);
     138      if (status) {
     139        fprintf (stderr, "ERROR running relastro_client\n");
     140        exit (2);
     141      }
     142    } else {
     143      // launch the job on the remote machine (no handshake)
     144      int errorInfo = 0;
     145      int pid = rconnect ("ssh", table->hosts[i].hostname, command, table->hosts[i].stdio, &errorInfo, FALSE);
     146      if (!pid) {
     147        if (DEBUG) fprintf (stderr, "failure to start %s (error %d)\n", table->hosts[i].hostname, errorInfo);
     148        exit (1);
     149      }
     150      table->hosts[i].pid = pid; // save for future reference
     151    }
     152  }
     153
     154  if (PARALLEL_MANUAL) {
     155    fprintf (stderr, "run the relastro_client commands above.  when these are done, hit return\n");
     156    getchar();
     157  }
     158  if (!PARALLEL_MANUAL && !PARALLEL_SERIAL) {
     159    HostTableWaitJobs (table, __FILE__, __LINE__);
     160  }
     161
     162  return TRUE;
     163}     

  • branches/eam_branches/ipp-20111122/Ohana/src/relastro/src/UpdateObjects.c

    r33422 r33448  
    3838} 
    3939
     40// This function operates on both Measure and MeasureTiny.  In the big stages, this should
     41// be called with just MeasureTiny set and Measure == NULL
    4042int UpdateObjects (Catalog *catalog, int Ncatalog) {
    4143
    … …  
    9395      m = catalog[i].average[j].measureOffset;
    9496      MeasureTiny *measure = &catalog[i].measureT[m];
     97      Measure *measureBig = catalog[i].measure ? &catalog[i].measure[m] : NULL;
     98      // when we update the output measure values, we need to do it here
    9599
    96100      Tmin = Tmax = (measure[0].t - T2000) / (86400*365.25);
    … …  
    102106
    103107        // does the measurement pass the supplied filtering constraints?
    104         if (!MeasFilterTest(&measure[k], FALSE)) {
     108        if (!MeasFilterTestTiny(&measure[k], FALSE)) {
    105109          measure[k].dbFlags &= ~ID_MEAS_USED_OBJ;
     110          if (measureBig) { measureBig[k].dbFlags &= ~ID_MEAS_USED_OBJ; }
    106111          continue;
    107112        }
    … …  
    110115        if (FlagOutlier && (measure[k].dbFlags & ID_MEAS_POOR_ASTROM)) {
    111116          measure[k].dbFlags &= ~ID_MEAS_USED_OBJ;
     117          if (measureBig) { measureBig[k].dbFlags &= ~ID_MEAS_USED_OBJ; }
    112118          continue;
    113119        }
    114120
    115121        measure[k].dbFlags |= ID_MEAS_USED_OBJ;
     122        if (measureBig) { measureBig[k].dbFlags |= ID_MEAS_USED_OBJ; }
    116123
    117124        R[N] = getMeanR (&measure[k], &catalog[i].average[j], &catalog[i].secfilt[j*Nsecfilt]);
    … …  
    127134
    128135        // dX, dY : error in arcsec --
    129         dX[N] = GetAstromError (&measure[k], ERROR_MODE_RA);
    130         dY[N] = GetAstromError (&measure[k], ERROR_MODE_DEC);
     136        dX[N] = GetAstromErrorTiny (&measure[k], ERROR_MODE_RA);
     137        dY[N] = GetAstromErrorTiny (&measure[k], ERROR_MODE_DEC);
    131138
    132139        // add systematic error in quadrature, if desired
    … …  
    293300
    294301      // the measure fields must be updated before the average fields
    295       m = catalog[i].average[j].measureOffset;
    296       MeasureTiny *measure = &catalog[i].measureT[m];
    297       for (k = 0; k < catalog[i].average[j].Nmeasure; k++, m++) {
     302      for (k = 0; k < catalog[i].average[j].Nmeasure; k++) {
    298303        setMeanR (fit.Ro, &measure[k], &catalog[i].average[j], &catalog[i].secfilt[j*Nsecfilt]);
    299304        setMeanD (fit.Do, &measure[k], &catalog[i].average[j], &catalog[i].secfilt[j*Nsecfilt]);
     305        if (measureBig) {
     306          setMeanR_Big (fit.Ro, &measureBig[k], &catalog[i].average[j], &catalog[i].secfilt[j*Nsecfilt]);
     307          setMeanD_Big (fit.Do, &measureBig[k], &catalog[i].average[j], &catalog[i].secfilt[j*Nsecfilt]);
     308        }
    300309      }     
    301310

  • branches/eam_branches/ipp-20111122/Ohana/src/relastro/src/args.c

    r32695 r33448  
    11# include "relastro.h"
    22void usage (void);
     3void usage_client (void);
    34void usage_merge_source (void);
    45void usage_merge_source_id (char *name);
    … …  
    4243  if ((N = get_argument (argc, argv, "-update-objects"))) {
    4344    remove_argument (N, &argc, argv);
    44     FIT_TARGET = TARGET_OBJECTS;
     45    FIT_TARGET = TARGET_UPDATE_OBJECTS;
    4546  }
    4647
    … …  
    9293    RADIUS = atof(argv[N]);
    9394    remove_argument (N, &argc, argv);
    94     HIGH_SPEED_DIR = strcreate(argv[N]);
     95    HIGH_SPEED_DIR = abspath(argv[N], MAX_PATH_LENGTH);
    9596    remove_argument (N, &argc, argv);
    9697  }
    … …  
    118119
    119120  if (FIT_TARGET == TARGET_NONE) usage();
    120 
    121   UserCatalog = FALSE;
    122   if ((N = get_argument (argc, argv, "-catalog"))) {
    123     remove_argument (N, &argc, argv);
    124     UserCatalogRA = atof(argv[N]);
    125     remove_argument (N, &argc, argv);
    126     UserCatalogDEC = atof(argv[N]);
    127     remove_argument (N, &argc, argv);
    128     UserCatalog = TRUE;
    129   }
    130121
    131122  /* specify portion of the sky : allow default of all sky? */
    … …  
    144135    UserPatch.Dmax = atof (argv[N]);
    145136    remove_argument (N, &argc, argv);
    146   } else {
    147     if (!UserCatalog) {
    148       usage ();
     137  }
     138  if ((N = get_argument (argc, argv, "-catalog"))) {
     139    remove_argument (N, &argc, argv);
     140    UserPatch.Rmin = atof (argv[N]);
     141    UserPatch.Rmax = UserPatch.Rmin + 0.001;
     142    remove_argument (N, &argc, argv);
     143    UserPatch.Dmin = atof (argv[N]);
     144    UserPatch.Dmax = UserPatch.Dmin + 0.001;
     145    remove_argument (N, &argc, argv);
     146  }
     147
     148  // XXX for the moment, make this selection manual.  it needs to be automatic
     149  // based on the state of the SkyTable
     150  PARALLEL = FALSE;
     151  if ((N = get_argument (argc, argv, "-parallel"))) {
     152    PARALLEL = TRUE;
     153    remove_argument (N, &argc, argv);
     154  }
     155  // this is a test mode : rather than launching the remote jobs and waiting for completion,
     156  // relastro will simply list the remote command and wait for the user to signal completion
     157  PARALLEL_MANUAL = FALSE;
     158  if ((N = get_argument (argc, argv, "-parallel-manual"))) {
     159    PARALLEL = TRUE; // -parallel-manual implies -parallel
     160    PARALLEL_MANUAL = TRUE;
     161    remove_argument (N, &argc, argv);
     162  }
     163  // this is a test mode : rather than launching the relastro_client jobs remotely, they are
     164  // run in serial via 'system'
     165  PARALLEL_SERIAL = FALSE;
     166  if ((N = get_argument (argc, argv, "-parallel-serial"))) {
     167    if (PARALLEL_MANUAL) {
     168      fprintf (stderr, "ERROR: cannot mix -parallel-manual and -parallel-serial\n");
     169      exit (1);
    149170    }
     171    PARALLEL = TRUE; // -parallel-serial implies -parallel
     172    PARALLEL_SERIAL = TRUE;
     173    remove_argument (N, &argc, argv);
    150174  }
    151175
    … …  
    290314  } 
    291315
    292   AreaSelect = FALSE;
    293   if ((N = get_argument (argc, argv, "-area"))) {
    294     remove_argument (N, &argc, argv);
    295     AreaXmin = atof (argv[N]);
    296     remove_argument (N, &argc, argv);
    297     AreaXmax = atof (argv[N]);
    298     remove_argument (N, &argc, argv);
    299     AreaYmin = atof (argv[N]);
    300     remove_argument (N, &argc, argv);
    301     AreaYmax = atof (argv[N]);
    302     remove_argument (N, &argc, argv);
    303     AreaSelect = TRUE;
    304   }
    305 
    306316  ImagSelect = FALSE;
    307317  if ((N = get_argument (argc, argv, "-instmag"))) {
    … …  
    336346  }
    337347
    338   /* XXX drop this? */
    339   DophotSelect = FALSE;
    340   if ((N = get_argument (argc, argv, "-dophot"))) {
    341     remove_argument (N, &argc, argv);
    342     DophotValue = atof (argv[N]);
    343     remove_argument (N, &argc, argv);
    344     DophotSelect = TRUE;
    345   }
    346 
    347348  NLOOP = 4;
    348349  if ((N = get_argument (argc, argv, "-nloop"))) {
    … …  
    356357}
    357358
     359int args_client (int argc, char **argv) {
     360
     361  int N;
     362  double trange;
     363
     364  /* possible operations */
     365  FIT_TARGET = TARGET_NONE;
     366  FIT_MODE = FIT_AVERAGE;
     367
     368  // by definition, the client is not parallel
     369  PARALLEL = FALSE;
     370  PARALLEL_MANUAL = FALSE;
     371  PARALLEL_SERIAL = FALSE;
     372
     373  BCATALOG = NULL;
     374
     375  HOST_ID = 0;
     376  if ((N = get_argument (argc, argv, "-hostID"))) {
     377    remove_argument (N, &argc, argv);
     378    HOST_ID = atoi (argv[N]);
     379    remove_argument (N, &argc, argv);
     380  }
     381  if (!HOST_ID) usage_client();
     382
     383  HOSTDIR = NULL;
     384  if ((N = get_argument (argc, argv, "-hostdir"))) {
     385    remove_argument (N, &argc, argv);
     386    HOSTDIR = strcreate (argv[N]);
     387    remove_argument (N, &argc, argv);
     388  }
     389  if (!HOSTDIR) usage_client();
     390
     391  if ((N = get_argument (argc, argv, "-load-objects"))) {
     392    remove_argument (N, &argc, argv);
     393    BCATALOG = strcreate(argv[N]);
     394    remove_argument (N, &argc, argv);
     395    FIT_TARGET = TARGET_LOAD_OBJECTS;
     396  }
     397
     398  if ((N = get_argument (argc, argv, "-update-objects"))) {
     399    remove_argument (N, &argc, argv);
     400    FIT_TARGET = TARGET_UPDATE_OBJECTS;
     401  }
     402
     403  if ((N = get_argument (argc, argv, "-update-offsets"))) {
     404    remove_argument (N, &argc, argv);
     405    FIT_TARGET = TARGET_UPDATE_OFFSETS;
     406  }
     407
     408  // check for object fitting modes
     409  if ((N = get_argument (argc, argv, "-pm"))) {
     410    remove_argument (N, &argc, argv);
     411    FIT_MODE = FIT_PM_ONLY;
     412  }
     413  if ((N = get_argument (argc, argv, "-par"))) {
     414    remove_argument (N, &argc, argv);
     415    FIT_MODE = FIT_PAR_ONLY;
     416  }
     417  if ((N = get_argument (argc, argv, "-pmpar"))) {
     418    remove_argument (N, &argc, argv);
     419    FIT_MODE = FIT_PM_AND_PAR;
     420  }
     421
     422  if ((N = get_argument (argc, argv, "-high-speed"))) {
     423    // XXX include a parallax / no-parallax option
     424    if (N >= argc - 3) usage_client();
     425    FIT_TARGET = TARGET_HIGH_SPEED;
     426    remove_argument (N, &argc, argv);
     427    PHOTCODE_A_LIST = strcreate(argv[N]);
     428    remove_argument (N, &argc, argv);
     429    PHOTCODE_B_LIST = strcreate(argv[N]);
     430    remove_argument (N, &argc, argv);
     431    RADIUS = atof(argv[N]);
     432    remove_argument (N, &argc, argv);
     433    HIGH_SPEED_DIR = strcreate(argv[N]);
     434    remove_argument (N, &argc, argv);
     435  }
     436
     437  FlagOutlier = FALSE;
     438  if ((N = get_argument (argc, argv, "-clip"))) {
     439    remove_argument (N, &argc, argv);
     440    CLIP_THRESH = atof (argv[N]);
     441    remove_argument (N, &argc, argv);
     442    FlagOutlier = TRUE;
     443  }
     444
     445  if (FIT_TARGET == TARGET_NONE) usage_client();
     446
     447  /* specify portion of the sky : allow default of all sky? */
     448  UserPatch.Rmin = 0;
     449  UserPatch.Rmax = 360;
     450  UserPatch.Dmin = -90;
     451  UserPatch.Dmax = +90;
     452  if ((N = get_argument (argc, argv, "-region"))) {
     453    remove_argument (N, &argc, argv);
     454    UserPatch.Rmin = atof (argv[N]);
     455    remove_argument (N, &argc, argv);
     456    UserPatch.Rmax = atof (argv[N]);
     457    remove_argument (N, &argc, argv);
     458    UserPatch.Dmin = atof (argv[N]);
     459    remove_argument (N, &argc, argv);
     460    UserPatch.Dmax = atof (argv[N]);
     461    remove_argument (N, &argc, argv);
     462  }
     463  if ((N = get_argument (argc, argv, "-catalog"))) {
     464    remove_argument (N, &argc, argv);
     465    UserPatch.Rmin = atof (argv[N]);
     466    UserPatch.Rmax = UserPatch.Rmin + 0.001;
     467    remove_argument (N, &argc, argv);
     468    UserPatch.Dmin = atof (argv[N]);
     469    UserPatch.Dmax = UserPatch.Dmin + 0.001;
     470    remove_argument (N, &argc, argv);
     471  }
     472
     473  USE_BASIC_CHECK = FALSE;
     474  if ((N = get_argument (argc, argv, "-basic-image-search"))) {
     475    remove_argument (N, &argc, argv);
     476    USE_BASIC_CHECK = TRUE;
     477  }
     478
     479  MaxDensityUse = FALSE;
     480  if ((N = get_argument (argc, argv, "-max-density"))) {
     481    remove_argument (N, &argc, argv);
     482    MaxDensityValue = atof(argv[N]);
     483    remove_argument (N, &argc, argv);
     484    MaxDensityUse = TRUE;
     485  }
     486
     487  /* define time */
     488  TimeSelect = FALSE;
     489  if ((N = get_argument (argc, argv, "-time"))) {
     490    remove_argument (N, &argc, argv);
     491    if (!ohana_str_to_time (argv[N], &TSTART)) {
     492      fprintf (stderr, "ERROR: syntax error\n");
     493      return (FALSE);
     494    }
     495    remove_argument (N, &argc, argv);
     496    if (!ohana_str_to_dtime (argv[N], &trange)) {
     497      if (!ohana_str_to_time (argv[N], &TSTOP)) {
     498        fprintf (stderr, "ERROR: syntax error\n");
     499        return (FALSE);
     500      }
     501    } else {
     502      if (trange < 0) {
     503        trange = fabs (trange);
     504        TSTOP = TSTART;
     505        TSTART -= trange;
     506      } else {
     507        TSTOP = TSTART + trange;
     508      }
     509    }
     510    remove_argument (N, &argc, argv);
     511    TimeSelect = TRUE;
     512  }
     513
     514  PHOTCODE_KEEP_LIST = NULL;
     515  if ((N = get_argument (argc, argv, "+photcode"))) {
     516    remove_argument (N, &argc, argv);
     517    PHOTCODE_KEEP_LIST = strcreate(argv[N]);
     518    remove_argument (N, &argc, argv);
     519  }
     520
     521  PHOTCODE_SKIP_LIST = NULL;
     522  if ((N = get_argument (argc, argv, "-photcode"))) {
     523    remove_argument (N, &argc, argv);
     524    PHOTCODE_SKIP_LIST = strcreate(argv[N]);
     525    remove_argument (N, &argc, argv);
     526  }
     527
     528  VERBOSE = VERBOSE2 = FALSE;
     529  if ((N = get_argument (argc, argv, "-v"))) {
     530    VERBOSE = TRUE;
     531    remove_argument (N, &argc, argv);
     532  }
     533  if ((N = get_argument (argc, argv, "-vv"))) {
     534    VERBOSE = VERBOSE2 = TRUE;
     535    remove_argument (N, &argc, argv);
     536  }
     537
     538  // by default, require > 10pts to clip
     539  strcpy (STATMODE, "CHI_INNER_80_WTMEAN");
     540  if ((N = get_argument (argc, argv, "-statmode"))) {
     541    remove_argument (N, &argc, argv);
     542    strcpy (STATMODE, argv[N]);
     543    remove_argument (N, &argc, argv);
     544  }
     545
     546  RESET = FALSE;
     547  if ((N = get_argument (argc, argv, "-reset"))) {
     548    remove_argument (N, &argc, argv);
     549    RESET = TRUE;
     550  }
     551
     552  // do not hold the lock on the Image table
     553  UPDATE = FALSE;
     554
     555  /* XXX update this */
     556  MIN_ERROR = 0.001;
     557  if ((N = get_argument (argc, argv, "-minerror"))) {
     558    remove_argument (N, &argc, argv);
     559    MIN_ERROR = atof (argv[N]);
     560    remove_argument (N, &argc, argv);
     561    /* require MIN_ERROR > 0 */
     562  } 
     563
     564  ImagSelect = FALSE;
     565  if ((N = get_argument (argc, argv, "-instmag"))) {
     566    remove_argument (N, &argc, argv);
     567    ImagMin = atof (argv[N]);
     568    remove_argument (N, &argc, argv);
     569    ImagMax = atof (argv[N]);
     570    remove_argument (N, &argc, argv);
     571    ImagSelect = TRUE;
     572  }
     573 
     574  // for now, make the default to ignore the photflags
     575  // XXX make it true by default instead?
     576  PhotFlagSelect = FALSE;
     577  if ((N = get_argument (argc, argv, "+photflags"))) {
     578    remove_argument (N, &argc, argv);
     579    PhotFlagSelect = TRUE;
     580  }
     581
     582  PhotFlagBad = 0;
     583  if ((N = get_argument (argc, argv, "-photflagbad"))) {
     584    remove_argument (N, &argc, argv);
     585    PhotFlagBad = atoi (argv[N]);
     586    remove_argument (N, &argc, argv);
     587  }
     588
     589  PhotFlagPoor = 0;
     590  if ((N = get_argument (argc, argv, "-photflagpoor"))) {
     591    remove_argument (N, &argc, argv);
     592    PhotFlagPoor = atoi (argv[N]);
     593    remove_argument (N, &argc, argv);
     594  }
     595
     596  if (argc != 1) usage_client ();
     597  return TRUE;
     598}
     599
    358600void usage () {
    359   fprintf (stderr, "ERROR: USAGE: relastro -region RA RA DEC DEC\n");
    360   fprintf (stderr, "       OR:    relastro -catalog (ra) (dec)\n\n");
    361   fprintf (stderr, "       OR:    relastro -high-speed [options]\n\n");
     601  fprintf (stderr, "ERROR: USAGE: relastro -update-simple [options]\n");
     602  fprintf (stderr, "       OR:    relastro -update-chips [options]\n");
     603  fprintf (stderr, "       OR:    relastro -update-mosaic [options]\n");
     604  fprintf (stderr, "       OR:    relastro -update-objects [options]\n");
     605  fprintf (stderr, "       OR:    relastro -high-speed [options]\n");
    362606  fprintf (stderr, "       OR:    relastro -merge-source [options]\n\n");
     607
    363608  fprintf (stderr, "  specify one of the following modes: \n");
    364609  fprintf (stderr, "  -update-objects\n");
    … …  
    369614  fprintf (stderr, "  -update-chips\n");
    370615  fprintf (stderr, "  -update-mosaics\n");
    371   fprintf (stderr, "  -high-speed (code[,code,code]) (code[,code,code]) (radius) (output catdir)\n\n");
    372   fprintf (stderr, "  additional options: \n");
     616  fprintf (stderr, "  -high-speed (code[,code,code]) (code[,code,code]) (radius) (output catdir)\n");
     617  fprintf (stderr, "  -merge-source (objID) (catID) into (objID) (catID)\n\n");
     618
     619  fprintf (stderr, " additional options: \n");
     620  fprintf (stderr, "  -region RA RA DEC DEC\n");
     621  fprintf (stderr, "  -catalog (ra) (dec)\n\n");
    373622  fprintf (stderr, "  -time (start)(stop)\n");
    374623  fprintf (stderr, "  +photcode (code)[,code,code...]\n");
    … …  
    383632  fprintf (stderr, "  -plrange\n");
    384633  fprintf (stderr, "  -minerror\n");
    385   fprintf (stderr, "  -area\n");
    386   fprintf (stderr, "  -area Xmin Xmax Ymin Ymax\n");
    387634  fprintf (stderr, "  -instmag min max\n\n");
    388635  fprintf (stderr, "  +photflags\n");
    … …  
    394641}
    395642
     643void usage_client () {
     644  fprintf (stderr, "ERROR: USAGE: relastro_client -load\n");
     645  fprintf (stderr, "       OR:    relastro_client -update-offsets\n");
     646  fprintf (stderr, "       OR:    relastro_client -update-objects\n");
     647  fprintf (stderr, "       OR:    relastro_client -high-speed\n\n");
     648
     649  fprintf (stderr, "  specify one of the following modes: \n");
     650  fprintf (stderr, "  -load : load the bright source detections for analysis\n");
     651  fprintf (stderr, "  -update-offsets : apply the updated image parameters\n");
     652  fprintf (stderr, "  -update-objects : calculate average astrometric properties of objects\n");
     653  fprintf (stderr, "  -high-speed (code[,code,code]) (code[,code,code]) (radius) (output catdir)\n");
     654  fprintf (stderr, "    support for the special high-speed mode\n");
     655
     656  fprintf (stderr, " additional options: \n");
     657  fprintf (stderr, "  -pm : calculate proper motions\n");
     658  fprintf (stderr, "  -par : calculate parallaxes\n");
     659  fprintf (stderr, "  -pmpar : calculate motions and parallaxes\n");
     660  fprintf (stderr, "  -region Rmin Rmax Dmin Dmax");
     661  fprintf (stderr, "  -catalog RA DEC");
     662  fprintf (stderr, "  -time (start)(stop)\n");
     663  fprintf (stderr, "  +photcode (code)[,code,code...]\n");
     664  fprintf (stderr, "  -photcode (code)[,code,code...]\n");
     665  fprintf (stderr, "  -statmode (mode)\n");
     666  fprintf (stderr, "  -reset\n");
     667  fprintf (stderr, "  -minerror\n");
     668  fprintf (stderr, "  -instmag min max\n\n");
     669  fprintf (stderr, "  +photflags\n");
     670  fprintf (stderr, "  -photflagbad\n");
     671  fprintf (stderr, "  -photflagpoor\n");
     672  fprintf (stderr, "  -v\n");
     673  fprintf (stderr, "  \n");
     674  exit (2);
     675}
     676
    396677void usage_merge_source_id (char *name) {
    397678

  • branches/eam_branches/ipp-20111122/Ohana/src/relastro/src/bcatalog.c

    r33422 r33448  
    7171      // allowed.
    7272
    73       CopyMeasureTiny (&subcatalog[0].measureT[Nmeasure], &catalog[0].measure[offset]);
     73      CopyMeasureToTiny (&subcatalog[0].measureT[Nmeasure], &catalog[0].measure[offset]);
    7474      // subcatalog[0].measure[Nmeasure] = catalog[0].measure[offset];
    7575      subcatalog[0].measureT[Nmeasure].dbFlags &= ~ID_MEAS_SKIP_ASTROM;
    … …  
    102102  subcatalog[0].Nmeasure = Nmeasure;
    103103  subcatalog[0].Nsecfilt = catalog[0].Nsecfilt;
    104   subcatalog[0].Nsecf_mem = Naverage * Nsecfilt;
     104  subcatalog[0].Nsecfilt = catalog[0].Nsecfilt;
     105  subcatalog[0].catID    = catalog[0].catID;
    105106  assert (Nsecfilt == catalog[0].Nsecfilt);
    106107

  • branches/eam_branches/ipp-20111122/Ohana/src/relastro/src/dvo_astrom_ops.c

    r12332 r33448  
    11# include "relastro.h"
    22
    3 double getMeanR (Measure *measure, Average *average, SecFilt *secfilt) {
     3double getMeanR (MeasureTiny *measure, Average *average, SecFilt *secfilt) {
    44
    55  double ra;
    … …  
    2121}
    2222
    23 double getMeanD (Measure *measure, Average *average, SecFilt *secfilt) {
     23double getMeanD (MeasureTiny *measure, Average *average, SecFilt *secfilt) {
    2424
    2525  double dec;
    … …  
    4141}
    4242
    43 int setMeanR (double ra_fit, Measure *measure, Average *average, SecFilt *secfilt) {
     43int setMeanR (double ra_fit, MeasureTiny *measure, Average *average, SecFilt *secfilt) {
    4444
    4545  /* math to get from new fitted position to new measure offset
    … …  
    6666}
    6767
    68 int setMeanD (double dec_fit, Measure *measure, Average *average, SecFilt *secfilt) {
     68int setMeanD (double dec_fit, MeasureTiny *measure, Average *average, SecFilt *secfilt) {
    6969
    7070  /* math to get from new fitted position to new measure offset
    … …  
    9090  return (TRUE);
    9191}
     92
     93int setMeanR_Big (double ra_fit, Measure *measure, Average *average, SecFilt *secfilt) {
     94
     95  if (!measure) return TRUE;
     96
     97  /* math to get from new fitted position to new measure offset
     98     ra_obs = average[0].R - measure[0].dR / 3600.0;
     99     measure[0].dR = (ra_fit - ra_obs) * 3600.0;
     100     measure[0].dR = (ra_fit - average[0].R + measure[0].dR/3600) * 3600.0
     101     measure[0].dR = (ra_fit - average[0].R) * 3600.0 + measure[0].dR;
     102  */
     103
     104  /* the measure carries the instantaneous mean position at the epoch t */
     105  measure[0].dR += (ra_fit - average[0].R) * 3600.0;
     106
     107  /* possible corrections to mean ra:
     108
     109  - proper-motion and parallax
     110  - abberation
     111  - precession and nutation, etc
     112  - refraction
     113  - DCR
     114
     115  */
     116
     117  return (TRUE);
     118}
     119
     120int setMeanD_Big (double dec_fit, Measure *measure, Average *average, SecFilt *secfilt) {
     121
     122  if (!measure) return TRUE;
     123
     124  /* math to get from new fitted position to new measure offset
     125     dec_obs = average[0].D - measure[0].dD / 3600.0;
     126     measure[0].dD = (dec_fit - dec_obs) * 3600.0;
     127     measure[0].dD = (dec_fit - average[0].D + measure[0].dD/3600) * 3600.0
     128     measure[0].dD = (dec_fit - average[0].D) * 3600.0 + measure[0].dD;
     129  */
     130
     131  /* the measure carries the instantaneous mean position at the epoch t */
     132  measure[0].dD += (dec_fit - average[0].D) * 3600.0;
     133
     134  /* possible corrections to mean ra:
     135
     136  - proper-motion and parallax
     137  - abberation
     138  - precession and nutation, etc
     139  - refraction
     140  - DCR
     141
     142  */
     143
     144  return (TRUE);
     145}

  • branches/eam_branches/ipp-20111122/Ohana/src/relastro/src/high_speed_catalogs.c

    r33396 r33448  
    11# include "relastro.h"
    22
    3 int high_speed_catalogs () {
     3int high_speed_catalogs_parallel (SkyList *sky);
     4
     5int high_speed_catalogs (SkyTable *sky, SkyList *skylist, int hostID, char *hostpath) {
    46
    57  int i;
    6 
    7   SkyTable *sky = NULL;
    8   SkyList *skylist = NULL;
    98  Catalog catalog;
    109
    … …  
    1211  dvo_set_catdir(CATDIR);
    1312
    14   // load the current sky table (layout of all SkyRegions)
    15   sky = SkyTableLoadOptimal (CATDIR, SKY_TABLE, GSCFILE, FALSE, SKY_DEPTH, VERBOSE);
    16   SkyTableSetFilenames (sky, CATDIR, "cpt");
    17  
    18   // determine the populated SkyRegions overlapping the requested area (default depth)
    19   if (UserCatalog) {
    20     skylist = SkyRegionByPoint (sky, -1, UserCatalogRA, UserCatalogDEC);
    21     setupAreaSelection(skylist->regions[0]);
    22   } else {
    23     skylist = SkyListByPatch (sky, -1, &UserPatch);
    24     setupAreaSelection(&UserPatch);
    25   }
    26 
    2713  initializeConstraints();
    2814
    2915  // XXX need to decide how to determine PARALLEL mode...
    3016  if (PARALLEL & !hostID) {
    31     high_speed_catalogs (skylist);
    32     return TRUE;
     17    high_speed_catalogs_parallel (skylist);
     18    goto finish;
    3319  }
    3420
    … …  
    6450  }
    6551 
     52finish: {
     53    char photcodeFile[1024];
     54    snprintf (photcodeFile, 1024, "%s/Photcodes.dat", HIGH_SPEED_DIR);
     55    if (!SavePhotcodesFITS (photcodeFile)) {
     56      fprintf (stderr, "error saving photcode table %s\n", photcodeFile);
     57      exit (1);
     58    }
     59
     60    // need to copy across the Images, SkyTable, and Photcode tables:
     61    char *skyfile = SkyTableFilename (HIGH_SPEED_DIR);
     62    SkyTableSave (sky, skyfile);
     63
     64    char line[2048];
     65    snprintf (line, 2048, "cp %s/Images.dat %s/Images.dat", CATDIR, HIGH_SPEED_DIR);
     66    system (line);
     67  }
     68
    6669  return (TRUE);
    6770}
    … …  
    6972// CATDIR is supplied globally
    7073# define DEBUG 1
    71 int high_speed_catalogs_parallel (SkyList *sky) {
     74int high_speed_catalogs_parallel (SkyList *skylist) {
    7275
    7376  // launch the setphot_client jobs to the parallel hosts
    7477
    7578  // load the list of hosts
    76   HostTable *table = HostTableLoad (CATDIR, sky->hosts);
     79  HostTable *table = HostTableLoad (CATDIR, skylist->hosts);
    7780  if (!table) {
    78     fprintf (stderr, "ERROR: failure reading Host Table %s for database %s\n", sky->hosts, CATDIR);
     81    fprintf (stderr, "ERROR: failure reading Host Table %s for database %s\n", skylist->hosts, CATDIR);
    7982    exit (1);
    8083  }   
    … …  
    8891    table->hosts[i].pathname = tmppath;
    8992
    90     char catalogFile[512];
    91     snprintf (catalogFile, 512, "%s/relastro.catalog.subset.dat", table->hosts[i].pathname);
    92 
    93     // options / arguments that can affect relastro_client -update-objects:
     93    // options / arguments that can affect relastro_client -high-speed
    9494
    9595    char command[1024];
    96     snprintf (command, 1024, "relastro_client -high-speed -hostID %d -D CATDIR %s -hostdir %s -region %f %f %f %f",
    97               table->hosts[i].hostID, CATDIR, table->hosts[i].pathname, UserPatch.Rmin, UserPatch.Rmax, UserPatch.Dmin, UserPatch.Dmax);
     96    snprintf (command, 1024, "relastro_client -high-speed %s %s %f %s -hostID %d -D CATDIR %s -hostdir %s -region %f %f %f %f",
     97              PHOTCODE_A_LIST, PHOTCODE_B_LIST, RADIUS, HIGH_SPEED_DIR, table->hosts[i].hostID, CATDIR, table->hosts[i].pathname, UserPatch.Rmin, UserPatch.Rmax, UserPatch.Dmin, UserPatch.Dmax);
    9898
    9999    char tmpline[1024];
    100     if (VERBOSE)       { snprintf (tmpline, 1024, "%s -v",              command);                    strcpy (command, tmpline); }
    101     if (VERBOSE2)      { snprintf (tmpline, 1024, "%s -vv",             command);                    strcpy (command, tmpline); }
    102     if (RESET)         { snprintf (tmpline, 1024, "%s -reset",          command);                    strcpy (command, tmpline); }
    103     if (KEEP_UBERCAL)  { snprintf (tmpline, 1024, "%s -keep-ubercal",   command);                    strcpy (command, tmpline); }
     100    if (FIT_MODE == FIT_PM_ONLY)         { snprintf (tmpline, 1024, "%s -pm",    command);              strcpy (command, tmpline); }
     101    if (FIT_MODE == FIT_PAR_ONLY)        { snprintf (tmpline, 1024, "%s -par",   command);              strcpy (command, tmpline); }
     102    if (FIT_MODE == FIT_PM_AND_PAR)      { snprintf (tmpline, 1024, "%s -pmpar", command);              strcpy (command, tmpline); }
     103
     104    if (VERBOSE)       { snprintf (tmpline, 1024, "%s -v",              command);                       strcpy (command, tmpline); }
     105    if (VERBOSE2)      { snprintf (tmpline, 1024, "%s -vv",             command);                       strcpy (command, tmpline); }
     106    if (RESET)         { snprintf (tmpline, 1024, "%s -reset",          command);                       strcpy (command, tmpline); }
     107    if (ImagSelect)    { snprintf (tmpline, 1024, "%s -instmag %f %f",  command, ImagMin, ImagMax);     strcpy (command, tmpline); }
     108    if (MaxDensityUse) { snprintf (tmpline, 1024, "%s -max-density %f", command, MaxDensityValue);      strcpy (command, tmpline); }
     109   
     110    if (USE_BASIC_CHECK) { snprintf (tmpline, 1024, "%s -basic-image-search", command);                 strcpy (command, tmpline); }
     111
     112    if (FlagOutlier)   { snprintf (tmpline, 1024, "%s -clip %d",        command, CLIP_THRESH);          strcpy (command, tmpline); }
     113   
     114    if (PHOTCODE_KEEP_LIST) { snprintf (tmpline, 1024, "%s +photcode %s", command, PHOTCODE_KEEP_LIST); strcpy (command, tmpline); }
     115    if (PHOTCODE_SKIP_LIST) { snprintf (tmpline, 1024, "%s -photcode %s", command, PHOTCODE_SKIP_LIST); strcpy (command, tmpline); }
     116    if (PhotFlagSelect)    { snprintf (tmpline, 1024, "%s +photflags",   command);                      strcpy (command, tmpline); }
     117    if (PhotFlagBad)       { snprintf (tmpline, 1024, "%s +photflagbad %d", command, PhotFlagBad);      strcpy (command, tmpline); }
     118    if (PhotFlagPoor)      { snprintf (tmpline, 1024, "%s +photflagpoor %d", command, PhotFlagPoor);    strcpy (command, tmpline); }
     119    // XXX note that the above pass in the flag as decimal -- also note that args.c cannot handle 0xHEX values
     120
     121    if (WHERE_A[0]) { snprintf (tmpline, 1024, "%s -D WHERE_A \"%s\"", command, WHERE_A);               strcpy (command, tmpline); }
     122    if (WHERE_B[0]) { snprintf (tmpline, 1024, "%s -D WHERE_B \"%s\"", command, WHERE_B);               strcpy (command, tmpline); }
     123
     124    if (TimeSelect) {
     125      char *tstart = ohana_sec_to_date (TSTART);
     126      char *tstop  = ohana_sec_to_date (TSTOP);
     127      snprintf (tmpline, 1024, "%s -time %s %s", command, tstart, tstop);
     128      free (tstart);
     129      free (tstop);
     130      strcpy (command, tmpline);
     131    }
    104132
    105133    fprintf (stderr, "command: %s\n", command);

  • branches/eam_branches/ipp-20111122/Ohana/src/relastro/src/high_speed_objects.c

    r31635 r33448  
    66int high_speed_objects (SkyRegion *region, Catalog *catalog) {
    77
    8   // XXX seems silly to require region here just to find the catalog bounds / center
    9 
    10   off_t i, j, m, J, ni, nj, *N1, Nslow, Ninvalid, NgroupA, NgroupB, NgroupAbad, NgroupBbad,l,i1,Noff,Nmeasure,Naverage, NAVERAGE, NMEASURE;
    11   int *slowMoving, *groupA, *groupB, status, Nmatch, Nmatchmeas,Nepoch,nv[2],Nmatchmeasobj;
     8  off_t i, j, m, J, ni, nj, *N1;
     9  off_t Nslow, Ninvalid, NgroupA, NgroupB, NgroupAbad, NgroupBbad;
     10  off_t l, i1, NAVERAGE, NMEASURE;
     11  int *slowMoving, *groupA, *groupB, status, Nmatch, Nmatchmeas, Nepoch, nv[2], Nmatchmeasobj;
    1212  int foundA, foundB;
    1313  double *X1, *Y1;
    1414  double dX, dY, dR, RADIUS2;
    1515  Coords tcoords;
    16   Catalog catalog1;
     16  Catalog catalogOut;
    1717
    1818  int Nsecfilt;
    1919  char filename[1024];
    20   Noff = strlen(CATDIR);
    21   sprintf (filename, "%s%s", HIGH_SPEED_DIR, &catalog[0].filename[Noff]);
    22   fprintf(stderr, "%s\n",filename);
    23   dvo_catalog_init(&catalog1, TRUE); /*initialise new catalogue*/
    24   catalog1.filename = strcreate(filename);
     20  snprintf (filename, 1024, "%s/%s.cpt", HIGH_SPEED_DIR, region[0].name);
     21  fprintf (stderr, "%s\n",filename);
     22
     23  dvo_catalog_init (&catalogOut, TRUE); /*initialise new catalogue*/
     24  catalogOut.filename = strcreate(filename);
    2525
    2626  SIGMA_LIM = 0.0;
    2727
    2828  Nsecfilt = GetPhotcodeNsecfilt();
    29   Naverage=catalog[0].Naverage;
    30   Nmeasure=catalog[0].Nmeasure;
    31   catalog1.filename = filename; // based on the input name, need to keep everything below the catdir portion
    32   catalog1.Nsecfilt  = Nsecfilt;
    33   // always load all of the data (if any exists)
    34   catalog1.catflags = LOAD_AVES | LOAD_MEAS | LOAD_MISS | LOAD_SECF;
     29  // off_t Naverage = catalog[0].Naverage;
     30  // off_t Nmeasure = catalog[0].Nmeasure;
     31  catalogOut.filename = filename; // based on the input name, need to keep everything below the catdir portion
     32  catalogOut.Nsecfilt = Nsecfilt;
     33  catalogOut.catflags = LOAD_AVES | LOAD_MEAS | LOAD_MISS | LOAD_SECF; // load all data
    3534 
    36   catalog1.catformat = dvo_catalog_catformat (CATFORMAT);  // set the default catformat from config data
    37   catalog1.catmode   = dvo_catalog_catmode (CATMODE);      // set the default catmode from config data
    38 
    39   if (!dvo_catalog_open (&catalog1, region, VERBOSE,"w")) {
     35  catalogOut.catformat = dvo_catalog_catformat (CATFORMAT);  // set the default catformat from config data
     36  catalogOut.catmode   = dvo_catalog_catmode (CATMODE);      // set the default catmode from config data
     37
     38  if (!dvo_catalog_open (&catalogOut, region, VERBOSE,"w")) {
    4039    fprintf (stderr, "ERROR: failure to open catalog file %s\n",
    4140             filename);
    … …  
    4544  NAVERAGE = 1000;
    4645  NMEASURE = 10000;
    47   REALLOCATE (catalog1.average, Average, NAVERAGE);
    48   REALLOCATE (catalog1.measure, Measure, NMEASURE);
    49   REALLOCATE(catalog1.secfilt, SecFilt, NAVERAGE*Nsecfilt);
     46  REALLOCATE (catalogOut.average, Average, NAVERAGE);
     47  REALLOCATE (catalogOut.measure, Measure, NMEASURE);
     48  REALLOCATE (catalogOut.secfilt, SecFilt, NAVERAGE*Nsecfilt);
    5049
    5150  // high-speed between different surveys (easier case):
    … …  
    114113            // fprintf (stderr, "skip "OFF_T_FMT" (group A)\n",  i);
    115114        }
    116 
    117     } else if (foundB && !foundA) {
    118 
     115        continue;
     116    }
     117    if (foundB && !foundA) {
    119118        if (applyConstraintsB(catalog, i)) {
    120119            groupB[i] = TRUE;
    … …  
    124123            NgroupBbad++;
    125124        }
     125        continue;
    126126    }
    127127    // this object does not have a detection matching the contraints from either gorupA or groupB -- skip as if slow
    … …  
    209209      nv[1]=nj;
    210210
    211       catalog1.average[Nmatch]=catalog[0].average[nv[0]];
     211      catalogOut.average[Nmatch]=catalog[0].average[nv[0]];
    212212      /*Loop over index values and set measurements in new catalogue*/
    213213      Nmatchmeasobj=0;
    214       catalog1.average[Nmatch].measureOffset=Nmatchmeas;
     214      catalogOut.average[Nmatch].measureOffset=Nmatchmeas;
    215215      for(l=0;l<Nepoch;l++) {
    216216          m = catalog[0].average[nv[l]].measureOffset;
    217217          for(i1=0;i1<catalog[0].average[nv[l]].Nmeasure;i1++) {
    218               catalog1.measure[Nmatchmeas]=catalog[0].measure[m+i1];
     218              catalogOut.measure[Nmatchmeas]=catalog[0].measure[m+i1];
    219219              /*Set offset RA and Dec wrt correct average value*/
    220               catalog1.measure[Nmatchmeas].dR=catalog[0].measure[m+i1].dR+3600.0*(catalog[0].average[nv[0]].R-catalog[0].average[nv[l]].R);
    221               catalog1.measure[Nmatchmeas].dD=catalog[0].measure[m+i1].dD+3600.0*(catalog[0].average[nv[0]].D-catalog[0].average[nv[l]].D);
    222               catalog1.measure[Nmatchmeas].averef = Nmatch;
     220              catalogOut.measure[Nmatchmeas].dR=catalog[0].measure[m+i1].dR+3600.0*(catalog[0].average[nv[0]].R-catalog[0].average[nv[l]].R);
     221              catalogOut.measure[Nmatchmeas].dD=catalog[0].measure[m+i1].dD+3600.0*(catalog[0].average[nv[0]].D-catalog[0].average[nv[l]].D);
     222              catalogOut.measure[Nmatchmeas].averef = Nmatch;
    223223              Nmatchmeasobj++;
    224224              Nmatchmeas++;
    … …  
    226226              if (Nmatchmeas == NMEASURE - 1) {
    227227                NMEASURE += 10000;
    228                 REALLOCATE (catalog1.measure, Measure, NMEASURE);
     228                REALLOCATE (catalogOut.measure, Measure, NMEASURE);
    229229              }
    230230          }
    231231      }
    232       catalog1.average[Nmatch].Nmeasure=Nmatchmeasobj;
     232      catalogOut.average[Nmatch].Nmeasure = Nmatchmeasobj;
    233233      Nmatch ++;
    234234
    235235      if (Nmatch == NAVERAGE - 1) {
    236236        NAVERAGE += 1000;
    237         REALLOCATE (catalog1.average, Average, NAVERAGE);
    238         REALLOCATE(catalog1.secfilt, SecFilt, NAVERAGE*Nsecfilt);
     237        REALLOCATE (catalogOut.average, Average, NAVERAGE);
     238        REALLOCATE(catalogOut.secfilt, SecFilt, NAVERAGE*Nsecfilt);
    239239      }
    240240    }
    241241    i++;
    242242  }
    243   catalog1.Naverage=Nmatch;
    244   catalog1.Nmeasure=Nmatchmeas;
    245   catalog1.Nsecfilt=Nsecfilt;
    246   catalog1.Nsecf_mem=Nmatch*Nsecfilt;
    247 
    248   UpdateObjects(&catalog1, 1);
     243  catalogOut.Naverage=Nmatch;
     244  catalogOut.Nmeasure=Nmatchmeas;
     245  catalogOut.Nsecfilt=Nsecfilt;
     246  catalogOut.Nsecf_mem=Nmatch*Nsecfilt;
     247
     248  populate_tiny_values (&catalogOut, DVO_TV_MEASURE);
     249  UpdateObjects (&catalogOut, 1);
    249250
    250251  fprintf (stderr, "found %d matches\n", Nmatch);
    251   dvo_catalog_save (&catalog1, VERBOSE);
    252   dvo_catalog_unlock (&catalog1);
    253   dvo_catalog_free (&catalog1);
     252  dvo_catalog_save (&catalogOut, VERBOSE);
     253  dvo_catalog_unlock (&catalogOut);
     254  dvo_catalog_free (&catalogOut);
    254255  free (slowMoving);
    255256  free (groupA);

  • branches/eam_branches/ipp-20111122/Ohana/src/relastro/src/high_speed_utils.c

    r33099 r33448  
    1717static dbStack *stackA = NULL;
    1818static int NstackA = 0;
     19
    1920static dbField *fieldsB = NULL;
    2021static int NfieldsB = 0;
    2122static dbStack *stackB = NULL;
    2223static int NstackB = 0;
     24
    2325static dbValue *valuesA;
    2426// static int NvaluesA;
    … …  
    6163int setupConstraints(char *whereString, dbField **pFields, int *pNfields, dbStack **pStack, int *pNstack) {
    6264
     65  // an empty where is allowed (no further filtering)
     66  if (!whereString[0]) return TRUE;
     67
    6368    // split the string into words
    6469    int nWords;
    … …  
    8994
    9095int initializeConstraints() {
     96
    9197  // if constaints haven't been set by the user use 2mass
    92     if (!WHERE_A[0]) {
    93         strcpy(WHERE_A, DEFAULT_WHERE_A);
    94     }
    95     if (!WHERE_B[0]) {
    96         strcpy(WHERE_B, DEFAULT_WHERE_B);
    97     }
     98  // XXX  if (!WHERE_A[0]) {
     99  // XXX      strcpy(WHERE_A, DEFAULT_WHERE_A);
     100  // XXX  }
     101  // XXX  if (!WHERE_B[0]) {
     102  // XXX      strcpy(WHERE_B, DEFAULT_WHERE_B);
     103  // XXX  }
     104
     105    fprintf (stderr, "where A: %s\n", WHERE_A);
     106    fprintf (stderr, "where B: %s\n", WHERE_B);
     107
     108    // pull TIMEFORMAT & TIMEREF from the Config system
     109    dbExtractAveragesInit();
    98110
    99111    if (! setupConstraints(WHERE_A, &fieldsA, &NfieldsA, &stackA, &NstackA)) {
    … …  
    102114        exit(1);
    103115    }
    104     ALLOCATE(valuesA, dbValue, NfieldsA);
     116    if (NfieldsA) {
     117      ALLOCATE(valuesA, dbValue, NfieldsA);
     118    }
     119
    105120    if (! setupConstraints(WHERE_B, &fieldsB, &NfieldsB, &stackB, &NstackB)) {
    106121        print_error();
    … …  
    108123        exit(1);
    109124    }
    110     ALLOCATE(valuesB, dbValue, NfieldsB);
     125    if (NfieldsB) {
     126      ALLOCATE(valuesB, dbValue, NfieldsB);
     127    }
     128
    111129    return TRUE;
    112130}
    … …  
    116134int applyConstraintsA(Catalog *catalog, off_t i) {
    117135    off_t m = catalog[0].average[i].measureOffset;
    118     off_t n;
    119136    int Nsecfilt = GetPhotcodeNsecfilt();
    120137   
    121     dbExtractAveragesInit();
     138    // an empty WHERE matches all objects
     139    if (!NfieldsA) return TRUE;
     140
     141    off_t n;
    122142    for (n = 0; n < NfieldsA; n++) {
    123143        valuesA[n] = dbExtractAverages (&catalog[0].average[i], &catalog[0].secfilt[i*Nsecfilt], &catalog[0].measure[m], &fieldsA[n]);
    … …  
    130150    int Nsecfilt = GetPhotcodeNsecfilt();
    131151   
    132     dbExtractAveragesInit();
     152    // an empty WHERE matches all objects
     153    if (!NfieldsB) return TRUE;
    133154
    134155    off_t n;

  • branches/eam_branches/ipp-20111122/Ohana/src/relastro/src/initialize.c

    r32346 r33448  
    123123}
    124124
     125void initialize_client (int argc, char **argv) {
     126
     127  int NPHOTCODES;
     128  char *codename, *ptr, *list;
     129
     130  ptr = NULL;
     131
     132  ConfigInit (&argc, argv);
     133  args_client (argc, argv);
     134
     135  /* build a list of accepted photcodes. these will be used by bcatalog to accept or
     136     reject loaded data */
     137
     138  NphotcodesKeep = 0;
     139  photcodesKeep = NULL;
     140  if (PHOTCODE_KEEP_LIST != NULL) {
     141    NPHOTCODES = 10;
     142    ALLOCATE (photcodesKeep, PhotCode *, NPHOTCODES);
     143
     144    /* parse the comma-separated list of photcodesKeep */
     145    list = PHOTCODE_KEEP_LIST;
     146    while ((codename = strtok_r (list, ",", &ptr)) != NULL) {
     147      list = NULL; // pass NULL on successive strtok_r calls
     148      fprintf (stderr, "PHOTCODE_LIST: %s\n", PHOTCODE_KEEP_LIST);
     149      fprintf (stderr, "codename: %s\n", codename);
     150      if ((photcodesKeep[NphotcodesKeep] = GetPhotcodebyName (codename)) == NULL) {
     151        fprintf (stderr, "ERROR: photcode %s not found in photcode table\n", codename);
     152        exit (1);
     153      }
     154      NphotcodesKeep ++;
     155      CHECK_REALLOCATE (photcodesKeep, PhotCode *, NPHOTCODES, NphotcodesKeep, 10);
     156    }
     157  }
     158
     159  NphotcodesSkip = 0;
     160  photcodesSkip = NULL;
     161  if (PHOTCODE_SKIP_LIST != NULL) {
     162    NPHOTCODES = 10;
     163    ALLOCATE (photcodesSkip, PhotCode *, NPHOTCODES);
     164
     165    /* parse the comma-separated list of photcodesSkip */
     166    list = PHOTCODE_SKIP_LIST;
     167    while ((codename = strtok_r (list, ",", &ptr)) != NULL) {
     168      list = NULL; // pass NULL on successive strtok_r calls
     169      fprintf (stderr, "PHOTCODE_LIST: %s\n", PHOTCODE_SKIP_LIST);
     170      fprintf (stderr, "codename: %s\n", codename);
     171      if ((photcodesSkip[NphotcodesSkip] = GetPhotcodebyName (codename)) == NULL) {
     172        fprintf (stderr, "ERROR: photcode %s not found in photcode table\n", codename);
     173        exit (1);
     174      }
     175      NphotcodesSkip ++;
     176      CHECK_REALLOCATE (photcodesSkip, PhotCode *, NPHOTCODES, NphotcodesSkip, 10);
     177    }
     178  }
     179
     180  NphotcodesGroupA = 0;
     181  photcodesGroupA = NULL;
     182  if (PHOTCODE_A_LIST != NULL) {
     183    NPHOTCODES = 10;
     184    ALLOCATE (photcodesGroupA, PhotCode *, NPHOTCODES);
     185
     186    /* parse the comma-separated list of photcodesGroupA */
     187    fprintf (stderr, "PHOTCODE_A_LIST: %s\n", PHOTCODE_A_LIST);
     188    list = PHOTCODE_A_LIST;
     189    while ((codename = strtok_r (list, ",", &ptr)) != NULL) {
     190      list = NULL; // pass NULL on successive strtok_r calls
     191      fprintf (stderr, "codename: %s\n", codename);
     192      if ((photcodesGroupA[NphotcodesGroupA] = GetPhotcodebyName (codename)) == NULL) {
     193        fprintf (stderr, "ERROR: photcode %s not found in photcode table\n", codename);
     194        exit (1);
     195      }
     196      NphotcodesGroupA ++;
     197      CHECK_REALLOCATE (photcodesGroupA, PhotCode *, NPHOTCODES, NphotcodesGroupA, 10);
     198    }
     199  }
     200
     201  NphotcodesGroupB = 0;
     202  photcodesGroupB = NULL;
     203  if (PHOTCODE_B_LIST != NULL) {
     204    NPHOTCODES = 10;
     205    ALLOCATE (photcodesGroupB, PhotCode *, NPHOTCODES);
     206
     207    fprintf (stderr, "PHOTCODE_B_LIST: %s\n", PHOTCODE_B_LIST);
     208    /* parse the comma-separated list of photcodesGroupB */
     209    list = PHOTCODE_B_LIST;
     210    while ((codename = strtok_r (list, ",", &ptr)) != NULL) {
     211      list = NULL; // pass NULL on successive strtok_r calls
     212      fprintf (stderr, "codename: %s\n", codename);
     213      if ((photcodesGroupB[NphotcodesGroupB] = GetPhotcodebyName (codename)) == NULL) {
     214        fprintf (stderr, "ERROR: photcode %s not found in photcode table\n", codename);
     215        exit (1);
     216      }
     217      NphotcodesGroupB ++;
     218      CHECK_REALLOCATE (photcodesGroupB, PhotCode *, NPHOTCODES, NphotcodesGroupB, 10);
     219    }
     220  }
     221
     222  initstats (STATMODE);
     223
     224  // IMAGE_BAD = ID_IMAGE_ASTROM_POOR | ID_IMAGE_ASTROM_FEW | ID_IMAGE_ASTROM_SKIP;
     225  // STAR_BAD  = ID_STAR_POOR | ID_STAR_FEW;
     226  // MEAS_BAD  = ID_MEAS_NOCAL | ID_MEAS_POOR_ASTROM | ID_MEAS_SKIP_ASTROM | ID_MEAS_AREA;
     227}
     228

  • branches/eam_branches/ipp-20111122/Ohana/src/relastro/src/load_catalogs.c

    r33396 r33448  
    11# include "relastro.h"
    22
    3 Catalog *load_catalogs (SkyList *skylist, int *Ncatalog, int subselect) {
     3Catalog *load_catalogs_parallel (SkyList *sky, int *Ncatalog);
     4
     5Catalog *load_catalogs (SkyList *skylist, int *Ncatalog, int subselect, int hostID, char *hostpath) {
    46
    57  int i, j, k, m, Nstar;
    … …  
    1618  ALLOCATE (catalog, Catalog, skylist[0].Nregions);
    1719
    18   initStarMaps();
    19 
    2020  // load data from each region file, only use bright stars
     21  int Ncat = 0;
    2122  for (i = 0; i < skylist[0].Nregions; i++) {
    2223
    … …  
    3839    pcatalog->Nsecfilt  = GetPhotcodeNsecfilt ();
    3940
     41    // loads Average, Measure, SecFilt
    4042    if (!dvo_catalog_open (pcatalog, skylist[0].regions[i], VERBOSE2, "w")) {
    4143      fprintf (stderr, "ERROR: failure reading catalog %s\n", pcatalog[0].filename);
    4244      exit (1);
    4345    }
    44     if (VERBOSE2 && !pcatalog[0].Naves_disk) fprintf (stderr, "no data in %s, skipping\n", pcatalog[0].filename);
    45 
    46     //outlier rejection
     46    if (!pcatalog[0].Naves_disk) {
     47      if (VERBOSE2) fprintf (stderr, "no data in %s, skipping\n", pcatalog[0].filename);
     48      dvo_catalog_unlock (pcatalog);
     49      dvo_catalog_free (pcatalog);
     50      continue;
     51    }
     52
     53    if (!pcatalog->sorted) {
     54      fprintf (stderr, "this database is not sorted.  please sort using addstar -resort\n");
     55      exit (1);
     56    }
     57
     58    // outlier rejection (operates on full tables (Average, Measure, Secfilt))
    4759    if (FlagOutlier) {
    4860      FlagOutliers(&tcatalog);
    4961    }
    5062
    51     updateStarMaps (&tcatalog);
    52 
    5363    // select only the brighter stars
    5464    if (subselect) {
    55       bcatalog (&catalog[i], &tcatalog);
     65      // results are in Average, Secfilt, and MeasureTiny
     66      // Ncat tracks the actually used catalogs
     67      bcatalog (&catalog[Ncat], &tcatalog);
    5668      dvo_catalog_unlock (&tcatalog);
    5769      dvo_catalog_free (&tcatalog);
    5870    } else {
    5971      if (RESET) {
    60         for (j = 0; j < catalog[i].Naverage; j++) {
    61           catalog[i].average[j].flags = 0;
    62           m = catalog[i].average[j].measureOffset;
    63           for (k = 0; k < catalog[i].average[j].Nmeasure; k++) {
    64             catalog[i].measure[m+k].dbFlags = 0;
     72        for (j = 0; j < catalog[Ncat].Naverage; j++) {
     73          catalog[Ncat].average[j].flags = 0;
     74          m = catalog[Ncat].average[j].measureOffset;
     75          for (k = 0; k < catalog[Ncat].average[j].Nmeasure; k++) {
     76            catalog[Ncat].measure[m+k].dbFlags = 0;
    6577          }
    6678        }
    6779      }
    6880    }
     81    Ncat ++;
    6982  }
    7083
    7184  // XXX keep this test?
    7285  Nstar = 0;
    73   for (i = 0; i < skylist[0].Nregions; i++) {
     86  for (i = 0; i < Ncat; i++) {
    7487    Nstar += catalog[i].Naverage;
    7588  }
    … …  
    7891  }
    7992
    80   // XXX consider only returning the populated catalogs
    81   *Ncatalog = skylist[0].Nregions;
     93  // only return the populated catalogs
     94  REALLOCATE (catalog, Catalog, Ncat);
     95  *Ncatalog = Ncat;
    8296  return (catalog);
    8397}
    … …  
    92106# define DEBUG 1
    93107Catalog *load_catalogs_parallel (SkyList *sky, int *Ncatalog) {
    94 
    95   int Nsecfilt  = GetPhotcodeNsecfilt ();               // set the desired number in case we need to create the catalog
    96 
    97   // launch the setphot_client jobs to the parallel hosts
    98108
    99109  // load the list of hosts
    … …  
    117127    // options / arguments that can affect relastro_client -load:
    118128    // VERBOSE, VERBOSE2
    119     // KEEP_UBERCAL
    120129    // RESET (-reset)
    121130    // TimeSelect -time
    122     // DophotSelect
    123131    // (note that psfQual is applied rigidly at 0.85, as is the galaxy test)
    124     // MAG_LIM
    125     // SIGMA_LIM
    126132    // ImagSelect, ImagMin, ImagMax
    127133    // MaxDensityUse, MaxDensityValue
    128134
    129135    char command[1024];
    130     snprintf (command, 1024, "relastro_client %s -load %s -hostID %d -D CATDIR %s -hostdir %s -region %f %f %f %f -D CAMERA %s -D MAG_LIM %f -D SIGMA_LIM %f",
    131               PhotcodeList, catalogFile, table->hosts[i].hostID, CATDIR, table->hosts[i].pathname, UserPatch.Rmin, UserPatch.Rmax, UserPatch.Dmin, UserPatch.Dmax, CAMERA, MAG_LIM, SIGMA_LIM);
     136    snprintf (command, 1024, "relastro_client -load-objects %s -hostID %d -D CATDIR %s -hostdir %s -region %f %f %f %f -statmode %s -minerror %f",
     137              catalogFile, table->hosts[i].hostID, CATDIR, table->hosts[i].pathname, UserPatch.Rmin, UserPatch.Rmax, UserPatch.Dmin, UserPatch.Dmax, STATMODE, MIN_ERROR);
    132138
    133139    char tmpline[1024];
     140    if (FIT_MODE == FIT_PM_ONLY)         { snprintf (tmpline, 1024, "%s -pm",    command);           strcpy (command, tmpline); }
     141    if (FIT_MODE == FIT_PAR_ONLY)        { snprintf (tmpline, 1024, "%s -par",   command);           strcpy (command, tmpline); }
     142    if (FIT_MODE == FIT_PM_AND_PAR)      { snprintf (tmpline, 1024, "%s -pmpar", command);           strcpy (command, tmpline); }
     143
    134144    if (VERBOSE)       { snprintf (tmpline, 1024, "%s -v",              command);                    strcpy (command, tmpline); }
    135145    if (VERBOSE2)      { snprintf (tmpline, 1024, "%s -vv",             command);                    strcpy (command, tmpline); }
    136146    if (RESET)         { snprintf (tmpline, 1024, "%s -reset",          command);                    strcpy (command, tmpline); }
    137     if (KEEP_UBERCAL)  { snprintf (tmpline, 1024, "%s -keep-ubercal",   command);                    strcpy (command, tmpline); }
    138     if (DophotSelect)  { snprintf (tmpline, 1024, "%s -dophot %d",      command, DophotValue);       strcpy (command, tmpline); }
    139147    if (ImagSelect)    { snprintf (tmpline, 1024, "%s -instmag %f %f",  command, ImagMin, ImagMax);  strcpy (command, tmpline); }
    140148    if (MaxDensityUse) { snprintf (tmpline, 1024, "%s -max-density %f", command, MaxDensityValue);   strcpy (command, tmpline); }
     149   
     150    if (FlagOutlier)   { snprintf (tmpline, 1024, "%s -clip %d",        command, CLIP_THRESH);       strcpy (command, tmpline); }
     151   
     152    if (PHOTCODE_KEEP_LIST) { snprintf (tmpline, 1024, "%s +photcode %s", command, PHOTCODE_KEEP_LIST); strcpy (command, tmpline); }
     153    if (PHOTCODE_SKIP_LIST) { snprintf (tmpline, 1024, "%s -photcode %s", command, PHOTCODE_SKIP_LIST); strcpy (command, tmpline); }
     154    if (PhotFlagSelect)    { snprintf (tmpline, 1024, "%s +photflags",   command);                     strcpy (command, tmpline); }
     155    if (PhotFlagBad)       { snprintf (tmpline, 1024, "%s +photflagbad %d", command, PhotFlagBad);     strcpy (command, tmpline); }
     156    if (PhotFlagPoor)      { snprintf (tmpline, 1024, "%s +photflagpoor %d", command, PhotFlagPoor);   strcpy (command, tmpline); }
     157    // XXX note that the above pass in the flag as decimal -- also note that args.c cannot handle 0xHEX values
     158
    141159    if (TimeSelect) {
    142160      char *tstart = ohana_sec_to_date (TSTART);
    … …  
    183201
    184202  // set up an initial array of catalogs
     203  int Nsecfilt  = GetPhotcodeNsecfilt ();
    185204  CatalogSplitter *catalogs = BrightCatalogSplitInit (Nsecfilt);
    186205

  • branches/eam_branches/ipp-20111122/Ohana/src/relastro/src/load_images.c

    r32346 r33448  
    77  fprintf (stderr, MSG, __VA_ARGS__); }
    88
    9 SkyList *load_images (FITS_DB *db, SkyRegion *region) {
     9int load_images (FITS_DB *db, SkyList *skylist) {
    1010
    1111  Image     *image, *subset;
    … …  
    1414  struct timeval start, stop;
    1515
    16   SkyTable *sky = NULL;
    17   SkyList *skylist = NULL;
    18 
    1916  gettimeofday (&start, (void *) NULL);
    20 
    21   // load the current sky table (layout of all SkyRegions)
    22   sky = SkyTableLoadOptimal (CATDIR, SKY_TABLE, GSCFILE, TRUE, SKY_DEPTH, VERBOSE);
    23   SkyTableSetFilenames (sky, CATDIR, "cpt");
    24  
    25   // determine the populated SkyRegions overlapping the requested area
    26   if (UserCatalog) {
    27     skylist = SkyRegionByPoint (sky, -1, UserCatalogRA, UserCatalogDEC);
    28   } else {
    29     skylist = SkyListByPatch (sky, -1, region);
    30   }
    31   MARKTIME("  setup sky: %f sec\n", dtime);
    3217
    3318  // convert database table to internal structure
    … …  
    3722      exit (2);
    3823  }
    39 
    4024  MARKTIME("  convert image table: %f sec\n", dtime);
    4125
    … …  
    5337  if (!UPDATE) dvo_image_unlock (db);
    5438
    55   return (skylist);
     39  return TRUE;
    5640}
    5741
    … …  
    7155    memcpy (vtable[0].buffer[i], &image[i], Nx);
    7256  }
    73   return (TRUE);
     57  return TRUE;
    7458}

  • branches/eam_branches/ipp-20111122/Ohana/src/relastro/src/relastro.c

    r33396 r33448  
    33int main (int argc, char **argv) {
    44
    5   int i, status, Ncatalog;
    6   Catalog *catalog;
    7   FITS_DB db;
    8   struct timeval start, stop;
    9   SkyList *skylist = NULL;
    10 
    11   gettimeofday (&start, (void *) NULL);
    12 
    135  /* get configuration info, args */
    146  initialize (argc, argv);
    157
    16   /* the object analysis is a separate process iterating over catalogs */
    17   if (FIT_TARGET == TARGET_OBJECTS) {
    18     relastro_objects ();
    19     exit (0);
    20   }
     8  SkyTable *sky = SkyTableLoadOptimal (CATDIR, SKY_TABLE, GSCFILE, TRUE, SKY_DEPTH, VERBOSE);
     9  SkyTableSetFilenames (sky, CATDIR, "cpt");
     10  SkyList *skylist = SkyListByPatch (sky, -1, &UserPatch);
    2111
    22   /* high-speed is a 2pt cross-correlation process for linking moving objects (high PM) */
    23   if (FIT_TARGET == TARGET_HIGH_SPEED) {
    24     high_speed_catalogs ();
    25     exit (0);
    26   }
     12  switch (FIT_TARGET) {
     13    case TARGET_UPDATE_OBJECTS:
     14      /* the object analysis is a separate process iterating over catalogs */
     15      relastro_objects (skylist, 0, NULL);
     16      exit (0);
    2717
    28   /* a special method to manually merge unlinked detections of sources togther (not parallel) */
    29   if (FIT_TARGET == TARGET_MERGE_SOURCE) {
    30     relastro_merge_source ();
    31     exit (0);
    32   }
     18    case TARGET_HIGH_SPEED:
     19      /* high-speed is a 2pt cross-correlation process for linking moving objects (high PM) */
     20      high_speed_catalogs (sky, skylist, 0, NULL);
     21      exit (0);
    3322
    34   /* register database handle with shutdown procedure */
    35   set_db (&db);
     23    case TARGET_MERGE_SOURCE:
     24      /* a special method to manually merge unlinked detections of sources togther (not parallel) */
     25      relastro_merge_source (sky);
     26      exit (0);
    3627
    37   /* lock and load the image db table */
    38   status = dvo_image_lock (&db, ImageCat, 60.0, (UPDATE ? LCK_XCLD : LCK_SOFT));
    39   if (!status) Shutdown ("ERROR: failure to lock image catalog %s", db.filename);
    40   if (db.dbstate == LCK_EMPTY) Shutdown ("ERROR: No images in catalog %s (1)", db.filename);
    41   if (!dvo_image_load (&db, VERBOSE, FALSE)) Shutdown ("can't read image catalog %s", db.filename);
    42   MARKTIME("load image data: %f sec\n", dtime);
    43 
    44   /* load regions and images based on specified sky patch (default depth) */
    45   skylist = load_images (&db, &UserPatch);
    46   MARKTIME("load images: %f sec\n", dtime);
    47 
    48   initCoords();
    49 
    50   /* load catalog data from region files : subselect high-quality measurements */
    51   // XXX pass in the image table
    52   // XXX who carries the image grid?
    53   catalog = load_catalogs (skylist, &Ncatalog, TRUE);
    54   MARKTIME("load catalog data: %f sec\n", dtime);
    55 
    56   /* match measurements with images */
    57   initImageBins (catalog, Ncatalog, TRUE);
    58   MARKTIME("make image bins: %f sec\n", dtime);
    59 
    60   findImages (catalog, Ncatalog, TRUE);
    61   MARKTIME("set up image indexes: %f sec\n", dtime);
    62 
    63   if (PLOTSTUFF) {
    64     // plot_star_coords (catalog, Ncatalog);
    65     // plot_mosaic_fields (catalog);
    66   }
    67 
    68   // set test points based on the starmap
    69   createStarMapPoints();
    70 
    71   /* major modes */
    72   switch (FIT_TARGET) {
    7328    case TARGET_SIMPLE:
    74       for (i = 0; i < NLOOP; i++) {
    75         UpdateObjects (catalog, Ncatalog);
    76         UpdateSimple (catalog, Ncatalog);
    77       }
    78       break;
    79 
    8029    case TARGET_CHIPS:
    81       for (i = 0; i < NLOOP; i++) {
    82         UpdateObjects (catalog, Ncatalog);
    83         UpdateChips (catalog, Ncatalog);
    84         MARKTIME("update chips: %f sec\n", dtime);
    85       }
    86       // create summary plots of the process
    87       // relastroVisualSummaryChips();
    88       break;
    89 
    9030    case TARGET_MOSAICS:
    91       for (i = 0; i < NLOOP; i++) {
    92         UpdateObjects (catalog, Ncatalog);
    93         UpdateMosaic (catalog, Ncatalog);
    94       }
    95       break;
     31      relastro_images (skylist);
     32      exit (0);
    9633
    9734    default:
    98       fprintf (stderr, "programming error at %s:%d", __FILE__, __LINE__);
    99       exit (2);
     35      fprintf (stderr, "impossible!\n");
     36      abort();
    10037  }
    101 
    102   if (!UPDATE) exit (0);
    103 
    104   // free the image / measurement pointers
    105   freeImageBins (1);
    106 
    107   reload_images (&db);
    108    
    109   // iterate over catalogs to make detection coordinates consistant
    110   UpdateObjectOffsets (skylist);
    111 
    112   // iterate over catalogs to make detection coordinates consistant
    113   // FixProblemImages (skylist);
    114 
    115   // save the updated image parameters
    116   dvo_image_update (&db, VERBOSE);
    117   dvo_image_unlock (&db);
    118 
    119   exit (0);
     38  exit (1);
    12039}

  • branches/eam_branches/ipp-20111122/Ohana/src/relastro/src/relastro_client.c

    r33396 r33448  
    44// owned by that host. 
    55
    6 // there are two modes:
     6// there are four modes:
    77
    8 // relastro_client -load : extract the bright catalog subset from the client's tables and
    9 //                        save to a local FITS table for relastro to load and used
     8// relastro_client -load-objects : extract the bright catalog subset from the client's
     9//                                 tables and save to a local FITS table for relastro to
     10//                                 load and used
    1011
    11 // (relastro loads the bright catalog subsets from all clients, then uses this to determine
    12 // the per-image zero points (and potentially the flat-field corrections)
     12// relastro_client -update-offsets : load image table containing the updated astrometry
     13//                                   information and determine per-object astrometry
    1314
    14 // relastro_client -update : load image table containing the update zero point information
    15 //                          and apply this to the client's tables
     15// relastro_client -update-objects : load image table containing the updated astrometry
     16//                                   information and determine per-object astrometry
     17
     18// relastro loads the bright catalog subsets from all clients, then determines the
     19// astrometry for chips/mosaics, while iteratively improving the per-object mean positions
    1620
    1721int main (int argc, char **argv) {
    … …  
    2024  initialize_client (argc, argv);
    2125
    22   // load the current sky table (layout of all SkyRegions)
    23   SkyTable *sky = SkyTableLoadOptimal (CATDIR, NULL, NULL, TRUE, -1, VERBOSE);
    24   if (!sky) {
    25       fprintf (stderr, "ERROR running loading sky table from %s\n", CATDIR);
    26       exit (2);
    27   }
     26  SkyTable *sky = SkyTableLoadOptimal (CATDIR, SKY_TABLE, GSCFILE, TRUE, SKY_DEPTH, VERBOSE);
    2827  SkyTableSetFilenames (sky, CATDIR, "cpt");
     28  SkyList *skylist = SkyListByPatch (sky, -1, &UserPatch);
    2929
    30   SkyList *skylist = SkyListByPatch (sky, -1, &UserPatch);
    31   if (!skylist) {
    32       fprintf (stderr, "ERROR setting up skylist for %s\n", CATDIR);
    33       exit (2);
    34   }
    35  
    36   switch (MODE) {
    37     case MODE_LOAD: {
     30  switch (FIT_TARGET) {
     31
     32    case TARGET_LOAD_OBJECTS: {
     33      // USAGE: relastro_client -load-objects
    3834      int Ncatalog;
    39       Catalog *catalog = load_catalogs (skylist, &Ncatalog, HOST_ID, HOSTDIR);
     35      Catalog *catalog = load_catalogs (skylist, &Ncatalog, TRUE, HOST_ID, HOSTDIR);
    4036      if (!catalog) {
    4137          fprintf (stderr, "ERROR loading catalogs from %s\n", CATDIR);
    … …  
    5046    }
    5147     
    52     case MODE_UPDATE: {
    53       // load the image subset table from the specified location
    54       off_t Nimage;
    55       ImageSubset *image = ImageSubsetLoad (IMAGES, &Nimage);
    56       if (!image) {
    57           fprintf (stderr, "ERROR loading image subset %s\n", CATDIR);
    58           exit (2);
    59       }
    60      
    61       // save the available image information in the static array in ImageOps.c
    62       initImagesSubset (image, NULL, Nimage);
    63 
    64       // load the flat-field correction table from CATDIR
    65       char flatcorrFile[1024];
    66       snprintf (flatcorrFile, 1024, "%s/flatcorr.fits", CATDIR);
    67       FlatCorrectionTable *flatcorr = FlatCorrectionLoad (flatcorrFile, VERBOSE);
    68 
    69       reload_catalogs (skylist, flatcorr, HOST_ID, HOSTDIR);
     48    case TARGET_UPDATE_OBJECTS: {
     49      // USAGE: relastro_client -update-objects
     50      relastro_objects (skylist, HOST_ID, HOSTDIR);
    7051      break;
    7152    }
    7253
    73     case MODE_UPDATE_OBJECTS: {
    74       relastro_objects (HOST_ID, HOSTDIR);
     54    case TARGET_HIGH_SPEED: {
     55      // USAGE: relastro_client -high-speed
     56      high_speed_catalogs (sky, skylist, HOST_ID, HOSTDIR);
     57      break;
     58    }
     59
     60      // XXX loading the images is fairly costly -- see if we can do an image subset
     61    case TARGET_UPDATE_OFFSETS: {
     62      FITS_DB db;
     63     
     64      set_db (&db);
     65
     66      /* lock and load the image db table */
     67      int status = dvo_image_lock (&db, ImageCat, 60.0, LCK_SOFT);
     68      if (!status) Shutdown ("ERROR: failure to lock image catalog %s", db.filename);
     69      if (db.dbstate == LCK_EMPTY) Shutdown ("ERROR: No images in catalog %s (1)", db.filename);
     70      if (!dvo_image_load (&db, VERBOSE, FALSE)) Shutdown ("can't read image catalog %s", db.filename);
     71
     72      /* load regions and images based on specified sky patch (default depth) */
     73      load_images (&db, skylist);
     74     
     75      // // load the image subset table from the specified location
     76      // off_t Nimage;
     77      // ImageSubset *image = ImageSubsetLoad (IMAGES, &Nimage);
     78      // if (!image) {
     79      //          fprintf (stderr, "ERROR loading image subset %s\n", CATDIR);
     80      //          exit (2);
     81      // }
     82      //
     83      // // save the available image information in the static array in ImageOps.c
     84      // initImagesSubset (image, NULL, Nimage);
     85      //
     86      // // load the flat-field correction table from CATDIR
     87      // char flatcorrFile[1024];
     88      // snprintf (flatcorrFile, 1024, "%s/flatcorr.fits", CATDIR);
     89      // FlatCorrectionTable *flatcorr = FlatCorrectionLoad (flatcorrFile, VERBOSE);
     90
     91      UpdateObjectOffsets (skylist, HOST_ID, HOSTDIR);
     92     
    7593      break;
    7694    }

  • branches/eam_branches/ipp-20111122/Ohana/src/relastro/src/relastro_merge_source.c

    r32346 r33448  
    11# include "relastro.h"
    22
    3 int relastro_merge_source () {
     3int relastro_merge_source (SkyTable *sky) {
    44
    55  int i, m;
    6   SkyTable *sky = NULL;
    76  Catalog catalog_src, catalog_dst;
    87
    … …  
    1312  dvo_set_catdir(CATDIR);
    1413
    15   // load the current sky table (layout of all SkyRegions)
    16   sky = SkyTableLoadOptimal (CATDIR, NULL, NULL, FALSE, SKY_DEPTH, VERBOSE);
    17   SkyTableSetFilenames (sky, CATDIR, "cpt");
    18  
    1914  // we are merging the detections of the src into the dst : OBJ_ID_SRC,CAT_ID_SRC -> OBJ_ID_DST,CAT_ID_DST
    2015

  • branches/eam_branches/ipp-20111122/Ohana/src/relastro/src/relastro_objects.c

    r33396 r33448  
    33int relastro_objects_parallel (SkyList *sky);
    44
    5 int relastro_objects (int hostID, char *hostpath) {
     5int relastro_objects (SkyList *skylist, int hostID, char *hostpath) {
    66
    77  int i, j, k, m;
    88
    9   SkyTable *sky = NULL;
    10   SkyList *skylist = NULL;
    119  Catalog catalog;
    12 
    13   // load the current sky table (layout of all SkyRegions)
    14   sky = SkyTableLoadOptimal (CATDIR, SKY_TABLE, GSCFILE, TRUE, SKY_DEPTH, VERBOSE);
    15   SkyTableSetFilenames (sky, CATDIR, "cpt");
    16  
    17   // determine the populated SkyRegions overlapping the requested area (default depth)
    18   if (UserCatalog) {
    19     skylist = SkyRegionByPoint (sky, -1, UserCatalogRA, UserCatalogDEC);
    20   } else {
    21     skylist = SkyListByPatch (sky, -1, &UserPatch);
    22   }
    2310
    2411  // XXX need to decide how to determine PARALLEL mode...
    2512  if (PARALLEL & !hostID) {
    26     relphot_objects_parallel (skylist);
     13    relastro_objects_parallel (skylist);
    2714    return TRUE;
    2815  }
    … …  
    4229    catalog.Nsecfilt  = GetPhotcodeNsecfilt ();
    4330
     31    // loads Average, Measure, SecFilt
    4432    if (!dvo_catalog_open (&catalog, skylist[0].regions[i], VERBOSE, "w")) {
    4533      fprintf (stderr, "ERROR: failure reading catalog %s\n", catalog.filename);
    … …  
    5341    }
    5442
     43    // operates on full tables (Average, Measure, Secfilt)
    5544    if (FlagOutlier) {
    5645      FlagOutliers(&catalog);
    … …  
    6857    }
    6958
     59    populate_tiny_values(&catalog, DVO_TV_MEASURE);
     60
    7061    UpdateObjects (&catalog, 1);
     62
     63    free_tiny_values(&catalog);
    7164
    7265    if (!UPDATE) {
    … …  
    10396    table->hosts[i].pathname = tmppath;
    10497
    105     char catalogFile[512];
    106     snprintf (catalogFile, 512, "%s/relastro.catalog.subset.dat", table->hosts[i].pathname);
    107 
    10898    // options / arguments that can affect relastro_client -update-objects:
    10999    // VERBOSE, VERBOSE2
    110100    // RESET (-reset)
    111101    // TimeSelect -time
    112     // DophotSelect
    113102    // (note that psfQual is applied rigidly at 0.85, as is the galaxy test)
    114     // MAG_LIM
    115     // SIGMA_LIM
    116103    // ImagSelect, ImagMin, ImagMax
    117104    // MaxDensityUse, MaxDensityValue
    … …  
    127114
    128115    char tmpline[1024];
     116    if (FIT_MODE == FIT_PM_ONLY)         { snprintf (tmpline, 1024, "%s -pm",    command);           strcpy (command, tmpline); }
     117    if (FIT_MODE == FIT_PAR_ONLY)        { snprintf (tmpline, 1024, "%s -par",   command);           strcpy (command, tmpline); }
     118    if (FIT_MODE == FIT_PM_AND_PAR)      { snprintf (tmpline, 1024, "%s -pmpar", command);           strcpy (command, tmpline); }
     119
    129120    if (VERBOSE)       { snprintf (tmpline, 1024, "%s -v",              command);                    strcpy (command, tmpline); }
    130121    if (VERBOSE2)      { snprintf (tmpline, 1024, "%s -vv",             command);                    strcpy (command, tmpline); }
    131122    if (RESET)         { snprintf (tmpline, 1024, "%s -reset",          command);                    strcpy (command, tmpline); }
    132     if (KEEP_UBERCAL)  { snprintf (tmpline, 1024, "%s -keep-ubercal",   command);                    strcpy (command, tmpline); }
     123
     124    if (ImagSelect)    { snprintf (tmpline, 1024, "%s -instmag %f %f",  command, ImagMin, ImagMax);  strcpy (command, tmpline); }
     125    if (MaxDensityUse) { snprintf (tmpline, 1024, "%s -max-density %f", command, MaxDensityValue);   strcpy (command, tmpline); }
     126   
     127    if (FlagOutlier)   { snprintf (tmpline, 1024, "%s -clip %d",        command, CLIP_THRESH);       strcpy (command, tmpline); }
     128   
     129    if (PHOTCODE_KEEP_LIST) { snprintf (tmpline, 1024, "%s +photcode %s", command, PHOTCODE_KEEP_LIST); strcpy (command, tmpline); }
     130    if (PHOTCODE_SKIP_LIST) { snprintf (tmpline, 1024, "%s -photcode %s", command, PHOTCODE_SKIP_LIST); strcpy (command, tmpline); }
     131    if (PhotFlagSelect)    { snprintf (tmpline, 1024, "%s +photflags",   command);                     strcpy (command, tmpline); }
     132    if (PhotFlagBad)       { snprintf (tmpline, 1024, "%s +photflagbad %d", command, PhotFlagBad);     strcpy (command, tmpline); }
     133    if (PhotFlagPoor)      { snprintf (tmpline, 1024, "%s +photflagpoor %d", command, PhotFlagPoor);   strcpy (command, tmpline); }
     134    // XXX note that the above pass in the flag as decimal -- also note that args.c cannot handle 0xHEX values
     135
     136    if (TimeSelect) {
     137      char *tstart = ohana_sec_to_date (TSTART);
     138      char *tstop  = ohana_sec_to_date (TSTOP);
     139      snprintf (tmpline, 1024, "%s -time %s %s", command, tstart, tstop);
     140      free (tstart);
     141      free (tstop);
     142      strcpy (command, tmpline);
     143    }
    133144
    134145    fprintf (stderr, "command: %s\n", command);
Note: See TracChangeset for help on using the changeset viewer.