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 <string.h>
7 #include <stdlib.h>
8 #include <stdio.h>
9 #include <math.h>
10 #include <assert.h>
11
12 #include "index.h"
13 #include "quadfile.h"
14 #include "kdtree.h"
15 #include "starutil.h"
16 #include "mathutil.h"
17 #include "fitsioutils.h"
18 #include "bl.h"
19 #include "starkd.h"
20 #include "boilerplate.h"
21 #include "rdlist.h"
22
23 #define OPTIONS "hr:R"
24
25
print_help(char * progname)26 void print_help(char* progname)
27 {
28 BOILERPLATE_HELP_HEADER(stderr);
29 fprintf(stderr, "Usage: %s\n"
30 " -r <rdls-output-file>\n"
31 " [-R]: add quad radius columns to the RDLS file.\n"
32 " [-h]: help\n"
33 " <base-name> [<base-name> ...]\n\n"
34 "Reads index (.quad and .skdt or merged index) files.\n"
35 "Writes an RDLS containing the quad centers (midpoints of AB), one field per input file.\n\n",
36 progname);
37 }
38
main(int argc,char ** args)39 int main(int argc, char** args) {
40 int argchar;
41 char* basename;
42 char* outfn = NULL;
43 index_t* index;
44 quadfile* qf;
45 rdlist_t* rdls;
46 startree_t* skdt = NULL;
47 anbool addradius = FALSE;
48 int i;
49 int radcolumn = -1;
50
51 while ((argchar = getopt (argc, args, OPTIONS)) != -1)
52 switch (argchar) {
53 case 'R':
54 addradius = TRUE;
55 break;
56 case 'r':
57 outfn = optarg;
58 break;
59 case 'h':
60 print_help(args[0]);
61 exit(0);
62 }
63
64 if (!outfn || (optind == argc)) {
65 print_help(args[0]);
66 exit(-1);
67 }
68
69 rdls = rdlist_open_for_writing(outfn);
70 if (!rdls) {
71 fprintf(stderr, "Failed to open RDLS file %s for output.\n", outfn);
72 exit(-1);
73 }
74 if (rdlist_write_primary_header(rdls)) {
75 fprintf(stderr, "Failed to write RDLS header.\n");
76 exit(-1);
77 }
78
79 if (addradius) {
80 radcolumn = rdlist_add_tagalong_column(rdls, fitscolumn_double_type(),
81 1, fitscolumn_double_type(),
82 "QUADRADIUS", "deg");
83 }
84
85 for (; optind<argc; optind++) {
86 //int Nstars;
87 int dimquads;
88
89 basename = args[optind];
90 printf("Reading files with basename %s\n", basename);
91
92 index = index_load(basename, 0, NULL);
93 if (!index) {
94 fprintf(stderr, "Failed to read index with base name \"%s\"\n", basename);
95 exit(-1);
96 }
97
98 qf = index->quads;
99 skdt = index->starkd;
100 //Nstars = startree_N(skdt);
101
102 if (rdlist_write_header(rdls)) {
103 fprintf(stderr, "Failed to write new RDLS field header.\n");
104 exit(-1);
105 }
106
107 dimquads = quadfile_dimquads(qf);
108
109 printf("Reading quads...\n");
110 for (i=0; i<qf->numquads; i++) {
111 unsigned int stars[dimquads];
112 double axyz[3], bxyz[3];
113 double midab[3];
114 double radec[2];
115 if (!(i % 200000)) {
116 printf(".");
117 fflush(stdout);
118 }
119 quadfile_get_stars(qf, i, stars);
120 startree_get(skdt, stars[0], axyz);
121 startree_get(skdt, stars[1], bxyz);
122 star_midpoint(midab, axyz, bxyz);
123 xyzarr2radecdegarr(midab, radec);
124
125 if (rdlist_write_one_radec(rdls, radec[0], radec[1])) {
126 fprintf(stderr, "Failed to write a RA,Dec entry.\n");
127 exit(-1);
128 }
129
130 if (addradius) {
131 double rad = arcsec2deg(distsq2arcsec(distsq(midab, axyz, 3)));
132 if (rdlist_write_tagalong_column(rdls, radcolumn, i, 1, &rad, 0)) {
133 fprintf(stderr, "Failed to write quad radius.\n");
134 exit(-1);
135 }
136 }
137 }
138 printf("\n");
139
140 index_close(index);
141
142 if (rdlist_fix_header(rdls)) {
143 fprintf(stderr, "Failed to fix RDLS field header.\n");
144 exit(-1);
145 }
146
147 rdlist_next_field(rdls);
148 }
149
150 if (rdlist_fix_primary_header(rdls) ||
151 rdlist_close(rdls)) {
152 fprintf(stderr, "Failed to close RDLS file.\n");
153 exit(-1);
154 }
155
156 return 0;
157 }
158