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 <stdio.h>
7 #include <stdlib.h>
8 #include <assert.h>
9 
10 #include "solver.h"
11 #include "index.h"
12 #include "pquad.h"
13 #include "permutedsort.h"
14 #include "log.h"
15 
compare_n(const void * v1,const void * v2,int N)16 static int compare_n(const void* v1, const void* v2, int N) {
17     const int* u1 = v1;
18     const int* u2 = v2;
19     int i;
20     for (i=0; i<N; i++) {
21         if (u1[i] < u2[i]) return -1;
22         if (u1[i] > u2[i]) return 1;
23     }
24     return 0;
25 }
26 
compare_tri(const void * v1,const void * v2)27 static int compare_tri(const void* v1, const void* v2) {
28     return compare_n(v1, v2, 3);
29 }
compare_quad(const void * v1,const void * v2)30 static int compare_quad(const void* v1, const void* v2) {
31     return compare_n(v1, v2, 4);
32 }
compare_quint(const void * v1,const void * v2)33 static int compare_quint(const void* v1, const void* v2) {
34     return compare_n(v1, v2, 5);
35 }
36 
37 
38 
39 
40 bl* quadlist;
41 
test_try_all_codes(pquad * pq,unsigned int * fieldstars,int dimquad,solver_t * solver,double tol2)42 void test_try_all_codes(pquad* pq,
43                         unsigned int* fieldstars, int dimquad,
44                         solver_t* solver, double tol2) {
45     int sorted[dimquad];
46     int i;
47     fflush(NULL);
48     printf("test_try_all_codes: [");
49     for (i=0; i<dimquad; i++) {
50         printf("%s%i", (i?" ":""), fieldstars[i]);
51     }
52     printf("]");
53 
54     // sort AB and C[DE]...
55     memcpy(sorted, fieldstars, dimquad * sizeof(int));
56     qsort(sorted, 2, sizeof(int), compare_ints_asc);
57     qsort(sorted+2, dimquad-2, sizeof(int), compare_ints_asc);
58 
59     printf(" -> [");
60     for (i=0; i<dimquad; i++) {
61         printf("%s%i", (i?" ":""), sorted[i]);
62     }
63     printf("]\n");
64     fflush(NULL);
65 
66     bl_append(quadlist, sorted);
67 }
68 
field1()69 static starxy_t* field1() {
70     starxy_t* starxy;
71     double field[14];
72     int i=0, N;
73     // star0 A: (0,0)
74     field[i++] = 0.0;
75     field[i++] = 0.0;
76     // star1 B: (2,2)
77     field[i++] = 2.0;
78     field[i++] = 2.0;
79     // star2
80     field[i++] = -1.0;
81     field[i++] = 3.0;
82     // star3
83     field[i++] = 0.5;
84     field[i++] = 1.5;
85     // star4
86     field[i++] = 1.0;
87     field[i++] = 1.0;
88     // star5
89     field[i++] = 1.5;
90     field[i++] = 0.5;
91     // star6
92     field[i++] = 3.0;
93     field[i++] = -1.0;
94 
95     N = i/2;
96     starxy = starxy_new(N, FALSE, FALSE);
97     for (i=0; i<N; i++) {
98         starxy_setx(starxy, i, field[i*2+0]);
99         starxy_sety(starxy, i, field[i*2+1]);
100     }
101     return starxy;
102 }
103 
test1()104 void test1() {
105     int i;
106     solver_t* solver;
107     index_t index;
108     starxy_t* starxy;
109     int wanted[][4] = { { 0,1,3,4 },
110                         { 0,2,3,4 },
111                         { 1,2,3,4 },
112                         { 2,5,0,1 },
113                         { 2,5,0,3 },
114                         { 2,5,0,4 },
115                         { 2,5,1,3 },
116                         { 2,5,1,4 },
117                         { 2,5,3,4 },
118                         { 0,1,3,5 },
119                         { 0,1,4,5 },
120                         { 0,6,4,5 },
121                         { 1,6,4,5 },
122                         { 2,6,0,1 },
123                         { 2,6,0,3 },
124                         { 2,6,0,4 },
125                         { 2,6,0,5 },
126                         { 2,6,1,3 },
127                         { 2,6,1,4 },
128                         { 2,6,1,5 },
129                         { 2,6,3,4 },
130                         { 2,6,3,5 },
131                         { 2,6,4,5 },
132                         { 3,6,0,1 },
133                         { 3,6,0,4 },
134                         { 3,6,0,5 },
135                         { 3,6,1,4 },
136                         { 3,6,1,5 },
137                         { 3,6,4,5 },
138     };
139 
140     starxy = field1();
141 
142     quadlist = bl_new(16, 4*sizeof(uint));
143 
144     solver = solver_new();
145 
146     memset(&index, 0, sizeof(index_t));
147     index.index_scale_lower = 1;
148     index.index_scale_upper = 10;
149     index.dimquads = 4;
150 
151     solver->funits_lower = 0.1;
152     solver->funits_upper = 10;
153 
154     solver_add_index(solver, &index);
155     solver_set_field(solver, starxy);
156     solver_preprocess_field(solver);
157 
158     solver_run(solver);
159 
160     solver_free_field(solver);
161     solver_free(solver);
162 
163     //
164     assert(bl_size(quadlist) == (sizeof(wanted) / (4*sizeof(uint))));
165     for (i=0; i<bl_size(quadlist); i++) {
166         assert(compare_quad(bl_access(quadlist, i), wanted[i]) == 0);
167     }
168 
169     bl_free(quadlist);
170 }
171 
test2()172 void test2() {
173     int i;
174     solver_t* solver;
175     index_t index;
176     starxy_t* starxy;
177     int wanted[][3] = { { 0, 1, 3 },
178                         { 0, 1, 4 },
179                         { 0, 1, 5 },
180                         { 0, 2, 3 },
181                         { 0, 2, 4 },
182                         { 0, 3, 4 },
183                         { 0, 5, 4 },
184                         { 0, 6, 4 },
185                         { 0, 6, 5 },
186                         { 1, 2, 3 },
187                         { 1, 2, 4 },
188                         { 1, 3, 4 },
189                         { 1, 5, 4 },
190                         { 1, 6, 4 },
191                         { 1, 6, 5 },
192                         { 2, 4, 3 },
193                         { 2, 5, 0 },
194                         { 2, 5, 1 },
195                         { 2, 5, 3 },
196                         { 2, 5, 4 },
197                         { 2, 6, 0 },
198                         { 2, 6, 1 },
199                         { 2, 6, 3 },
200                         { 2, 6, 4 },
201                         { 2, 6, 5 },
202                         { 3, 5, 4 },
203                         { 3, 6, 0 },
204                         { 3, 6, 1 },
205                         { 3, 6, 4 },
206                         { 3, 6, 5 },
207                         { 4, 6, 5 },
208     };
209 
210     starxy = field1();
211     quadlist = bl_new(16, 3*sizeof(uint));
212     solver = solver_new();
213     memset(&index, 0, sizeof(index_t));
214     index.index_scale_lower = 1;
215     index.index_scale_upper = 10;
216     index.dimquads = 3;
217 
218     solver->funits_lower = 0.1;
219     solver->funits_upper = 10;
220 
221     solver_add_index(solver, &index);
222     solver_set_field(solver, starxy);
223     solver_preprocess_field(solver);
224 
225     solver_run(solver);
226 
227     solver_free_field(solver);
228     solver_free(solver);
229 
230     //
231     assert(bl_size(quadlist) == (sizeof(wanted) / (3*sizeof(uint))));
232     bl_sort(quadlist, compare_tri);
233     for (i=0; i<bl_size(quadlist); i++) {
234         assert(compare_tri(bl_access(quadlist, i), wanted[i]) == 0);
235     }
236     bl_free(quadlist);
237 }
238 
239 
240 char* OPTIONS = "v";
241 
main(int argc,char ** args)242 int main(int argc, char** args) {
243     int argchar;
244 
245     while ((argchar = getopt(argc, args, OPTIONS)) != -1)
246         switch (argchar) {
247         case 'v':
248             log_init(LOG_ALL+1);
249             break;
250         }
251 
252     test1();
253     test2();
254     return 0;
255 }
256 
257