Index: /branches/eam_branches/ipp-20100621/Ohana/src/dvomerge/include/dvomerge.h =================================================================== --- /branches/eam_branches/ipp-20100621/Ohana/src/dvomerge/include/dvomerge.h (revision 28992) +++ /branches/eam_branches/ipp-20100621/Ohana/src/dvomerge/include/dvomerge.h (revision 28993) @@ -25,4 +25,5 @@ float RADIUS; int SKY_DEPTH; +char *ALTERNATE_PHOTCODE_FILE; SkyRegion UserPatch; // used by MODE CAT @@ -66,6 +67,6 @@ int LoadCatalog PROTO((Catalog *catalog, SkyRegion *region, char *filename, char *mode)); -int merge_catalogs_new PROTO((SkyRegion *region, Catalog *output, Catalog *input)); -int merge_catalogs_old PROTO((SkyRegion *region, Catalog *output, Catalog *input, double RADIUS)); +int merge_catalogs_new PROTO((SkyRegion *region, Catalog *output, Catalog *input, int *secflitMap)); +int merge_catalogs_old PROTO((SkyRegion *region, Catalog *output, Catalog *input, double RADIUS, int *secflitMap)); off_t *init_measure_links PROTO((Average *average, off_t Naverage, Measure *measure, off_t Nmeasure)); Index: /branches/eam_branches/ipp-20100621/Ohana/src/dvomerge/src/args.c =================================================================== --- /branches/eam_branches/ipp-20100621/Ohana/src/dvomerge/src/args.c (revision 28992) +++ /branches/eam_branches/ipp-20100621/Ohana/src/dvomerge/src/args.c (revision 28993) @@ -10,4 +10,9 @@ if ((N = get_argument (*argc, argv, "-v"))) { VERBOSE = TRUE; + remove_argument (N, argc, argv); + } + if ((N = get_argument (*argc, argv, "-photcode-file"))) { + remove_argument (N, argc, argv); + ALTERNATE_PHOTCODE_FILE = strdup(argv[N]); remove_argument (N, argc, argv); } Index: /branches/eam_branches/ipp-20100621/Ohana/src/dvomerge/src/dvoconvert.c =================================================================== --- /branches/eam_branches/ipp-20100621/Ohana/src/dvomerge/src/dvoconvert.c (revision 28992) +++ /branches/eam_branches/ipp-20100621/Ohana/src/dvomerge/src/dvoconvert.c (revision 28993) @@ -85,5 +85,5 @@ continue; } - merge_catalogs_new (&outsky[0].regions[i], &outcatalog, &incatalog); + merge_catalogs_new (&outsky[0].regions[i], &outcatalog, &incatalog, NULL); dvo_catalog_unlock (&incatalog); dvo_catalog_free (&incatalog); Index: /branches/eam_branches/ipp-20100621/Ohana/src/dvomerge/src/dvomergeCreate.c =================================================================== --- /branches/eam_branches/ipp-20100621/Ohana/src/dvomerge/src/dvomergeCreate.c (revision 28992) +++ /branches/eam_branches/ipp-20100621/Ohana/src/dvomerge/src/dvomergeCreate.c (revision 28993) @@ -10,5 +10,8 @@ SkyList *inlist; Catalog incatalog, outcatalog; - char filename[256], *input1, *input2, *output; + char filename[256]; + char filename1[256]; + char *codesFilename; + char *input1, *input2, *output; IDmapType IDmap1, IDmap2; @@ -22,11 +25,65 @@ dvomergeImagesCreate (&IDmap1, input1, &IDmap2, input2, output); - // the first input defines the photcode table & db layout - sprintf (filename, "%s/Photcodes.dat", input1); + PhotCodeData *input1Photcodes = NULL; + PhotCodeData *input2Photcodes = NULL; + PhotCodeData *outputPhotcodes = NULL; + int *secfiltMap1 = NULL; + int *secfiltMap2 = NULL; + + // Read the input1 photcodes + SetPhotcodeTable(NULL); + sprintf (filename1, "%s/Photcodes.dat", input1); + if (!LoadPhotcodes (filename1, NULL, FALSE)) { + fprintf (stderr, "error loading photcode table %s\n", filename1); + exit (1); + } + input1Photcodes = GetPhotcodeTable(); + + // Read the input2 photcodes + SetPhotcodeTable(NULL); + sprintf (filename, "%s/Photcodes.dat", input2); if (!LoadPhotcodes (filename, NULL, FALSE)) { fprintf (stderr, "error loading photcode table %s\n", filename); exit (1); } - // save the photcodes in the output catdir + input2Photcodes = GetPhotcodeTable(); + + if (ALTERNATE_PHOTCODE_FILE) { + SetPhotcodeTable(NULL); + // First attempt to read file as fits table ... + if (!LoadPhotcodes(ALTERNATE_PHOTCODE_FILE, NULL, FALSE)) { + // and if that fails try it as a text file + if (!LoadPhotcodesText(ALTERNATE_PHOTCODE_FILE)) { + fprintf (stderr, "failed to read photcode file %s\n", ALTERNATE_PHOTCODE_FILE); + exit (1); + } + } + outputPhotcodes = GetPhotcodeTable(); + + secfiltMap1 = GetSecFiltMap(outputPhotcodes, input1Photcodes); + if (!secfiltMap1) { + fprintf (stderr, "failed to map %s secfilt photcodes to alternate photcodes table %s\n", input1, + ALTERNATE_PHOTCODE_FILE); + exit (1); + } + + codesFilename = ALTERNATE_PHOTCODE_FILE; + + } else { + // Use the first input file's photcodes + outputPhotcodes = input1Photcodes; + codesFilename = filename1; + // secfitMap1 can remain NULL since input1 defines the set of codes + } + + secfiltMap2 = GetSecFiltMap(outputPhotcodes, input2Photcodes); + if (!secfiltMap2) { + fprintf (stderr, "falied to map %s secfilt photcodes to output photcodes table %s\n", input2, + codesFilename); + exit (1); + } + + // save the output photcodes in the output catdir + SetPhotcodeTable(outputPhotcodes); sprintf (filename, "%s/Photcodes.dat", output); if (!check_file_access (filename, TRUE, TRUE, VERBOSE)) { @@ -64,7 +121,13 @@ if (VERBOSE) fprintf (stderr, "output: %s\n", outsky[0].regions[i].name); + if (outputPhotcodes) { + SetPhotcodeTable(outputPhotcodes); + } // load / create output catalog LoadCatalog (&outcatalog, &outsky[0].regions[i], outsky[0].filename[i], "w"); + if (input1Photcodes) { + SetPhotcodeTable(input1Photcodes); + } // combine only tables at equal or larger depth @@ -84,9 +147,13 @@ } dvo_update_image_IDs (&IDmap1, &incatalog); - merge_catalogs_new (&outsky[0].regions[i], &outcatalog, &incatalog); + merge_catalogs_new (&outsky[0].regions[i], &outcatalog, &incatalog, secfiltMap1); dvo_catalog_unlock (&incatalog); dvo_catalog_free (&incatalog); } SkyListFree (inlist); + + if (input2Photcodes) { + SetPhotcodeTable(input2Photcodes); + } // load in all of the tables from input2 for this region @@ -95,5 +162,5 @@ if (VERBOSE) fprintf (stderr, "input 2: %s\n", inlist[0].regions[j][0].name); - // load input catalog (1) + // load input catalog (2) LoadCatalog (&incatalog, inlist[0].regions[j], inlist[0].filename[j], "r"); @@ -105,5 +172,5 @@ } dvo_update_image_IDs (&IDmap2, &incatalog); - merge_catalogs_old (&outsky[0].regions[i], &outcatalog, &incatalog, RADIUS); + merge_catalogs_old (&outsky[0].regions[i], &outcatalog, &incatalog, RADIUS, secfiltMap2); dvo_catalog_unlock (&incatalog); dvo_catalog_free (&incatalog); Index: /branches/eam_branches/ipp-20100621/Ohana/src/dvomerge/src/dvomergeUpdate.c =================================================================== --- /branches/eam_branches/ipp-20100621/Ohana/src/dvomerge/src/dvomergeUpdate.c (revision 28992) +++ /branches/eam_branches/ipp-20100621/Ohana/src/dvomerge/src/dvomergeUpdate.c (revision 28993) @@ -10,4 +10,7 @@ char filename[256], *input, *output; IDmapType IDmap; + PhotCodeData *inputPhotcodes; + PhotCodeData *outputPhotcodes; + int *secfiltMap = NULL; if (strcasecmp (argv[2], "into")) dvomerge_usage(); @@ -16,23 +19,44 @@ output = argv[3]; + if (ALTERNATE_PHOTCODE_FILE) { + fprintf (stderr, "cannot specify photcodes when merging into an existing catdir\n"); + exit (1); + } + + SetPhotcodeTable(NULL); + sprintf (filename, "%s/Photcodes.dat", input); + if (!LoadPhotcodes (filename, NULL, FALSE)) { + fprintf (stderr, "error reading input database directory %s\n", input); + exit (1); + } + inputPhotcodes = GetPhotcodeTable(); + // since we are merging the input db into the output db, the output defines the photcode // table & db layout but, this requires the output to exist. if it does not, instead use the // input. + + SetPhotcodeTable(NULL); sprintf (filename, "%s/Photcodes.dat", output); - if (!LoadPhotcodes (filename, NULL, FALSE)) { - sprintf (filename, "%s/Photcodes.dat", input); - if (!LoadPhotcodes (filename, NULL, FALSE)) { - fprintf (stderr, "error loading photcode table: tried %s/Photcodes.dat and %s/Photcodes.dat\n", output, input); + if (LoadPhotcodes (filename, NULL, FALSE)) { + + outputPhotcodes = GetPhotcodeTable(); + + secfiltMap = GetSecFiltMap(outputPhotcodes, inputPhotcodes); + if (!secfiltMap) { + fprintf (stderr, "failed to map input secfilt photcodes to output photcodes table\n"); exit (1); } + } else { + fprintf(stderr, "%s not found using photcodes from %s/Photcodes.dat\n", filename, input); + outputPhotcodes = inputPhotcodes; if (!check_dir_access (output, VERBOSE)) { fprintf (stderr, "error creating output database directory %s\n", output); exit (1); - } - sprintf (filename, "%s/Photcodes.dat", output); + } + SetPhotcodeTable(inputPhotcodes); if (!SavePhotcodesFITS (filename)) { fprintf (stderr, "error saving photcode table in %s/Photcodes.dat\n", output); exit (1); - } + } } @@ -66,4 +90,5 @@ if (VERBOSE) fprintf (stderr, "input: %s\n", inlist[0].regions[i][0].name); + SetPhotcodeTable(inputPhotcodes); // load / create output catalog (if catalog does not exist, it will be created) LoadCatalog (&incatalog, &inlist[0].regions[i][0], inlist[0].filename[i], "r"); @@ -86,8 +111,9 @@ // load input catalog + SetPhotcodeTable(outputPhotcodes); LoadCatalog (&outcatalog, outlist[0].regions[j], outlist[0].filename[j], "w"); dvo_update_image_IDs (&IDmap, &incatalog); - merge_catalogs_old (&outsky[0].regions[j], &outcatalog, &incatalog, RADIUS); + merge_catalogs_old (&outsky[0].regions[j], &outcatalog, &incatalog, RADIUS, secfiltMap); outcatalog.catflags = LOAD_AVES | LOAD_MEAS | LOAD_MISS | LOAD_SECF; Index: /branches/eam_branches/ipp-20100621/Ohana/src/dvomerge/src/merge_catalogs_new.c =================================================================== --- /branches/eam_branches/ipp-20100621/Ohana/src/dvomerge/src/merge_catalogs_new.c (revision 28992) +++ /branches/eam_branches/ipp-20100621/Ohana/src/dvomerge/src/merge_catalogs_new.c (revision 28993) @@ -7,13 +7,17 @@ // merge the input data into the output catalog // input entries always define new objects -int merge_catalogs_new (SkyRegion *region, Catalog *output, Catalog *input) { + +#define notyet 1 + +int merge_catalogs_new (SkyRegion *region, Catalog *output, Catalog *input, int *secfiltMap) { off_t i, j, offset; - off_t NAVERAGE, NMEASURE, Naverage, Nmeasure, Nsecfilt, Nm; + off_t NAVERAGE, NMEASURE, Naverage, Nmeasure, NsecfiltIn, NsecfiltOut, Nm; Naverage = output[0].Naverage; Nmeasure = output[0].Nmeasure; - Nsecfilt = output[0].Nsecfilt; - assert (output[0].Nsecfilt == input[0].Nsecfilt); + + NsecfiltOut = output[0].Nsecfilt; + NsecfiltIn = input[0].Nsecfilt; NAVERAGE = Naverage + 100; @@ -21,5 +25,5 @@ REALLOCATE (output[0].average, Average, NAVERAGE); - REALLOCATE (output[0].secfilt, SecFilt, NAVERAGE*Nsecfilt); + REALLOCATE (output[0].secfilt, SecFilt, NAVERAGE*NsecfiltOut); REALLOCATE (output[0].measure, Measure, NMEASURE); @@ -31,7 +35,15 @@ output[0].average[Naverage].measureOffset = Nmeasure; - // XXX input[0].Nsecfilt == output[0].Nsecfilt ??? - for (j = 0; j < Nsecfilt; j++) { - output[0].secfilt[Nsecfilt*Naverage+j] = input[0].secfilt[Nsecfilt*i+j]; + if (secfiltMap) { + // the output secfilt photcodes may be in a different place in the output table + assert (output[0].Nsecfilt >= input[0].Nsecfilt); + for (j = 0; j < NsecfiltIn; j++) { + output[0].secfilt[NsecfiltOut*Naverage + secfiltMap[j]] = input[0].secfilt[NsecfiltIn*i+j]; + } + } else { + assert (output[0].Nsecfilt == input[0].Nsecfilt); + for (j = 0; j < NsecfiltIn; j++) { + output[0].secfilt[NsecfiltOut*Naverage+j] = input[0].secfilt[NsecfiltIn*i+j]; + } } @@ -55,13 +67,13 @@ NAVERAGE += 50; REALLOCATE (output[0].average, Average, NAVERAGE); - REALLOCATE (output[0].secfilt, SecFilt, NAVERAGE*Nsecfilt); + REALLOCATE (output[0].secfilt, SecFilt, NAVERAGE*NsecfiltOut); } } REALLOCATE (output[0].average, Average, MAX (Naverage, 1)); REALLOCATE (output[0].measure, Measure, MAX (Nmeasure, 1)); - REALLOCATE (output[0].secfilt, SecFilt, Nsecfilt*MAX (Naverage, 1)); + REALLOCATE (output[0].secfilt, SecFilt, NsecfiltOut*MAX (Naverage, 1)); output[0].Naverage = Naverage; output[0].Nmeasure = Nmeasure; - output[0].Nsecf_mem = Naverage * Nsecfilt; + output[0].Nsecf_mem = Naverage * NsecfiltOut; // If we are using dvomergeCreate to split an existing catalog, then the max objID in Index: /branches/eam_branches/ipp-20100621/Ohana/src/dvomerge/src/merge_catalogs_old.c =================================================================== --- /branches/eam_branches/ipp-20100621/Ohana/src/dvomerge/src/merge_catalogs_old.c (revision 28992) +++ /branches/eam_branches/ipp-20100621/Ohana/src/dvomerge/src/merge_catalogs_old.c (revision 28993) @@ -16,5 +16,5 @@ // merge the input data into the output catalog // input entries are matched by position to the new objects -int merge_catalogs_old (SkyRegion *region, Catalog *output, Catalog *input, double RADIUS) { +int merge_catalogs_old (SkyRegion *region, Catalog *output, Catalog *input, double RADIUS, int *secfiltMap) { off_t i, j, k, Nin, offset, J, Jmin, status, Nstars; @@ -24,5 +24,6 @@ off_t *N1, *N2, *next_meas; off_t Nave, NAVE, Nmeas, NMEAS, Nmatch; - int Nsecfilt; + int NsecfiltIn; + int NsecfiltOut; unsigned int objID, catID; Coords tcoords; @@ -31,5 +32,7 @@ // gettimeofday (&start, (void *) NULL); - Nsecfilt = output[0].Nsecfilt; + NsecfiltOut = output[0].Nsecfilt; + NsecfiltIn = input[0].Nsecfilt; + assert (secfiltMap || (NsecfiltOut == NsecfiltIn)); /** allocate local arrays (stars) **/ @@ -228,8 +231,20 @@ // update the average properties to reflect the incoming entries: // if the original value is NAN but the input value is not, accept the input: - for (k = 0; k < Nsecfilt; k++) { - if ( isfinite(output[0].secfilt[n*Nsecfilt+k].M)) continue; - if (!isfinite( input[0].secfilt[N*Nsecfilt+k].M)) continue; - output[0].secfilt[n*Nsecfilt+k] = input[0].secfilt[N*Nsecfilt+k]; + if (!secfiltMap) { + // secfilt tables of input and output are the same + for (k = 0; k < NsecfiltIn; k++) { + if ( isfinite(output[0].secfilt[n*NsecfiltIn+k].M)) continue; + if (!isfinite( input[0].secfilt[N*NsecfiltIn+k].M)) continue; + output[0].secfilt[n*NsecfiltIn+k] = input[0].secfilt[N*NsecfiltIn+k]; + } + } else { + for (k = 0; k < NsecfiltIn; k++) { + // index for this entry in output's secfilt list + int outputIndex = n * NsecfiltOut + secfiltMap[k]; + + if ( isfinite(output[0].secfilt[outputIndex].M)) continue; + if (!isfinite( input[0].secfilt[N*NsecfiltIn+k].M)) continue; + output[0].secfilt[outputIndex] = input[0].secfilt[N*NsecfiltIn+k]; + } } @@ -255,5 +270,5 @@ NAVE = Nave + 1000; REALLOCATE (output[0].average, Average, NAVE); - REALLOCATE (output[0].secfilt, SecFilt, NAVE*output[0].Nsecfilt); + REALLOCATE (output[0].secfilt, SecFilt, NAVE*NsecfiltOut); } @@ -302,21 +317,28 @@ objID ++; - for (j = 0; j < Nsecfilt; j++) { - if (isfinite(input[0].secfilt[N*Nsecfilt+j].M)) { - output[0].secfilt[Nave*Nsecfilt+j].M = input[0].secfilt[N*Nsecfilt+j].M; - output[0].secfilt[Nave*Nsecfilt+j].dM = input[0].secfilt[N*Nsecfilt+j].dM; - output[0].secfilt[Nave*Nsecfilt+j].Xm = input[0].secfilt[N*Nsecfilt+j].Xm; - output[0].secfilt[Nave*Nsecfilt+j].M_20 = input[0].secfilt[N*Nsecfilt+j].M_20; - output[0].secfilt[Nave*Nsecfilt+j].M_80 = input[0].secfilt[N*Nsecfilt+j].M_80; - output[0].secfilt[Nave*Nsecfilt+j].Ncode = input[0].secfilt[N*Nsecfilt+j].Ncode; - output[0].secfilt[Nave*Nsecfilt+j].Nused = input[0].secfilt[N*Nsecfilt+j].Nused; + for (j = 0; j < NsecfiltIn; j++) { + int outputIndex; + if (secfiltMap) { + outputIndex = (Nave * NsecfiltOut) + secfiltMap[j]; } else { - output[0].secfilt[Nave*Nsecfilt+j].M = NAN; - output[0].secfilt[Nave*Nsecfilt+j].dM = NAN; - output[0].secfilt[Nave*Nsecfilt+j].Xm = NAN_S_SHORT; - output[0].secfilt[Nave*Nsecfilt+j].M_20 = NAN_S_SHORT; - output[0].secfilt[Nave*Nsecfilt+j].M_80 = NAN_S_SHORT; - output[0].secfilt[Nave*Nsecfilt+j].Ncode = 0; - output[0].secfilt[Nave*Nsecfilt+j].Nused = 0; + outputIndex = (Nave * NsecfiltIn) + j; + } + + if (isfinite(input[0].secfilt[N*NsecfiltIn+j].M)) { + output[0].secfilt[outputIndex].M = input[0].secfilt[N*NsecfiltIn+j].M; + output[0].secfilt[outputIndex].dM = input[0].secfilt[N*NsecfiltIn+j].dM; + output[0].secfilt[outputIndex].Xm = input[0].secfilt[N*NsecfiltIn+j].Xm; + output[0].secfilt[outputIndex].M_20 = input[0].secfilt[N*NsecfiltIn+j].M_20; + output[0].secfilt[outputIndex].M_80 = input[0].secfilt[N*NsecfiltIn+j].M_80; + output[0].secfilt[outputIndex].Ncode = input[0].secfilt[N*NsecfiltIn+j].Ncode; + output[0].secfilt[outputIndex].Nused = input[0].secfilt[N*NsecfiltIn+j].Nused; + } else { + output[0].secfilt[outputIndex].M = NAN; + output[0].secfilt[outputIndex].dM = NAN; + output[0].secfilt[outputIndex].Xm = NAN_S_SHORT; + output[0].secfilt[outputIndex].M_20 = NAN_S_SHORT; + output[0].secfilt[outputIndex].M_80 = NAN_S_SHORT; + output[0].secfilt[outputIndex].Ncode = 0; + output[0].secfilt[outputIndex].Nused = 0; } } @@ -376,5 +398,5 @@ output[0].Naverage = Nave; output[0].Nmeasure = Nmeas; - output[0].Nsecf_mem = Nave*Nsecfilt; + output[0].Nsecf_mem = Nave*NsecfiltOut; if (VERBOSE) fprintf (stderr, "Nstars, Nave, Nmeas: "OFF_T_FMT" "OFF_T_FMT" "OFF_T_FMT", ("OFF_T_FMT" matches)\n", Nstars, Nave, Nmeas, Nmatch);