1 /*
2  * Copyright (c) 2012-2021, The OSKAR Developers.
3  * See the LICENSE file at the top-level directory of this distribution.
4  */
5 
6 #include <gtest/gtest.h>
7 
8 #include "mem/oskar_mem.h"
9 #include "telescope/oskar_telescope.h"
10 #include "utility/oskar_dir.h"
11 #include "utility/oskar_get_error_string.h"
12 
13 #include <cstdio>
14 #include <cstdlib>
15 #include <vector>
16 
17 using std::vector;
18 
19 static void generate_noisy_telescope(const char* dir, int num_stations,
20         const vector<double>& freqs, const vector<double>& noise);
21 
TEST(telescope_model_load_save,test_0_level)22 TEST(telescope_model_load_save, test_0_level)
23 {
24     int err = 0;
25     const char* tm = "temp_test_telescope_0_level";
26     double longitude_rad = 0.1;
27     double latitude_rad = 0.5;
28     double altitude_m = 1.0;
29 
30     {
31         FILE* f = 0;
32         int num_stations = 10;
33         char* path = 0;
34 
35         // Create a telescope model directory.
36         oskar_dir_mkpath(tm);
37 
38         // Write position file.
39         path = oskar_dir_get_path(tm, "position.txt");
40         f = fopen(path, "w");
41         fprintf(f, "0.0, 0.0\n");
42         fclose(f);
43         free(path);
44 
45         // Write the top-level layout file.
46         path = oskar_dir_get_path(tm, "layout.txt");
47         f = fopen(path, "w");
48         for (int i = 0; i < num_stations; ++i)
49         {
50             fprintf(f, "%.1f, %.1f, %.1f\n", i * 10.0, i * 20.0, i * 30.0);
51         }
52         fclose(f);
53         free(path);
54     }
55 
56     // Load it back again.
57     {
58         oskar_Telescope* telescope = oskar_telescope_create(OSKAR_DOUBLE,
59                         OSKAR_CPU, 0, &err);
60         oskar_telescope_set_position(telescope,
61                 longitude_rad, latitude_rad, altitude_m);
62         oskar_telescope_set_enable_numerical_patterns(telescope, 0);
63         oskar_telescope_load(telescope, tm, NULL, &err);
64         ASSERT_EQ(0, err) << oskar_get_error_string(err);
65         oskar_telescope_free(telescope, &err);
66     }
67 
68     // Remove test directory.
69     oskar_dir_remove(tm);
70 }
71 
72 
TEST(telescope_model_load_save,test_1_level)73 TEST(telescope_model_load_save, test_1_level)
74 {
75     int err = 0;
76     const char* tm = "temp_test_telescope_1_level";
77 
78     // Create a telescope model.
79     int num_stations = 10;
80     int num_elements = 20;
81     double longitude_rad = 0.1;
82     double latitude_rad = 0.5;
83     double altitude_m = 1.0;
84     oskar_Telescope* telescope = oskar_telescope_create(OSKAR_SINGLE,
85             OSKAR_CPU, num_stations, &err);
86     oskar_telescope_resize_station_array(telescope, num_stations, &err);
87     oskar_telescope_set_unique_stations(telescope, 1, &err);
88     oskar_telescope_set_position(telescope, longitude_rad,
89             latitude_rad, altitude_m);
90 
91     // Populate a telescope model.
92     for (int i = 0; i < num_stations; ++i)
93     {
94         double xyz[3];
95         xyz[0] = 1.0 * i;
96         xyz[1] = 2.0 * i;
97         xyz[2] = 3.0 * i;
98         oskar_Station* st = oskar_telescope_station(telescope, i);
99         oskar_telescope_set_station_coords(telescope, i, xyz,
100                 xyz, xyz, xyz, xyz, &err);
101         ASSERT_EQ(0, err) << oskar_get_error_string(err);
102         oskar_station_resize(st, num_elements, &err);
103         ASSERT_EQ(0, err) << oskar_get_error_string(err);
104 
105         for (int j = 0; j < num_elements; ++j)
106         {
107             xyz[0] = 10.0 * i + 1.0 * j;
108             xyz[1] = 20.0 * i + 1.0 * j;
109             xyz[2] = 30.0 * i + 1.0 * j;
110             oskar_station_set_element_coords(st, 0, j, xyz, xyz, &err);
111             ASSERT_EQ(0, err) << oskar_get_error_string(err);
112         }
113     }
114 
115     // Save the telescope model.
116     oskar_telescope_save(telescope, tm, &err);
117     ASSERT_EQ(0, err) << oskar_get_error_string(err);
118 
119     // Load it back again.
120     oskar_Telescope* telescope2 = oskar_telescope_create(OSKAR_SINGLE,
121                     OSKAR_CPU, 0, &err);
122     oskar_telescope_set_position(telescope2,
123             longitude_rad, latitude_rad, altitude_m);
124     oskar_telescope_set_enable_numerical_patterns(telescope, 0);
125     oskar_telescope_load(telescope2, tm, NULL, &err);
126     ASSERT_EQ(0, err) << oskar_get_error_string(err);
127 
128     // Check the contents.
129     ASSERT_EQ(oskar_telescope_num_stations(telescope),
130             oskar_telescope_num_stations(telescope2));
131     EXPECT_NEAR(oskar_telescope_lon_rad(telescope),
132             oskar_telescope_lon_rad(telescope2), 1e-10);
133     EXPECT_NEAR(oskar_telescope_lat_rad(telescope),
134             oskar_telescope_lat_rad(telescope2), 1e-10);
135     EXPECT_NEAR(oskar_telescope_alt_metres(telescope),
136             oskar_telescope_alt_metres(telescope2), 1e-10);
137 
138     double max_ = 0.0, avg_ = 0.0;
139     for (int dim = 0; dim < 3; dim++)
140     {
141         oskar_mem_evaluate_relative_error(
142                 oskar_telescope_station_true_enu_metres(telescope, dim),
143                 oskar_telescope_station_true_enu_metres(telescope2, dim),
144                 0, &max_, &avg_, 0, &err);
145         ASSERT_EQ(0, err) << oskar_get_error_string(err);
146         EXPECT_LT(max_, 1e-5);
147         EXPECT_LT(avg_, 1e-5);
148     }
149     for (int i = 0; i < num_stations; ++i)
150     {
151         oskar_Station* st1 = oskar_telescope_station(telescope, i);
152         oskar_Station* st2 = oskar_telescope_station(telescope2, i);
153         for (int dim = 0; dim < 3; dim++)
154         {
155             oskar_mem_evaluate_relative_error(
156                     oskar_station_element_measured_enu_metres(st1, 0, dim),
157                     oskar_station_element_measured_enu_metres(st2, 0, dim),
158                     0, &max_, &avg_, 0, &err);
159             ASSERT_EQ(0, err) << oskar_get_error_string(err);
160             EXPECT_LT(max_, 1e-5);
161             EXPECT_LT(avg_, 1e-5);
162         }
163     }
164 
165     // Remove test directory.
166     oskar_dir_remove(tm);
167 
168     // Free models.
169     oskar_telescope_free(telescope, &err);
170     oskar_telescope_free(telescope2, &err);
171 }
172 
173 
TEST(telescope_model_load_save,test_2_level)174 TEST(telescope_model_load_save, test_2_level)
175 {
176     int err = 0;
177     const char* tm = "temp_test_telescope_2_level";
178 
179     // Create a telescope model.
180     int num_stations = 3;
181     int num_tiles = 4;
182     int num_elements = 8;
183     double longitude_rad = 0.1;
184     double latitude_rad = 0.5;
185     double altitude_m = 1.0;
186     oskar_Telescope* telescope = oskar_telescope_create(OSKAR_SINGLE,
187             OSKAR_CPU, num_stations, &err);
188     oskar_telescope_resize_station_array(telescope, num_stations, &err);
189     oskar_telescope_set_unique_stations(telescope, 1, &err);
190     oskar_telescope_set_position(telescope, longitude_rad,
191             latitude_rad, altitude_m);
192 
193     // Populate a telescope model.
194     for (int i = 0; i < num_stations; ++i)
195     {
196         double xyz[3];
197         xyz[0] = 1.0 * i;
198         xyz[1] = 2.0 * i;
199         xyz[2] = 3.0 * i;
200         oskar_Station* st = oskar_telescope_station(telescope, i);
201         oskar_telescope_set_station_coords(telescope, i, xyz,
202                 xyz, xyz, xyz, xyz, &err);
203         ASSERT_EQ(0, err) << oskar_get_error_string(err);
204         oskar_station_resize(st, num_tiles, &err);
205         ASSERT_EQ(0, err) << oskar_get_error_string(err);
206         oskar_station_create_child_stations(st, &err);
207         ASSERT_EQ(0, err) << oskar_get_error_string(err);
208 
209         for (int j = 0; j < num_tiles; ++j)
210         {
211             xyz[0] = 10.0 * i + 1.0 * j;
212             xyz[1] = 20.0 * i + 1.0 * j;
213             xyz[2] = 30.0 * i + 1.0 * j;
214             oskar_station_set_element_coords(st, 0, j, xyz, xyz, &err);
215             ASSERT_EQ(0, err) << oskar_get_error_string(err);
216             oskar_station_resize(oskar_station_child(st, j), num_elements, &err);
217             ASSERT_EQ(0, err) << oskar_get_error_string(err);
218 
219             for (int k = 0; k < num_elements; ++k)
220             {
221                 xyz[0] = 100.0 * i + 10.0 * j + 1.0 * k;
222                 xyz[1] = 200.0 * i + 10.0 * j + 1.0 * k;
223                 xyz[2] = 300.0 * i + 10.0 * j + 1.0 * k;
224                 oskar_station_set_element_coords(oskar_station_child(st, j), 0,
225                         k, xyz, xyz, &err);
226                 ASSERT_EQ(0, err) << oskar_get_error_string(err);
227             }
228         }
229     }
230 
231     // Save the telescope model.
232     oskar_telescope_save(telescope, tm, &err);
233     ASSERT_EQ(0, err) << oskar_get_error_string(err);
234 
235     // Load it back again.
236     oskar_Telescope* telescope2 = oskar_telescope_create(OSKAR_SINGLE,
237                     OSKAR_CPU, 0, &err);
238     oskar_telescope_set_position(telescope2,
239             longitude_rad, latitude_rad, altitude_m);
240     oskar_telescope_set_enable_numerical_patterns(telescope, 0);
241     ASSERT_EQ(0, err) << oskar_get_error_string(err);
242     oskar_telescope_load(telescope2, tm, NULL, &err);
243     ASSERT_EQ(0, err) << oskar_get_error_string(err);
244 
245     // Check the contents.
246     ASSERT_EQ(oskar_telescope_num_stations(telescope),
247             oskar_telescope_num_stations(telescope2));
248     EXPECT_NEAR(oskar_telescope_lon_rad(telescope),
249             oskar_telescope_lon_rad(telescope2), 1e-10);
250     EXPECT_NEAR(oskar_telescope_lat_rad(telescope),
251             oskar_telescope_lat_rad(telescope2), 1e-10);
252     EXPECT_NEAR(oskar_telescope_alt_metres(telescope),
253             oskar_telescope_alt_metres(telescope2), 1e-10);
254 
255     double max_ = 0.0, avg_ = 0.0;
256     for (int dim = 0; dim < 3; dim++)
257     {
258         oskar_mem_evaluate_relative_error(
259                 oskar_telescope_station_true_enu_metres(telescope, dim),
260                 oskar_telescope_station_true_enu_metres(telescope2, dim),
261                 0, &max_, &avg_, 0, &err);
262         ASSERT_EQ(0, err) << oskar_get_error_string(err);
263         EXPECT_LT(max_, 1e-5);
264         EXPECT_LT(avg_, 1e-5);
265     }
266     for (int i = 0; i < num_stations; ++i)
267     {
268         oskar_Station* s1 = oskar_telescope_station(telescope, i);
269         oskar_Station* s2 = oskar_telescope_station(telescope2, i);
270         for (int dim = 0; dim < 3; dim++)
271         {
272             oskar_mem_evaluate_relative_error(
273                     oskar_station_element_measured_enu_metres(s1, 0, dim),
274                     oskar_station_element_measured_enu_metres(s2, 0, dim),
275                     0, &max_, &avg_, 0, &err);
276             ASSERT_EQ(0, err) << oskar_get_error_string(err);
277             EXPECT_LT(max_, 1e-5);
278             EXPECT_LT(avg_, 1e-5);
279         }
280 
281         for (int j = 0; j < num_tiles; ++j)
282         {
283             oskar_Station *c1 = oskar_station_child(s1, j);
284             oskar_Station *c2 = oskar_station_child(s2, j);
285             for (int dim = 0; dim < 3; dim++)
286             {
287                 oskar_mem_evaluate_relative_error(
288                         oskar_station_element_measured_enu_metres(c1, 0, dim),
289                         oskar_station_element_measured_enu_metres(c2, 0, dim),
290                         0, &max_, &avg_, 0, &err);
291                 ASSERT_EQ(0, err) << oskar_get_error_string(err);
292                 EXPECT_LT(max_, 1e-5);
293                 EXPECT_LT(avg_, 1e-5);
294             }
295         }
296     }
297 
298     // Copy the telescope model to a new structure.
299     oskar_Telescope* telescope3 = oskar_telescope_create_copy(telescope2,
300             OSKAR_CPU, &err);
301 
302     // Check the contents.
303     ASSERT_EQ(oskar_telescope_num_stations(telescope),
304             oskar_telescope_num_stations(telescope3));
305     EXPECT_NEAR(oskar_telescope_lon_rad(telescope),
306             oskar_telescope_lon_rad(telescope3), 1e-10);
307     EXPECT_NEAR(oskar_telescope_lat_rad(telescope),
308             oskar_telescope_lat_rad(telescope3), 1e-10);
309     EXPECT_NEAR(oskar_telescope_alt_metres(telescope),
310             oskar_telescope_alt_metres(telescope3), 1e-10);
311 
312     for (int dim = 0; dim < 3; dim++)
313     {
314         oskar_mem_evaluate_relative_error(
315                 oskar_telescope_station_true_enu_metres(telescope, dim),
316                 oskar_telescope_station_true_enu_metres(telescope3, dim),
317                 0, &max_, &avg_, 0, &err);
318         ASSERT_EQ(0, err) << oskar_get_error_string(err);
319         EXPECT_LT(max_, 1e-5);
320         EXPECT_LT(avg_, 1e-5);
321     }
322     for (int i = 0; i < num_stations; ++i)
323     {
324         oskar_Station* s1 = oskar_telescope_station(telescope, i);
325         oskar_Station* s3 = oskar_telescope_station(telescope3, i);
326         for (int dim = 0; dim < 3; dim++)
327         {
328             oskar_mem_evaluate_relative_error(
329                     oskar_station_element_measured_enu_metres(s1, 0, dim),
330                     oskar_station_element_measured_enu_metres(s3, 0, dim),
331                     0, &max_, &avg_, 0, &err);
332             ASSERT_EQ(0, err) << oskar_get_error_string(err);
333             EXPECT_LT(max_, 1e-5);
334             EXPECT_LT(avg_, 1e-5);
335         }
336 
337         for (int j = 0; j < num_tiles; ++j)
338         {
339             oskar_Station *c1 = oskar_station_child(s1, j);
340             oskar_Station *c3 = oskar_station_child(s3, j);
341             for (int dim = 0; dim < 3; dim++)
342             {
343                 oskar_mem_evaluate_relative_error(
344                         oskar_station_element_measured_enu_metres(c1, 0, dim),
345                         oskar_station_element_measured_enu_metres(c3, 0, dim),
346                         0, &max_, &avg_, 0, &err);
347                 ASSERT_EQ(0, err) << oskar_get_error_string(err);
348                 EXPECT_LT(max_, 1e-5);
349                 EXPECT_LT(avg_, 1e-5);
350             }
351         }
352     }
353 
354     // Free models.
355     oskar_telescope_free(telescope, &err);
356     oskar_telescope_free(telescope2, &err);
357     oskar_telescope_free(telescope3, &err);
358 
359     // Remove test directory.
360     oskar_dir_remove(tm);
361 }
362 
363 //
364 // TODO: check combinations of telescope model loading and overrides...
365 //
366 
TEST(telescope_model_load_save,test_load_telescope_noise_rms)367 TEST(telescope_model_load_save, test_load_telescope_noise_rms)
368 {
369     // Test cases that should be considered.
370     // -- stddev file at various depths
371     // -- number of noise values vs number of stddev
372     // -- different modes of getting stddev and freq data.
373 
374     const char* root = "temp_test_noise_rms";
375     int err = 0;
376     int num_stations = 2;
377     int num_values = 5;
378     int type = OSKAR_DOUBLE;
379 
380     // Generate the telescope model.
381     vector<double> stddev(num_values), freq_values(num_values);
382     for (int i = 0; i < num_values; ++i)
383     {
384         stddev[i] = i * 0.25 + 0.5;
385         freq_values[i] = 20.0e6 + i * 10.0e6;
386     }
387     generate_noisy_telescope(root, num_stations, freq_values, stddev);
388 
389     // Load it back again.
390     oskar_Telescope* telescope = oskar_telescope_create(type,
391             OSKAR_CPU, 0, &err);
392     oskar_telescope_set_enable_numerical_patterns(telescope, 0);
393     oskar_telescope_set_enable_noise(telescope, true, 1);
394     oskar_telescope_load(telescope, root, NULL, &err);
395     ASSERT_EQ(0, err) << oskar_get_error_string(err);
396     ASSERT_EQ(num_stations, oskar_telescope_num_stations(telescope));
397 
398     // Check the loaded std.dev. values
399     for (int i = 0; i < oskar_telescope_num_station_models(telescope); ++i)
400     {
401         oskar_Station* s = oskar_telescope_station(telescope, i);
402         oskar_Mem *freq = 0, *rms = 0;
403         freq = oskar_station_noise_freq_hz(s);
404         rms = oskar_station_noise_rms_jy(s);
405         ASSERT_EQ(num_values, (int)oskar_mem_length(rms));
406         ASSERT_EQ(num_values, (int)oskar_mem_length(freq));
407         if (type == OSKAR_DOUBLE)
408         {
409             double* r = oskar_mem_double(rms, &err);
410             double* f = oskar_mem_double(freq, &err);
411             for (int j = 0; j < num_values; ++j)
412             {
413                 EXPECT_NEAR(stddev[j], r[j], 1.0e-6);
414                 EXPECT_NEAR(freq_values[j], f[j], 1.0e-6);
415             }
416         }
417         else
418         {
419             float* r = oskar_mem_float(rms, &err);
420             float* f = oskar_mem_float(freq, &err);
421             for (int j = 0; j < num_values; ++j)
422             {
423                 EXPECT_NEAR(stddev[j], r[j], 1.0e-5);
424                 EXPECT_NEAR(freq_values[j], f[j], 1.0e-5);
425             }
426         }
427     }
428 
429     oskar_dir_remove(root);
430     oskar_telescope_free(telescope, &err);
431 }
432 
433 
generate_noisy_telescope(const char * dir,int num_stations,const vector<double> & freqs,const vector<double> & noise)434 static void generate_noisy_telescope(const char* dir, int num_stations,
435         const vector<double>& freqs, const vector<double>& noise)
436 {
437     FILE* f = 0;
438     char* path = 0;
439 
440     // Create a telescope model directory.
441     if (oskar_dir_exists(dir)) oskar_dir_remove(dir);
442     oskar_dir_mkpath(dir);
443 
444     // Write position file.
445     path = oskar_dir_get_path(dir, "position.txt");
446     f = fopen(path, "w");
447     fprintf(f, "0,0\n");
448     fclose(f);
449     free(path);
450 
451     // Write the layout file.
452     path = oskar_dir_get_path(dir, "layout.txt");
453     f = fopen(path, "w");
454     for (int i = 0; i < num_stations; ++i) fprintf(f, "0,0\n");
455     fclose(f);
456     free(path);
457 
458     // Write frequency file.
459     if (!freqs.empty())
460     {
461         path = oskar_dir_get_path(dir, "noise_frequencies.txt");
462         f = fopen(path, "w");
463         for (size_t i = 0; i < freqs.size(); ++i)
464         {
465             fprintf(f, "%.10f\n", freqs[i]);
466         }
467         fclose(f);
468         free(path);
469     }
470 
471     // Write RMS noise values.
472     if (!noise.empty())
473     {
474         path = oskar_dir_get_path(dir, "rms.txt");
475         f = fopen(path, "w");
476         for (size_t i = 0; i < noise.size(); ++i)
477         {
478             fprintf(f, "%.10f\n", noise[i]);
479         }
480         fclose(f);
481         free(path);
482     }
483 }
484