Index: /trunk/Ohana/src/libohana/src/dvo_convert_elixir.c =================================================================== --- /trunk/Ohana/src/libohana/src/dvo_convert_elixir.c (revision 4830) +++ /trunk/Ohana/src/libohana/src/dvo_convert_elixir.c (revision 4830) @@ -0,0 +1,108 @@ +# include +# include + +/* convert loneos-style measures to internal measures */ +Measure *MeasureElixirToInternal (MeasureElixir *in, int Nvalues) { + + int i; + Measure *out; + + ALLOCATE (out, Measure, Nvalues); + + for (i = 0; i < Nvalues; i++) { + out[i].dR_PS = in[i].dR * 0.01; + out[i].dD_PS = in[i].dD * 0.01; + out[i].M = in[i].M; + out[i].Mcal = in[i].Mcal; + out[i].Mgal = in[i].Mgal; + out[i].airmass = in[i].airmass; + out[i].FWx = in[i].FWx; + out[i].dM = in[i].dM; + out[i].fwy = in[i].fwy; + out[i].theta = in[i].theta; + out[i].dophot = in[i].dophot; + out[i].source = in[i].source; + out[i].t = in[i].t; + out[i].averef = in[i].averef; + out[i].dt = in[i].dt; + out[i].flags = in[i].flags; + } + return (out); +} + +/* convert internal measures to loneos-style measures */ +MeasureElixir *MeasureInternalToElixir (Measure *in, int Nvalues) { + + int i; + MeasureElixir *out; + + ALLOCATE (out, MeasureElixir, Nvalues); + + for (i = 0; i < Nvalues; i++) { + out[i].dR = in[i].dR_PS * 100.0; + out[i].dD = in[i].dD_PS * 100.0; + out[i].M = in[i].M; + out[i].Mcal = in[i].Mcal; + out[i].Mgal = in[i].Mgal; + out[i].airmass = in[i].airmass; + out[i].FWx = in[i].FWx; + out[i].dM = in[i].dM; + out[i].fwy = in[i].fwy; + out[i].theta = in[i].theta; + out[i].dophot = in[i].dophot; + out[i].source = in[i].source; + out[i].t = in[i].t; + out[i].averef = in[i].averef; + out[i].dt = in[i].dt; + out[i].flags = in[i].flags; + } + return (out); +} + +Average *AverageElixirToInternal (AverageElixir *in, int Nvalues) { + + int i; + Average *out; + + ALLOCATE (out, Average, Nvalues); + + for (i = 0; i < Nvalues; i++) { + out[i].R_PS = in[i].R; + out[i].D_PS = in[i].D; + out[i].M = in[i].M; + out[i].Xp = in[i].Xp; + out[i].Xm = in[i].Xm; + out[i].Nm = in[i].Nm; + out[i].Nn = in[i].Nn; + out[i].code = in[i].code; + out[i].offset = in[i].offset; + out[i].missing = in[i].missing; + out[i].dM = in[i].dM; + out[i].Xg = in[i].Xg; + } + return (out); +} + +AverageElixir *AverageInternalToElixir (Average *in, int Nvalues) { + + int i; + AverageElixir *out; + + ALLOCATE (out, AverageElixir, Nvalues); + + for (i = 0; i < Nvalues; i++) { + out[i].R = in[i].R_PS; + out[i].D = in[i].D_PS; + out[i].M = in[i].M; + out[i].Xp = in[i].Xp; + out[i].Xm = in[i].Xm; + out[i].Nm = in[i].Nm; + out[i].Nn = in[i].Nn; + out[i].code = in[i].code; + out[i].offset = in[i].offset; + out[i].missing = in[i].missing; + out[i].dM = in[i].dM; + out[i].Xg = in[i].Xg; + } + return (out); +} Index: /trunk/Ohana/src/libohana/src/dvo_convert_loneos.c =================================================================== --- /trunk/Ohana/src/libohana/src/dvo_convert_loneos.c (revision 4830) +++ /trunk/Ohana/src/libohana/src/dvo_convert_loneos.c (revision 4830) @@ -0,0 +1,143 @@ +# include +# include + +/* convert loneos-style measures to internal measures */ +Measure *MeasureLoneosToInternal (MeasureLoneos *in, int Nvalues) { + + int i; + Measure *out; + + ALLOCATE (out, Measure, Nvalues); + + for (i = 0; i < Nvalues; i++) { + out[i].dR_PS = in[i].dR * 0.01; + out[i].dD_PS = in[i].dD * 0.01; + out[i].M = in[i].M; + out[i].Mcal = in[i].Mcal; + out[i].dM = in[i].dM; + out[i].dophot = in[i].dophot; + out[i].source = in[i].source; + out[i].t = in[i].t; + + /* flags and averef are now split */ + out[i].averef = in[i].averef & 0x00ffffff; + out[i].flags = (in[i].averef & 0xff000000) >> 24; + + /* these values don't exist in the Loneos format */ + out[i].Mgal = in[i].M; + out[i].dt = 0; + out[i].airmass = 0; + out[i].FWx = 0; + out[i].fwy = 0; + out[i].theta = 0; + } + return (out); +} + +/* convert internal measures to loneos-style measures */ +MeasureLoneos *MeasureInternalToLoneos (Measure *in, int Nvalues) { + + int i; + MeasureLoneos *out; + + ALLOCATE (out, MeasureLoneos, Nvalues); + + for (i = 0; i < Nvalues; i++) { + out[i].dR = in[i].dR_PS * 100.0; + out[i].dD = in[i].dD_PS * 100.0; + out[i].M = in[i].M; + out[i].dM = in[i].dM; + out[i].Mcal = in[i].Mcal; + out[i].dophot = in[i].dophot; + out[i].source = in[i].source; + out[i].t = in[i].t; + + /* flags and averef are merged in Loneos */ + out[i].averef = (in[i].averef & 0x00ffffff) | (in[i].flags << 24); + } + return (out); +} + +Average *AverageLoneosToInternal (AverageLoneos *in, int Nvalues) { + + int i; + Average *out; + + ALLOCATE (out, Average, Nvalues); + + for (i = 0; i < Nvalues; i++) { + out[i].R_PS = in[i].R; + out[i].D_PS = in[i].D; + out[i].M = in[i].M; + out[i].Xp = in[i].Xp; + out[i].Xm = in[i].Xm; + out[i].Nm = in[i].Nm; + out[i].Nn = in[i].Nn; + out[i].code = in[i].code; + out[i].offset = in[i].offset; + out[i].missing = in[i].missing; + + /* these don't exist in Loneos */ + out[i].dM = 0xffff; + out[i].Xg = 0xffff; + } + return (out); +} + +AverageLoneos *AverageInternalToLoneos (Average *in, int Nvalues) { + + int i; + AverageLoneos *out; + + ALLOCATE (out, AverageLoneos, Nvalues); + + for (i = 0; i < Nvalues; i++) { + out[i].R = in[i].R_PS; + out[i].D = in[i].D_PS; + out[i].M = in[i].M; + out[i].Xp = in[i].Xp; + out[i].Xm = in[i].Xm; + out[i].Nm = in[i].Nm; + out[i].Nn = in[i].Nn; + out[i].code = in[i].code; + out[i].offset = in[i].offset; + out[i].missing = in[i].missing; + } + return (out); +} + +/* +I need to add these to the dvo_convert structures: + +SecFilt *SecFiltLoneosToInternal (SecFiltLoneos *in, int Nvalues); +SecFiltLoneos *SecFiltInternalToLoneos (SecFilt *in, int Nvalues); + +SecFilt *SecFiltLoneosToInternal (SecFiltLoneos *in, int Nvalues) { + + int i; + SecFilt *out; + + ALLOCATE (out, SecFilt, Nvalues); + + for (i = 0; i < Nvalues; i++) { + out[i].M = in[i].M; + out[i].Xm = in[i].Xm; + out[i].dM = 0xffff; + } + return (out); +} + +SecFiltLoneos *SecFiltInternalToLoneos (SecFilt *in, int Nvalues) { + + int i; + SecFiltLoneos *out; + + ALLOCATE (out, SecFiltLoneos, Nvalues); + + for (i = 0; i < Nvalues; i++) { + out[i].M = in[i].M; + out[i].Xm = in[i].Xm; + } + return (out); +} +*/ Index: /trunk/Ohana/src/libohana/src/dvo_convert_panstarrs.c =================================================================== --- /trunk/Ohana/src/libohana/src/dvo_convert_panstarrs.c (revision 4830) +++ /trunk/Ohana/src/libohana/src/dvo_convert_panstarrs.c (revision 4830) @@ -0,0 +1,108 @@ +# include +# include + +/* convert loneos-style measures to internal measures */ +Measure *MeasurePanstarrsToInternal (MeasurePanstarrs *in, int Nvalues) { + + int i; + Measure *out; + + ALLOCATE (out, Measure, Nvalues); + + for (i = 0; i < Nvalues; i++) { + out[i].dR_PS = in[i].dR; + out[i].dD_PS = in[i].dD; + out[i].M = in[i].M; + out[i].Mcal = in[i].Mcal; + out[i].Mgal = in[i].Mgal; + out[i].airmass = in[i].airmass; + out[i].FWx = in[i].FWx; + out[i].dM = in[i].dM; + out[i].fwy = in[i].fwy; + out[i].theta = in[i].theta; + out[i].dophot = in[i].dophot; + out[i].source = in[i].source; + out[i].t = in[i].t; + out[i].averef = in[i].averef; + out[i].dt = in[i].dt; + out[i].flags = in[i].flags; + } + return (out); +} + +/* convert internal measures to loneos-style measures */ +MeasurePanstarrs *MeasureInternalToPanstarrs (Measure *in, int Nvalues) { + + int i; + MeasurePanstarrs *out; + + ALLOCATE (out, MeasurePanstarrs, Nvalues); + + for (i = 0; i < Nvalues; i++) { + out[i].dR = in[i].dR_PS; + out[i].dD = in[i].dD_PS; + out[i].M = in[i].M; + out[i].Mcal = in[i].Mcal; + out[i].Mgal = in[i].Mgal; + out[i].airmass = in[i].airmass; + out[i].FWx = in[i].FWx; + out[i].dM = in[i].dM; + out[i].fwy = in[i].fwy; + out[i].theta = in[i].theta; + out[i].dophot = in[i].dophot; + out[i].source = in[i].source; + out[i].t = in[i].t; + out[i].averef = in[i].averef; + out[i].dt = in[i].dt; + out[i].flags = in[i].flags; + } + return (out); +} + +Average *AveragePanstarrsToInternal (AveragePanstarrs *in, int Nvalues) { + + int i; + Average *out; + + ALLOCATE (out, Average, Nvalues); + + for (i = 0; i < Nvalues; i++) { + out[i].R_PS = in[i].R; + out[i].D_PS = in[i].D; + out[i].M = in[i].M; + out[i].Xp = in[i].Xp; + out[i].Xm = in[i].Xm; + out[i].Nm = in[i].Nm; + out[i].Nn = in[i].Nn; + out[i].code = in[i].code; + out[i].offset = in[i].offset; + out[i].missing = in[i].missing; + out[i].dM = in[i].dM; + out[i].Xg = in[i].Xg; + } + return (out); +} + +AveragePanstarrs *AverageInternalToPanstarrs (Average *in, int Nvalues) { + + int i; + AveragePanstarrs *out; + + ALLOCATE (out, AveragePanstarrs, Nvalues); + + for (i = 0; i < Nvalues; i++) { + out[i].R = in[i].R_PS; + out[i].D = in[i].D_PS; + out[i].M = in[i].M; + out[i].Xp = in[i].Xp; + out[i].Xm = in[i].Xm; + out[i].Nm = in[i].Nm; + out[i].Nn = in[i].Nn; + out[i].code = in[i].code; + out[i].offset = in[i].offset; + out[i].missing = in[i].missing; + out[i].dM = in[i].dM; + out[i].Xg = in[i].Xg; + } + return (out); +} Index: /trunk/Ohana/src/libohana/src/phot_catalog_mef.c =================================================================== --- /trunk/Ohana/src/libohana/src/phot_catalog_mef.c (revision 4830) +++ /trunk/Ohana/src/libohana/src/phot_catalog_mef.c (revision 4830) @@ -0,0 +1,233 @@ +# include +# include + +int load_catalog_mef (Catalog *catalog, char mode, int VERBOSE) { + + int Nitems, Nbytes, Naverage, Nmeasure, Nmissing, Nsecfilt; + FILE *f; + + Header header; + Matrix matrix; + FTable ftable; + + f = catalog[0].f; + fseek (f, 0, SEEK_SET); + ftable.header = &header; + + /* read PHU */ + if (!fits_load_header (f, &catalog[0].header)) { + if (VERBOSE) fprintf (stderr, "catalog file does not exist: %s\n", catalog[0].filename); + return (FALSE); + } + /* matrix should be empty */ + if (!fits_fread_matrix (f, &matrix, &catalog[0].header)) { + if (VERBOSE) fprintf (stderr, "can't read primary matrix"); + return (FALSE); + } + /* get the components from the header */ + if (!fits_scan (&catalog[0].header, "NSTARS", "%d", 1, &Naverage)) return (FALSE); + if (!fits_scan (&catalog[0].header, "NMEAS", "%d", 1, &Nmeasure)) return (FALSE); + if (!fits_scan (&catalog[0].header, "NMISS", "%d", 1, &Nmissing)) return (FALSE); + if (!fits_scan (&catalog[0].header, "NSECFILT", "%d", 1, &Nsecfilt)) Nsecfilt = 0; + + /* validate table mode */ + + /* read Average table header */ + if (!fits_fread_header (f, &header)) { + if (VERBOSE) fprintf (stderr, "can't read table average header"); + return (FALSE); + } + /* read Average table data (or skip) */ + if (mode & LOAD_AVES) { + if (!fits_fread_ftable_data (f, &ftable)) { + if (VERBOSE) fprintf (stderr, "can't read table average data"); + return (FALSE); + } + catalog[0].average = FtableToAverage (&ftable, &catalog[0].Naverage, catalog[0].average_format); + if (Naverage != catalog[0].Naverage) { + fprintf (stderr, "Warning: mismatch between Naverage in PHU and Table headers (%d vs %d)\n", Naverage, catalog[0].Naverage); + } + } else { + Nbytes = fits_matrix_size (&header); + fseek (f, Nbytes, SEEK_CUR); + } + fits_free_header (&header); + /** free the ftable or not? data is being used still..? **/ + + /* read Measure table header */ + if (!fits_fread_header (f, &header)) { + if (VERBOSE) fprintf (stderr, "can't read table measure header"); + return (FALSE); + } + /* read Measure table data */ + if (mode & LOAD_MEAS) { + if (!fits_fread_ftable_data (f, &ftable)) { + if (VERBOSE) fprintf (stderr, "can't read table measure data"); + return (FALSE); + } + catalog[0].measure = FtableToMeasure (&ftable, &catalog[0].Nmeasure, catalog[0].measure_format); + if (Nmeasure != catalog[0].Nmeasure) { + fprintf (stderr, "Warning: mismatch between Nmeasure in PHU and Table headers (%d vs %d)\n", Nmeasure, catalog[0].Nmeasure); + } + } else { + Nbytes = fits_matrix_size (&header); + fseek (f, Nbytes, SEEK_CUR); + } + + /* read Missing table header */ + if (!fits_fread_header (f, &header)) { + if (VERBOSE) fprintf (stderr, "can't read table missing header"); + return (FALSE); + } + /* read Missing table data */ + if (mode & LOAD_MISS) { + if (!fits_fread_ftable_data (f, &ftable)) { + if (VERBOSE) fprintf (stderr, "can't read table missing data"); + return (FALSE); + } + /* no conversions currently defined */ + catalog[0].missing = fits_table_get_Missing (&ftable, &catalog[0].Nmissing); + if (Nmissing != catalog[0].Nmissing) { + fprintf (stderr, "Warning: mismatch between Nmissing in PHU and Table headers (%d vs %d)\n", Nmissing, catalog[0].Nmissing); + } + } else { + Nbytes = fits_matrix_size (&header); + fseek (f, Nbytes, SEEK_CUR); + } + + /* read secfilt table header */ + if (!fits_fread_header (f, &header)) { + if (VERBOSE) fprintf (stderr, "can't read table secfilt header"); + return (FALSE); + } + /* read secfilt table data */ + if (mode & LOAD_SECF) { + if (!fits_fread_ftable_data (f, &ftable)) { + if (VERBOSE) fprintf (stderr, "can't read table secfilt data"); + return (FALSE); + } + /* no conversions currently defined */ + catalog[0].secfilt = fits_table_get_SecFilt (&ftable, &Nitems); + catalog[0].Nsecfilt = Nitems / catalog[0].Naverage; + if (Nsecfilt != catalog[0].Nsecfilt) { + fprintf (stderr, "Warning: mismatch between Nsecfilt in PHU and Table headers (%d vs %d)\n", Nsecfilt, catalog[0].Nsecfilt); + } + } else { + /* no real need to skip the data array here... */ + Nbytes = fits_matrix_size (&header); + fseek (f, Nbytes, SEEK_CUR); + } + return (TRUE); +} + +/* save_catalog_mef writes a complete new file from scratch */ + +int save_catalog_mef (Catalog *catalog, char VERBOSE) { + + FILE *f; + Matrix matrix; + Header header; + FTable ftable; + + if (catalog[0].Naverage == 0) { + if (VERBOSE) fprintf (stderr, "no stars in catalog, skipping\n"); + return (TRUE); + } + + /* make sure header is consistent with data */ + fits_modify (&catalog[0].header, "NSTARS", "%d", 1, catalog[0].Naverage); + fits_modify (&catalog[0].header, "NMEAS", "%d", 1, catalog[0].Nmeasure); + fits_modify (&catalog[0].header, "NMISS", "%d", 1, catalog[0].Nmissing); + fits_modify (&catalog[0].header, "NSECFILT", "%d", 1, catalog[0].Nsecfilt); + fits_modify (&catalog[0].header, "EXTEND", "%t", 1, TRUE); + + f = catalog[0].f; + /* rewind file pointers and truncate */ + fseek (f, 0, SEEK_SET); + ftruncate (fileno (catalog[0].f), 0); + ftable.header = &header; + + /* write table PHU header */ + if (!fits_fwrite_header (catalog[0].f, &catalog[0].header)) { + fprintf (stderr, "can't write primary header"); + return (FALSE); + } + + /* this is probably a NOP, do I have to keep it in? */ + fits_create_matrix (&catalog[0].header, &matrix); + if (!fits_fwrite_matrix (catalog[0].f, &matrix)) { + fprintf (stderr, "can't write primary matrix"); + return (FALSE); + } + fits_free_matrix (&matrix); + + /* write out Average table (convert to FITS table format) */ + AverageToFtable (&ftable, catalog[0].average, catalog[0].Naverage, catalog[0].average_format); + if (!fits_fwrite_Theader (catalog[0].f, &header)) { + fprintf (stderr, "can't write table header"); + return (FALSE); + } + if (!fits_fwrite_table (catalog[0].f, &ftable)) { + fprintf (stderr, "can't write table data"); + return (FALSE); + } + fits_free_table (&ftable); + fits_free_header (&header); + + /* write out Measure table (convert to FITS table format) */ + MeasureToFtable (&ftable, catalog[0].measure, catalog[0].Nmeasure, catalog[0].measure_format); + if (!fits_fwrite_Theader (catalog[0].f, &header)) { + fprintf (stderr, "can't write table header"); + return (FALSE); + } + if (!fits_fwrite_table (catalog[0].f, &ftable)) { + fprintf (stderr, "can't write table data"); + return (FALSE); + } + fits_free_table (&ftable); + fits_free_header (&header); + + /* write out Missing table (convert to FITS table format) */ + fits_table_set_Missing (&ftable, catalog[0].missing, catalog[0].Nmissing); + if (!fits_fwrite_Theader (catalog[0].f, &header)) { + fprintf (stderr, "can't write table header"); + return (FALSE); + } + if (!fits_fwrite_table (catalog[0].f, &ftable)) { + fprintf (stderr, "can't write table data"); + return (FALSE); + } + fits_free_table (&ftable); + fits_free_header (&header); + + /* write out Secfilt table (convert to FITS table format) */ + /*** note that Nsecfilt is NOT the number of rows! ***/ + fits_table_set_SecFilt (&ftable, catalog[0].secfilt, catalog[0].Naverage * catalog[0].Nsecfilt); + if (!fits_fwrite_Theader (catalog[0].f, &header)) { + fprintf (stderr, "can't write table header"); + return (FALSE); + } + if (!fits_fwrite_table (catalog[0].f, &ftable)) { + fprintf (stderr, "can't write table data"); + return (FALSE); + } + fits_free_table (&ftable); + fits_free_header (&header); + + return (TRUE); +} + +/* + catalog data is: + header (bitpix == 8) + matrix (empty) + average header + average table + measure header + measure table + missing header + missing table + secfilt header + secfilt table +*/ + Index: /trunk/Ohana/src/libohana/src/phot_catalog_raw.c =================================================================== --- /trunk/Ohana/src/libohana/src/phot_catalog_raw.c (revision 4830) +++ /trunk/Ohana/src/libohana/src/phot_catalog_raw.c (revision 4830) @@ -0,0 +1,254 @@ +# include +# include + +/* read data from raw-style catalog file; set data format values based on file */ + +int load_catalog_raw (Catalog *catalog, char mode, int VERBOSE) { + + int Nitems, nitems; + int i, Nmeas, Nmiss, size; + int NewMeasure, Nskip; + int AverageSize, MeasureSize, MissingSize, SecFiltSize; + FILE *f; + struct stat filestatus; + char format[80], telescope[80]; + + f = catalog[0].f; + fseek (f, 0, SEEK_SET); + + /* read catalog header */ + if (!fits_load_header (f, &catalog[0].header)) { + if (VERBOSE) fprintf (stderr, "catalog file does not exist: %s\n", catalog[0].filename); + return (FALSE); + } + + /* this is done by the load anyway... */ + fseek (f, catalog[0].header.size, SEEK_SET); + + /* get the components from the header */ + catalog[0].Naverage = catalog[0].Nmeasure = catalog[0].Nmissing = catalog[0].Nsecfilt = 0; + if (!fits_scan (&catalog[0].header, "NSTARS", "%d", 1, &catalog[0].Naverage)) return (FALSE); + if (!fits_scan (&catalog[0].header, "NMEAS", "%d", 1, &catalog[0].Nmeasure)) return (FALSE); + if (!fits_scan (&catalog[0].header, "NMISS", "%d", 1, &catalog[0].Nmissing)) return (FALSE); + if (!fits_scan (&catalog[0].header, "NSECFILT", "%d", 1, &catalog[0].Nsecfilt)) catalog[0].Nsecfilt = 0; + + /* determine average and measure format */ + strcpy (catalog[0].average_format, "NONE"); + strcpy (catalog[0].measure_format, "NONE"); + if (fits_scan (&catalog[0].header, "FORMAT", "%s", 1, format)) { + if (!strcmp (format, "INTERNAL")) { + sprintf (catalog[0].average_format, "DVO_AVERAGE"); + sprintf (catalog[0].measure_format, "DVO_MEASURE"); + } else { + sprintf (catalog[0].average_format, "DVO_AVERAGE_%s", format); + sprintf (catalog[0].measure_format, "DVO_MEASURE_%s", format); + } + goto got_format; + } + if (fits_scan (&catalog[0].header, "NEWMEAS", "%t", 1, &NewMeasure)) { + strcpy (catalog[0].average_format, "DVO_AVERAGE_ELIXIR"); + strcpy (catalog[0].measure_format, "DVO_MEASURE_ELIXIR"); + goto got_format; + } + if (fits_scan (&catalog[0].header, "TELESCOP", "%s", 1, telescope)) { + if (!strncmp (telescope, "LONEOS", strlen("LONEOS"))) { + strcpy (catalog[0].average_format, "DVO_AVERAGE_ELIXIR"); + strcpy (catalog[0].measure_format, "DVO_MEASURE_ELIXIR"); + goto got_format; + } + } + +got_format: + /* determine datatype sizes */ + AverageSize = 0; + if (!strcmp (catalog[0].average_format, "DVO_AVERAGE")) AverageSize = sizeof(Average); + if (!strcmp (catalog[0].average_format, "DVO_AVERAGE_ELIXIR")) AverageSize = sizeof(AverageElixir); + if (!strcmp (catalog[0].average_format, "DVO_AVERAGE_LONEOS")) AverageSize = sizeof(AverageLoneos); + if (!strcmp (catalog[0].average_format, "DVO_AVERAGE_PANSTARRS")) AverageSize = sizeof(AveragePanstarrs); + if (AverageSize == 0) { + fprintf (stderr, "average format is invalid: %s\n", catalog[0].average_format); + return (FALSE); + } + MeasureSize = 0; + if (!strcmp (catalog[0].measure_format, "DVO_MEASURE")) MeasureSize = sizeof(Measure); + if (!strcmp (catalog[0].measure_format, "DVO_MEASURE_ELIXIR")) MeasureSize = sizeof(MeasureElixir); + if (!strcmp (catalog[0].measure_format, "DVO_MEASURE_LONEOS")) MeasureSize = sizeof(MeasureLoneos); + if (!strcmp (catalog[0].measure_format, "DVO_MEASURE_PANSTARRS")) MeasureSize = sizeof(MeasurePanstarrs); + if (MeasureSize == 0) { + fprintf (stderr, "measure format is invalid: %s\n", catalog[0].measure_format); + return (FALSE); + } + MissingSize = sizeof (Missing); + SecFiltSize = sizeof (SecFilt); + + /* predicted file size - for double checking data validity */ + size = catalog[0].header.size; + size += AverageSize * catalog[0].Naverage; + size += MeasureSize * catalog[0].Nmeasure; + size += MissingSize * catalog[0].Nmissing; + size += SecFiltSize * catalog[0].Nsecfilt * catalog[0].Naverage; + + /* check that file size makes sense */ + if (stat (catalog[0].filename, &filestatus) == -1) { + if (VERBOSE) fprintf (stderr, "failed to get status of catalog\n"); + return (FALSE); + } + if (size > filestatus.st_size) { + if (VERBOSE) { + fprintf (stderr, "star catalog has inconsistent size\n"); + fprintf (stderr, "average: %d = %d bytes\n", catalog[0].Naverage, catalog[0].Naverage*AverageSize); + fprintf (stderr, "measure: %d = %d bytes\n", catalog[0].Nmeasure, catalog[0].Nmeasure*MeasureSize); + fprintf (stderr, "missing: %d = %d bytes\n", catalog[0].Nmissing, catalog[0].Nmissing*MissingSize); + fprintf (stderr, "secfilt: %d = %d bytes\n", catalog[0].Nsecfilt, catalog[0].Nsecfilt*SecFiltSize*catalog[0].Naverage); + fprintf (stderr, "expect: %d, found: %d\n", size, (int)filestatus.st_size); + } + return (FALSE); + } + if (size < filestatus.st_size) { + if (VERBOSE) fprintf (stderr, "warning: file larger than expected\n"); + } + + if (catalog[0].Naverage == 0) { + if (VERBOSE) fprintf (stderr, "no stars yet in catalog %s\n", catalog[0].filename); + return (TRUE); + } + + /* read and convert the averages (use a macro to clean this up?) */ + if (mode & LOAD_AVES) { + catalog[0].average = ReadRawAverage (catalog[0].f, catalog[0].Naverage, catalog[0].average_format); + } else { + /* skip over averages */ + Nskip = catalog[0].Naverage * AverageSize; + fseek (f, Nskip, SEEK_CUR); + } + + /* read and convert the measures (use a macro to clean this up?) */ + if (mode & LOAD_MEAS) { + catalog[0].measure = ReadRawMeasure (catalog[0].f, catalog[0].Nmeasure, catalog[0].measure_format); + } else { + /* skip over measures */ + Nskip = catalog[0].Nmeasure * MeasureSize; + fseek (f, Nskip, SEEK_CUR); + } + + /* read and convert missing */ + if (mode & LOAD_MISS) { + ALLOCATE (catalog[0].missing, Missing, MAX (catalog[0].Nmissing, 1)); + Nitems = catalog[0].Nmissing; + nitems = fread (catalog[0].missing, MissingSize, Nitems, f); + if (nitems != Nitems) { + if (VERBOSE) fprintf (stderr, "failed to read missing from catalog file %s (%d vs %d)\n", catalog[0].filename, nitems, Nitems); + return (FALSE); + } + fits_convert_Missing (catalog[0].missing, MissingSize, Nitems); + } else { + /* skip over missings */ + Nskip = catalog[0].Nmissing * MissingSize; + fseek (f, Nskip, SEEK_CUR); + } + + /* read and convert secfilt */ + if (mode & LOAD_SECF) { + Nitems = catalog[0].Naverage * catalog[0].Nsecfilt; + ALLOCATE (catalog[0].secfilt, SecFilt, MAX (Nitems, 1)); + nitems = fread (catalog[0].secfilt, SecFiltSize, Nitems, f); + if (nitems != Nitems) { + if (VERBOSE) fprintf (stderr, "failed to read secondary filters from catalog file %s (%d vs %d)\n", catalog[0].filename, nitems, Nitems); + return (FALSE); + } + fits_convert_SecFilt (catalog[0].secfilt, SecFiltSize, Nitems); + } else { + /* skip over secfilts */ + Nskip = catalog[0].Nsecfilt * SecFiltSize * catalog[0].Naverage; + fseek (f, Nskip, SEEK_CUR); + } + + if (VERBOSE) fprintf (stderr, "read %d stars from catalog file %s (%d measurements, %d missing, %d secondary filters)\n", + catalog[0].Naverage, catalog[0].filename, catalog[0].Nmeasure, catalog[0].Nmissing, catalog[0].Nsecfilt); + + /* check data integrity */ + if (mode & LOAD_AVES) { + for (i = Nmeas = Nmiss = 0; i < catalog[0].Naverage; i++) { + Nmeas += catalog[0].average[i].Nm; + Nmiss += catalog[0].average[i].Nn; + } + if ((Nmeas != catalog[0].Nmeasure) || (Nmiss != catalog[0].Nmissing)) { + if (VERBOSE) { + fprintf (stderr, "****** data in catalog %s is corrupt, sums don't check\n", catalog[0].filename); + fprintf (stderr, "****** Nmeas: %d, %d\n", Nmeas, catalog[0].Nmeasure); + fprintf (stderr, "****** Nmiss: %d, %d\n", Nmiss, catalog[0].Nmissing); + } + return (FALSE); + } + } + return (TRUE); +} + +int save_catalog_raw (Catalog *catalog, char VERBOSE) { + + int Nitems, nitems; + char *p; + FILE *f; + + if (catalog[0].Naverage == 0) { + if (VERBOSE) fprintf (stderr, "no stars in catalog, skipping\n"); + return (TRUE); + } + + /* make sure header is consistent with data */ + fits_modify (&catalog[0].header, "NSTARS", "%d", 1, catalog[0].Naverage); + fits_modify (&catalog[0].header, "NMEAS", "%d", 1, catalog[0].Nmeasure); + fits_modify (&catalog[0].header, "NMISS", "%d", 1, catalog[0].Nmissing); + fits_modify (&catalog[0].header, "NSECFILT", "%d", 1, catalog[0].Nsecfilt); + + /* specify the appropriate data format */ + if (!strcmp (catalog[0].average_format, "DVO_AVERAGE")) { + fits_modify (&catalog[0].header, "FORMAT", "%s", 1, "INTERNAL"); + } else { + p = catalog[0].average_format + strlen ("DVO_AVERAGE_"); + fits_modify (&catalog[0].header, "FORMAT", "%s", 1, p); + } + + /* rewind file pointers and truncate file */ + f = catalog[0].f; + fseek (f, 0, SEEK_SET); + ftruncate (fileno (catalog[0].f), 0); + + /* write header data (use fits_write_header?) */ + nitems = fwrite (catalog[0].header.buffer, 1, catalog[0].header.size, f); + if (nitems != catalog[0].header.size) { + if (VERBOSE) fprintf (stderr, "failed to write header\n"); + return (FALSE); + } + + /* write averages and measures */ + WriteRawAverage (f, catalog[0].average, catalog[0].Naverage, catalog[0].average_format); + WriteRawMeasure (f, catalog[0].measure, catalog[0].Nmeasure, catalog[0].measure_format); + + /* write missing data */ + Nitems = catalog[0].Nmissing; + fits_convert_Missing (catalog[0].missing, sizeof(Missing), Nitems); + nitems = fwrite (catalog[0].missing, sizeof(Missing), Nitems, f); + if (nitems != Nitems) { + if (VERBOSE) fprintf (stderr, "failed to write catalog file missing %s\n", catalog[0].filename); + return (FALSE); + } + + /* write secfilt data */ + Nitems = catalog[0].Naverage * catalog[0].Nsecfilt; + fits_convert_SecFilt (catalog[0].secfilt, sizeof(SecFilt), Nitems); + nitems = fwrite (catalog[0].secfilt, sizeof(SecFilt), Nitems, f); + if (nitems != Nitems) { + if (VERBOSE) fprintf (stderr, "failed to write catalog file secfilt %s\n", catalog[0].filename); + return (FALSE); + } + + return (TRUE); + +} + +/*** + some warnings: + - need to be wary of exit condition errors + - average_format and measure_format must be the same format (not checked) +***/ Index: /trunk/Ohana/src/libohana/src/time.c =================================================================== --- /trunk/Ohana/src/libohana/src/time.c (revision 4830) +++ /trunk/Ohana/src/libohana/src/time.c (revision 4830) @@ -0,0 +1,389 @@ +# include "ohana.h" + +/***** convert [-]00:00:00 to 0.0000 ****/ +int dms_to_ddd (double *Value, char *string) { + + int valid, neg, status; + double tmp, value; + char *p1, *p2, *px; + + valid = FALSE; + neg = FALSE; + stripwhite (string); + p1 = string; + px = string + strlen(string); + + if (string[0] == '-') { + valid = TRUE; + neg = TRUE; + p1 = &string[1]; + } + if (string[0] == '+') { + valid = TRUE; + neg = FALSE; + p1 = &string[1]; + } + if (isdigit(string[0])) { + valid = TRUE; + p1 = &string[0]; + } + if (!valid) { return (FALSE); } + + status = 1; + tmp = strtod (p1, &p2); + if (p2 == p1) return (FALSE); /* entry not a number: +fred */ + value = tmp; + if (p2 == px) goto escape; /* entry only number: +1.0 */ + p1 = p2 + 1; + + tmp = strtod (p1, &p2); + if (p2 == p1) goto escape; /* entry not a number: +1:fred */ + status = 2; + value += tmp / 60.0; + if (p2 == px) goto escape; /* entry only number: +1:1 */ + p1 = p2 + 1; + + tmp = strtod (p1, &p2); + if (p2 == p1) goto escape; /* entry not a number: +1:1:fred */ + value += tmp / 3600.0; + + escape: + if (neg) { + value *= -1; + } + *Value = value; + + return (status); +} + +/**********/ +int str_to_radec (double *ra, double *dec, char *str1, char *str2) { + + double Ra, Dec; + + *ra = *dec = 0; + switch (dms_to_ddd (&Ra, str1)) { + case 0: + fprintf (stderr, "syntax error in RA\n"); + return (FALSE); + case 1: + break; + case 2: + Ra = Ra * 15; + break; + } + switch (dms_to_ddd (&Dec, str2)) { + case 0: + fprintf (stderr, "syntax error in DEC\n"); + return (FALSE); + case 1: + case 2: + break; + } + *ra = Ra; + *dec = Dec; + return (TRUE); +} + +/**********/ +int chk_time (char *line) { + + char *p1, *p2; + double tmp; + int mode; + + p1 = line; + tmp = strtod (p1, &p2); + mode = TIME_DATE; + if (p2 == p1 + strlen (p1) - 1) { + if (*p2 == 'd') { + mode = TIME_DAYS; + } + if (*p2 == 'h') { + mode = TIME_HOURS; + } + if (*p2 == 'm') { + mode = TIME_MINUTES; + } + if (*p2 == 's') { + mode = TIME_SECONDS; + } + if (*p2 == 'j') { + mode = TIME_JD; + } + if (*p2 == 'J') { + mode = TIME_MJD; + } + } + return (mode); +} + +/**********/ +int str_to_time (char *line, time_t *second) { + + char *tmpline; + struct tm *gmt; + struct timeval now; + double jd; + time_t tsec; + + if (!strcasecmp (line, "NOW")) { + gettimeofday (&now, (struct timezone *) NULL); + *second = now.tv_sec; + return (TRUE); + } + + if (!strncasecmp (line, "TODAY", 5)) { + gettimeofday (&now, (struct timezone *) NULL); + if (line[5]) { /* line has extra data (ie, hh:mm:ss) */ + tsec = now.tv_sec; + ALLOCATE (tmpline, char, 64); + gmt = gmtime (&tsec); + sprintf (tmpline, "%04d/%02d/%02d,%s", 1900 + gmt[0].tm_year, gmt[0].tm_mon+1, gmt[0].tm_mday, &line[6]); + *second = date_to_sec (tmpline); + free (tmpline); + return (TRUE); + } else { + *second = 86400 * ((int)(now.tv_sec / 86400)); + return (TRUE); + } + } + + switch (chk_time (line)) { + case 0: + return (FALSE); + case TIME_DATE: + *second = date_to_sec (line); + return (TRUE); + case TIME_DAYS: + *second = strtod (line, 0) * 86400.0; + return (TRUE); + case TIME_HOURS: + *second = strtod (line, 0) * 3600.0; + return (TRUE); + case TIME_MINUTES: + *second = strtod (line, 0) * 60.0; + return (TRUE); + case TIME_SECONDS: + *second = strtod (line, 0); + return (TRUE); + case TIME_JD: + jd = strtod (line, 0); + *second = jd_to_sec (jd); + return (TRUE); + case TIME_MJD: + jd = strtod (line, 0); + *second = mjd_to_sec (jd); + return (TRUE); + } + return (FALSE); +} + + +/**********/ +int str_to_dtime (char *line, double *second) { + + switch (chk_time (line)) { + case 0: + case TIME_JD: + case TIME_MJD: + case TIME_DATE: + return (FALSE); + case TIME_DAYS: + *second = strtod (line, 0) * 86400.0; + return (TRUE); + case TIME_HOURS: + *second = strtod (line, 0) * 3600.0; + return (TRUE); + case TIME_MINUTES: + *second = strtod (line, 0) * 60.0; + return (TRUE); + case TIME_SECONDS: + *second = strtod (line, 0); + return (TRUE); + } + return (FALSE); +} + +/**********/ +double sec_to_jd (time_t second) { + + double jd; + + jd = second/86400.0 + 2440587.5; + return (jd); +} + +/**********/ +time_t jd_to_sec (double jd) { + + time_t second; + + second = (jd - 2440587.5)*86400; + return (second); +} + +/**********/ +double sec_to_mjd (time_t second) { + + double mjd; + + mjd = second/86400.0 + 40587.0; + return (mjd); +} + +/**********/ +time_t mjd_to_sec (double mjd) { + + time_t second; + + second = (mjd - 40587.0)*86400; + return (second); +} + +/**********/ +char *sec_to_date (time_t second) { + + struct tm *gmt; + char *line; + + ALLOCATE (line, char, 64); + gmt = gmtime (&second); + sprintf (line, "%04d/%02d/%02d,%02d:%02d:%02d", 1900 + gmt[0].tm_year, gmt[0].tm_mon+1, gmt[0].tm_mday, gmt[0].tm_hour, gmt[0].tm_min, gmt[0].tm_sec); + return (line); + +} + +/***** date in format yyyy/mm/dd,hh:mm:ss *****/ +time_t date_to_sec (char *date) { + + time_t second; + double tmp, jd; + struct tm now; + char *p1, *p2, *px; + + p1 = date; + px = date + strlen(date); + bzero (&now, sizeof(now)); + + tmp = strtod (p1, &p2); + if (p2 == p1) goto escape; + now.tm_year = tmp; + if (now.tm_year > 1000) now.tm_year -= 1900; + if (now.tm_year < 50) now.tm_year += 100; + if (p2 == px) goto escape; + p1 = p2 + 1; + + tmp = strtod (p1, &p2); + if (p2 == p1) goto escape; + now.tm_mon = tmp - 1; /* mon runs from 0 - 11 */ + if (p2 == px) goto escape; + p1 = p2 + 1; + + tmp = strtod (p1, &p2); + if (p2 == p1) goto escape; + now.tm_mday = tmp; + if (p2 == px) goto escape; + p1 = p2 + 1; + + tmp = strtod (p1, &p2); + if (p2 == p1) goto escape; + p1 = p2 + 1; + now.tm_hour = tmp; + if (p2 == px) goto escape; + + tmp = strtod (p1, &p2); + if (p2 == p1) goto escape; + now.tm_min = tmp; + if (p2 == px) goto escape; + p1 = p2 + 1; + + tmp = strtod (p1, &p2); + if (p2 == p1) goto escape; + now.tm_sec = tmp; + if (p2 == px) goto escape; + p1 = p2 + 1; + + escape: + jd = now.tm_mday - 32075 + (int)(1461*(1900 + now.tm_year + 4800 + (int)(((now.tm_mon+1)-14)/12))/4) + + (int)(367*((now.tm_mon+1) - 2 - (int)(((now.tm_mon+1) - 14)/12)*12)/12) + - (int)(3*(int)((1900 + now.tm_year + 4900 + (int)(((now.tm_mon+1) - 14)/12))/100)/4) - 0.5; + + second = (jd - 2440587.5)*86400 + 3600.0*now.tm_hour + now.tm_min*60.0 + now.tm_sec; + + return (second); +} + +/***** short date in format yymmdd *****/ +time_t short_date_to_sec (char *date) { + + time_t second; + double jd; + struct tm now; + + bzero (&now, sizeof(now)); + + sscanf (date, "%2d%2d%2d", &now.tm_year, &now.tm_mon, &now.tm_mday); + + if (now.tm_year > 51) now.tm_year += 0; + if (now.tm_year < 50) now.tm_year += 100; + now.tm_mon --; /* tm_mon runs from 0 - 11 */ + + jd = now.tm_mday - 32075 + (int)(1461*(1900 + now.tm_year + 4800 + (int)(((now.tm_mon+1)-14)/12))/4) + + (int)(367*((now.tm_mon+1) - 2 - (int)(((now.tm_mon+1) - 14)/12)*12)/12) + - (int)(3*(int)((1900 + now.tm_year + 4900 + (int)(((now.tm_mon+1) - 14)/12))/100)/4) - 0.5; + + second = (jd - 2440587.5)*86400 + 3600.0*now.tm_hour + now.tm_min*60.0 + now.tm_sec; + + return (second); +} + +/**********/ +int hstgsc_hms_to_deg (double *h0, double *h1, double *d0, double *d1, char *string) { + + int flag_d0, flag_d1, flag_h0, flag_h1; + double tmp; + + *d0 = *h0 = *d1 = *h1 = 0; + + flag_h0 = dparse (h0, 1, string); + flag_h1 = dparse (h1, 4, string); + flag_d0 = dparse (d0, 7, string); + flag_d1 = dparse (d1, 9, string); + *h0 *= flag_h0; + *h1 *= flag_h1; + *d0 *= flag_d0; + *d1 *= flag_d1; + + dparse (&tmp, 2, string); + *h0 += tmp/60.0; + dparse (&tmp, 3, string); + *h0 += tmp/3600.0; + + dparse (&tmp, 5, string); + *h1 += tmp/60.0; + dparse (&tmp, 6, string); + *h1 += tmp/3600.0; + + dparse (&tmp, 8, string); + *d0 += tmp/60.0; + + dparse (&tmp, 10, string); + *d1 += tmp/60.0; + + *h0 *= 15*flag_h0; + *h1 *= 15*flag_h1; + *d0 *= flag_d0; + *d1 *= flag_d1; + + return (TRUE); +} + +/* times may be in forms as: + * 20040200450s (N seconds since 1970.0) + * 2440900.232j (julian date) + * 99/02/23,03:22:18 (date string) + * (separators may be anything except space, +, -) + * 99:02:15:12:23:30 + * 99:02:15:12h23m30s + */