2 * Copyright (C) 2009-2011 Andy Spencer <andy753421@gmail.com>
4 * This program is free software: you can redistribute it and/or modify
5 * it under the terms of the GNU General Public License as published by
6 * the Free Software Foundation, either version 3 of the License, or
7 * (at your option) any later version.
9 * This program is distributed in the hope that it will be useful,
10 * but WITHOUT ANY WARRANTY; without even the implied warranty of
11 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12 * GNU General Public License for more details.
14 * You should have received a copy of the GNU General Public License
15 * along with this program. If not, see <http://www.gnu.org/licenses/>.
20 #include <gdk/gdkgl.h>
21 #include <gdk/gdkkeysyms.h>
23 #include <objects/grits-volume.h>
24 #include <objects/grits-callback.h>
25 #include <objects/marching.h>
34 static double dist = 0.75;
36 static gdouble distp(VolPoint *a,
37 gdouble bx, gdouble by, gdouble bz)
39 return 1/((a->c.x-bx)*(a->c.x-bx) +
40 (a->c.y-by)*(a->c.y-by) +
41 (a->c.z-bz)*(a->c.z-bz)) * 0.10;
42 //return 1-MIN(1,sqrt((a->c.x-bx)*(a->c.x-bx) +
43 // (a->c.y-by)*(a->c.y-by) +
44 // (a->c.z-bz)*(a->c.z-bz)));
47 static VolGrid *load_balls(float dist, int xs, int ys, int zs)
49 VolGrid *grid = vol_grid_new(xs, ys, zs);
50 for (int x = 0; x < xs; x++)
51 for (int y = 0; y < ys; y++)
52 for (int z = 0; z < zs; z++) {
53 VolPoint *point = vol_grid_get(grid, x, y, z);
54 point->c.x = ((double)x/(xs-1)*2-1);
55 point->c.y = ((double)y/(ys-1)*2-1);
56 point->c.z = ((double)z/(zs-1)*2-1);
58 distp(point, -dist, 0, 0) +
59 distp(point, dist, 0, 0);
68 /* Load the radar into a Grits Volume */
69 static void _cart_to_sphere(VolCoord *out, VolCoord *in)
71 gdouble angle = in->x;
74 gdouble lx = sin(angle);
75 gdouble ly = cos(angle);
76 gdouble lz = sin(tilt);
77 out->x = (ly*dist)/20000;
78 out->y = (lz*dist)/10000-0.5;
79 out->z = (lx*dist)/20000-1.5;
82 static VolGrid *load_radar(gchar *file, gchar *site)
85 RSL_read_these_sweeps("all", NULL);
86 Radar *rad = RSL_wsr88d_to_radar(file, site);
87 Volume *vol = RSL_get_volume(rad, DZ_INDEX);
88 RSL_sort_rays_in_volume(vol);
90 /* Count dimensions */
91 Sweep *sweep = vol->sweep[0];
92 Ray *ray = sweep->ray[0];
93 gint nsweeps = vol->h.nsweeps;
94 gint nrays = sweep->h.nrays/(1/sweep->h.beam_width)+1;
95 gint nbins = ray->h.nbins /(1000/ray->h.gate_size);
96 nbins = MIN(nbins, 100);
98 /* Convert to VolGrid */
99 VolGrid *grid = vol_grid_new(nrays, nbins, nsweeps);
102 gint si=0, ri=0, bi=0;
103 for (si = 0; si < nsweeps; si++) {
104 sweep = vol->sweep[si];
105 rs = 1.0/sweep->h.beam_width;
106 for (ri = 0; ri < nrays; ri++) {
107 /* TODO: missing rays, pick ri based on azimuth */
108 ray = sweep->ray[(ri*rs) % sweep->h.nrays];
109 bs = 1000/ray->h.gate_size;
110 for (bi = 0; bi < nbins; bi++) {
111 if (bi*bs >= ray->h.nbins)
113 val = ray->h.f(ray->range[bi*bs]);
114 if (val == BADVAL || val == RFVAL ||
115 val == APFLAG || val == NOECHO ||
116 val == NOTFOUND_H || val == NOTFOUND_V ||
119 VolPoint *point = vol_grid_get(grid, ri, bi, si);
121 point->c.x = deg2rad(ray->h.azimuth);
122 point->c.y = bi*bs*ray->h.gate_size + ray->h.range_bin1;
123 point->c.z = deg2rad(ray->h.elev);
126 /* Convert to spherical coords */
127 for (si = 0; si < nsweeps; si++)
128 for (ri = 0; ri < nrays; ri++)
129 for (bi = 0; bi < nbins; bi++) {
130 VolPoint *point = vol_grid_get(grid, ri, bi, si);
132 point->value = nan("");
134 _cart_to_sphere(&point->c, &point->c);
142 static gboolean key_press(GritsTester *tester, GdkEventKey *event, GritsVolume *volume)
144 if (event->keyval == GDK_v) grits_volume_set_level(volume, volume->level-0.5);
145 else if (event->keyval == GDK_V) grits_volume_set_level(volume, volume->level+0.5);
146 else if (event->keyval == GDK_d) dist += 0.5;
147 else if (event->keyval == GDK_D) dist -= 0.5;
154 int main(int argc, char **argv)
156 gtk_init(&argc, &argv);
157 GtkWidget *window = gtk_window_new(GTK_WINDOW_TOPLEVEL);
158 GritsTester *tester = grits_tester_new();
159 gtk_container_add(GTK_CONTAINER(window), GTK_WIDGET(tester));
160 gtk_widget_show_all(window);
163 VolGrid *balls_grid = load_balls(dist, 50, 50, 50);
164 GritsVolume *balls = grits_volume_new(balls_grid);
165 balls->proj = GRITS_VOLUME_CARTESIAN;
166 balls->disp = GRITS_VOLUME_SURFACE;
167 balls->color[0] = (0.8)*0xff;
168 balls->color[1] = (0.6)*0xff;
169 balls->color[2] = (0.2)*0xff;
170 balls->color[3] = (1.0)*0xff;
171 grits_volume_set_level(balls, 50);
172 grits_tester_add(tester, GRITS_OBJECT(balls));
173 g_signal_connect(tester, "key-press-event", G_CALLBACK(key_press), balls);
176 //char *file = "/home/andy/.cache/grits/nexrad/level2/KGWX/KGWX_20101130_0459.raw";
177 //char *file = "/home/andy/.cache/grits/nexrad/level2/KTLX/KTLX_19990503_2351.raw";
178 //VolGrid *radar_grid = load_radar(file, "KTLX");
179 //GritsVolume *radar = grits_volume_new(radar_grid);
180 //radar->proj = GRITS_VOLUME_SPHERICAL;
181 //radar->disp = GRITS_VOLUME_SURFACE;
182 //radar->color[0] = (0.8)*0xff;
183 //radar->color[1] = (0.6)*0xff;
184 //radar->color[2] = (0.2)*0xff;
185 //radar->color[3] = (1.0)*0xff;
186 //grits_volume_set_level(radar, 50);
187 //grits_tester_add(tester, GRITS_OBJECT(radar));
188 //g_signal_connect(tester, "key-press-event", G_CALLBACK(key_press), radar);