1 /*
2 # This file is part of the Astrometry.net suite.
3 # Licensed under a 3-clause BSD style license - see LICENSE
4 */
5
6 #include <math.h>
7 #include <stdio.h>
8 #include <stdlib.h>
9 #include <unistd.h>
10 #include <string.h>
11 #include <errno.h>
12 #include <sys/types.h>
13 #include <sys/mman.h>
14 #include <assert.h>
15
16 #include "nomad.h"
17 #include "nomad-fits.h"
18 #include "qfits_header.h"
19 #include "healpix.h"
20 #include "starutil.h"
21 #include "fitsioutils.h"
22 #include "boilerplate.h"
23
24 #define OPTIONS "ho:N:"
25
print_help(char * progname)26 void print_help(char* progname) {
27 BOILERPLATE_HELP_HEADER(stdout);
28 printf("\nUsage:\n"
29 " %s -o <output-filename-template> [eg, nomad_%%03i.fits]\n"
30 " [-N <healpix-nside>] (default = 9)\n"
31 " <input-file> [<input-file> ...]\n"
32 "\n"
33 "The output-filename-template should contain a \"printf\" sequence like \"%%03i\";\n"
34 "we use \"sprintf(filename, output-filename-template, healpix)\" to determine the filename\n"
35 "to be used for each healpix.\n\n"
36 "\nNOTE: WE ASSUME THE NOMAD FILES ARE GIVEN ON THE COMMAND LINE IN ORDER: 000/b0000.cat, 000/b0001.cat, etc.\n\n\n",
37 progname);
38 }
39
40
main(int argc,char ** args)41 int main(int argc, char** args) {
42 char* outfn = NULL;
43 int c;
44 int startoptind;
45 int nrecords, nfiles;
46 int Nside = 9;
47
48 nomad_fits** nomads;
49
50 int i, HP;
51 int slicecounts[1800];
52
53 while ((c = getopt(argc, args, OPTIONS)) != -1) {
54 switch (c) {
55 case '?':
56 case 'h':
57 print_help(args[0]);
58 exit(0);
59 case 'N':
60 Nside = atoi(optarg);
61 break;
62 case 'o':
63 outfn = optarg;
64 break;
65 }
66 }
67
68 if (!outfn || (optind == argc)) {
69 print_help(args[0]);
70 exit(-1);
71 }
72
73 if (Nside < 1) {
74 fprintf(stderr, "Nside must be >= 1.\n");
75 print_help(args[0]);
76 exit(-1);
77 }
78
79 HP = 12 * Nside * Nside;
80
81 printf("Nside = %i, using %i healpixes.\n", Nside, HP);
82
83 nomads = calloc(HP, sizeof(nomad_fits*));
84
85 memset(slicecounts, 0, 1800 * sizeof(uint));
86
87 nrecords = 0;
88 nfiles = 0;
89
90 printf("Reading NOMAD files... ");
91 fflush(stdout);
92
93 startoptind = optind;
94 for (; optind<argc; optind++) {
95 char* infn;
96 FILE* fid;
97 unsigned char* map;
98 size_t map_size;
99 int i;
100 int lastgrass;
101
102 infn = args[optind];
103 fid = fopen(infn, "rb");
104 if (!fid) {
105 fprintf(stderr, "Couldn't open input file %s: %s\n", infn, strerror(errno));
106 exit(-1);
107 }
108
109 if ((optind > startoptind) && ((optind - startoptind) % 100 == 0)) {
110 printf("\nReading file %i of %i: %s\n", optind - startoptind,
111 argc - startoptind, infn);
112 }
113
114 if (fseeko(fid, 0, SEEK_END)) {
115 fprintf(stderr, "Couldn't seek to end of input file %s: %s\n", infn, strerror(errno));
116 exit(-1);
117 }
118
119 map_size = ftello(fid);
120
121 fseeko(fid, 0, SEEK_SET);
122
123 map = mmap(NULL, map_size, PROT_READ, MAP_SHARED, fileno(fid), 0);
124 if (map == MAP_FAILED) {
125 fprintf(stderr, "Couldn't mmap input file %s: %s\n", infn, strerror(errno));
126 exit(-1);
127 }
128 fclose(fid);
129
130 if (map_size % NOMAD_RECORD_SIZE) {
131 fprintf(stderr, "Warning, input file %s has size %u which is not divisible into %i-byte records.\n",
132 infn, (unsigned int)map_size, NOMAD_RECORD_SIZE);
133 }
134
135 printf("File %i of %i: %s: %i records.\n", optind - startoptind, argc - startoptind, infn, (int)(map_size / NOMAD_RECORD_SIZE));
136
137 lastgrass = 0;
138 for (i=0; i<map_size; i+=NOMAD_RECORD_SIZE) {
139 nomad_entry entry;
140 int hp;
141 int slice;
142
143 if ((i * 80 / map_size) != lastgrass) {
144 printf(".");
145 fflush(stdout);
146 lastgrass = i * 80 / map_size;
147 }
148
149 if (nomad_parse_entry(&entry, map + i)) {
150 fprintf(stderr, "Failed to parse NOMAD entry: offset %i in file %s.\n",
151 i, infn);
152 exit(-1);
153 }
154
155 // compute the nomad_id based on its DEC slice and index.
156 slice = (int)(10.0 * (entry.dec + 90.0));
157 assert(slice < 1800);
158 assert((slicecounts[slice] & 0xffe00000) == 0);
159 entry.nomad_id = (slice << 21) | (slicecounts[slice]);
160 slicecounts[slice]++;
161
162 hp = radectohealpix(deg2rad(entry.ra), deg2rad(entry.dec), Nside);
163
164 if (!nomads[hp]) {
165 char fn[256];
166 sprintf(fn, outfn, hp);
167 nomads[hp] = nomad_fits_open_for_writing(fn);
168 if (!nomads[hp]) {
169 fprintf(stderr, "Failed to initialized FITS file %i (filename %s).\n", hp, fn);
170 exit(-1);
171 }
172
173 // header remarks...
174 fits_header_add_int(nomads[hp]->header, "HEALPIX", hp, "The healpix number of this catalog.");
175 fits_header_add_int(nomads[hp]->header, "NSIDE", Nside, "The healpix resolution.");
176 BOILERPLATE_ADD_FITS_HEADERS(nomads[hp]->header);
177 qfits_header_add(nomads[hp]->header, "HISTORY", "Created by the program \"nomadtofits\"", NULL, NULL);
178 qfits_header_add(nomads[hp]->header, "HISTORY", "nomadtofits command line:", NULL, NULL);
179 fits_add_args(nomads[hp]->header, args, argc);
180 qfits_header_add(nomads[hp]->header, "HISTORY", "(end of command line)", NULL, NULL);
181
182 if (nomad_fits_write_headers(nomads[hp])) {
183 fprintf(stderr, "Failed to write header for FITS file %s.\n", fn);
184 exit(-1);
185 }
186 }
187
188 if (nomad_fits_write_entry(nomads[hp], &entry)) {
189 fprintf(stderr, "Failed to write FITS entry.\n");
190 exit(-1);
191 }
192
193 nrecords++;
194 }
195
196 munmap(map, map_size);
197
198 nfiles++;
199 printf("\n");
200 }
201 printf("\n");
202
203 // close all the files...
204 for (i=0; i<HP; i++) {
205 if (!nomads[i])
206 continue;
207 if (nomad_fits_fix_headers(nomads[i]) ||
208 nomad_fits_close(nomads[i])) {
209 fprintf(stderr, "Failed to close file %i: %s\n", i, strerror(errno));
210 }
211 }
212
213 printf("Read %u files, %u records.\n", nfiles, nrecords);
214
215 free(nomads);
216
217 return 0;
218 }
219
220