]> Pileus Git - aweather/blobdiff - src/plugins/level2.c
Increase volume size and fix colormap bugs
[aweather] / src / plugins / level2.c
index 6e0a0666457a16ad63f7d801247d8df0340846b7..14f9e34bdfd260b8a0941523343b969d52b6c3d4 100644 (file)
 #include <math.h>
 #include <GL/gl.h>
 #include <glib/gstdio.h>
-#include <gis.h>
+#include <grits.h>
 #include <rsl.h>
 
 #include "level2.h"
 
+#define ISO_MIN 30
+#define ISO_MAX 80
 
 /**************************
  * Data loading functions *
@@ -46,18 +48,22 @@ static void _bscan_sweep(Sweep *sweep, AWeatherColormap *colormap,
        for (int ri = 0; ri < sweep->h.nrays; ri++) {
                Ray *ray  = sweep->ray[ri];
                for (int bi = 0; bi < ray->h.nbins; bi++) {
-                       /* copy RGBA into buffer */
-                       //guint val   = dz_f(ray->range[bi]);
-                       guint8 val   = (guint8)ray->h.f(ray->range[bi]);
                        guint  buf_i = (ri*max_bins+bi)*4;
-                       buf[buf_i+0] = colormap->data[val][0];
-                       buf[buf_i+1] = colormap->data[val][1];
-                       buf[buf_i+2] = colormap->data[val][2];
-                       buf[buf_i+3] = colormap->data[val][3]*0.75; // TESTING
-                       if (val == BADVAL     || val == RFVAL      || val == APFLAG ||
-                           val == NOTFOUND_H || val == NOTFOUND_V || val == NOECHO) {
+                       float  value = ray->h.f(ray->range[bi]);
+
+                       /* Check for bad values */
+                       if (value == BADVAL     || value == RFVAL      || value == APFLAG ||
+                           value == NOTFOUND_H || value == NOTFOUND_V || value == NOECHO) {
                                buf[buf_i+3] = 0x00; // transparent
+                               continue;
                        }
+
+                       /* Copy color to buffer */
+                       guint8 *data = colormap_get(colormap, value);
+                       buf[buf_i+0] = data[0];
+                       buf[buf_i+1] = data[1];
+                       buf[buf_i+2] = data[2];
+                       buf[buf_i+3] = data[3]*0.75; // TESTING
                }
        }
 
@@ -126,11 +132,78 @@ static gboolean _decompress_radar(const gchar *file, const gchar *raw)
        return TRUE;
 }
 
+/* Load the radar into a Grits Volume */
+static void _cart_to_sphere(VolCoord *out, VolCoord *in)
+{
+       gdouble angle = in->x;
+       gdouble dist  = in->y;
+       gdouble tilt  = in->z;
+       gdouble lx    = sin(angle);
+       gdouble ly    = cos(angle);
+       gdouble lz    = sin(tilt);
+       //out->x = (ly*dist)/20000;
+       //out->y = (lz*dist)/10000-0.5;
+       //out->z = (lx*dist)/20000-1.5;
+       out->x = (lx*dist);
+       out->y = (ly*dist);
+       out->z = (lz*dist);
+}
+
+static VolGrid *_load_grid(Volume *vol)
+{
+       g_debug("AWeatherLevel2: _load_grid");
+
+       Sweep *sweep   = vol->sweep[0];
+       Ray   *ray     = sweep->ray[0];
+       gint nsweeps   = vol->h.nsweeps;
+       gint nrays     = sweep->h.nrays/(1/sweep->h.beam_width)+1;
+       gint nbins     = ray->h.nbins  /(1000/ray->h.gate_size);
+       nbins = MIN(nbins, 150);
+
+       VolGrid  *grid = vol_grid_new(nrays, nbins, nsweeps);
+
+       gint rs, bs, val;
+       gint si=0, ri=0, bi=0;
+       for (si = 0; si < nsweeps; si++) {
+               sweep = vol->sweep[si];
+               rs    = 1.0/sweep->h.beam_width;
+       for (ri = 0; ri < nrays; ri++) {
+               /* TODO: missing rays, pick ri based on azmith */
+               ray   = sweep->ray[(ri*rs) % sweep->h.nrays];
+               bs    = 1000/ray->h.gate_size;
+       for (bi = 0; bi < nbins; bi++) {
+               if (bi*bs >= ray->h.nbins)
+                       break;
+               val   = ray->h.f(ray->range[bi*bs]);
+               if (val == BADVAL     || val == RFVAL      ||
+                   val == APFLAG     || val == NOECHO     ||
+                   val == NOTFOUND_H || val == NOTFOUND_V ||
+                   val > 80)
+                       val = 0;
+               VolPoint *point = vol_grid_get(grid, ri, bi, si);
+               point->value = val;
+               point->c.x = deg2rad(ray->h.azimuth);
+               point->c.y = bi*bs*ray->h.gate_size + ray->h.range_bin1;
+               point->c.z = deg2rad(ray->h.elev);
+       } } }
+
+       for (si = 0; si < nsweeps; si++)
+       for (ri = 0; ri < nrays; ri++)
+       for (bi = 0; bi < nbins; bi++) {
+               VolPoint *point = vol_grid_get(grid, ri, bi, si);
+               if (point->c.y == 0)
+                       point->value = nan("");
+               else
+                       _cart_to_sphere(&point->c, &point->c);
+       }
+       return grid;
+}
+
 
 /*********************
  * Drawing functions *
  *********************/
-static void _draw_radar(GisCallback *_self, gpointer _viewer)
+void aweather_level2_draw(GritsObject *_self, GritsOpenGL *opengl)
 {
        AWeatherLevel2 *self = AWEATHER_LEVEL2(_self);
        if (!self->sweep || !self->sweep_tex)
@@ -138,12 +211,6 @@ static void _draw_radar(GisCallback *_self, gpointer _viewer)
 
        /* Draw wsr88d */
        Sweep *sweep = self->sweep;
-       Radar_header *h = &self->radar->h;
-       gdouble lat  = (double)h->latd + (double)h->latm/60 + (double)h->lats/(60*60);
-       gdouble lon  = (double)h->lond + (double)h->lonm/60 + (double)h->lons/(60*60);
-       gdouble elev = h->height;
-       gis_viewer_center_position(self->viewer, lat, lon, elev);
-
        glDisable(GL_ALPHA_TEST);
        glDisable(GL_CULL_FACE);
        glDisable(GL_LIGHTING);
@@ -172,8 +239,8 @@ static void _draw_radar(GisCallback *_self, gpointer _viewer)
                double lx = sin(angle);
                double ly = cos(angle);
 
-               double near_dist = ray->h.range_bin1;
-               double far_dist  = ray->h.nbins*ray->h.gate_size + ray->h.range_bin1;
+               double near_dist = ray->h.range_bin1 - ((double)ray->h.gate_size/2.);
+               double far_dist  = near_dist + (double)ray->h.nbins*ray->h.gate_size;
 
                /* (find middle of bin) / scale for opengl */
                // near left
@@ -207,7 +274,7 @@ static gboolean _set_sweep_cb(gpointer _self)
        g_debug("AWeatherLevel2: _set_sweep_cb");
        AWeatherLevel2 *self = _self;
        _load_sweep_gl(self);
-       gtk_widget_queue_draw(GTK_WIDGET(self->viewer));
+       grits_object_queue_draw(_self);
        g_object_unref(self);
        return FALSE;
 }
@@ -224,31 +291,68 @@ void aweather_level2_set_sweep(AWeatherLevel2 *self,
 
        /* Find colormap */
        self->sweep_colors = NULL;
-       for (int i = 0; self->colormap[i].name; i++)
+       for (int i = 0; self->colormap[i].file; i++)
                if (self->colormap[i].type == type)
                        self->sweep_colors = &self->colormap[i];
-       if (!self->sweep_colors) return;
+       if (!self->sweep_colors) {
+               g_warning("AWeatherLevel2: set_sweep - missing colormap[%d]", type);
+               self->sweep_colors = &self->colormap[0];
+       }
 
        /* Load data */
        g_object_ref(self);
        g_idle_add(_set_sweep_cb, self);
 }
 
-AWeatherLevel2 *aweather_level2_new(GisViewer *viewer,
-               AWeatherColormap *colormap, Radar *radar)
+void aweather_level2_set_iso(AWeatherLevel2 *level2, gfloat level)
+{
+       g_debug("AWeatherLevel2: set_iso - %f", level);
+
+       if (!level2->volume) {
+               g_debug("AWeatherLevel2: set_iso - creating new volume");
+               Volume      *rvol = RSL_get_volume(level2->radar, DZ_INDEX);
+               VolGrid     *grid = _load_grid(rvol);
+               GritsVolume *vol  = grits_volume_new(grid);
+               vol->proj = GRITS_VOLUME_CARTESIAN;
+               vol->disp = GRITS_VOLUME_SURFACE;
+               GRITS_OBJECT(vol)->center = GRITS_OBJECT(level2)->center;
+               grits_viewer_add(GRITS_OBJECT(level2)->viewer,
+                               GRITS_OBJECT(vol), GRITS_LEVEL_WORLD, TRUE);
+               level2->volume = vol;
+       }
+       if (ISO_MIN < level && level < ISO_MAX) {
+               guint8 *data = colormap_get(&level2->colormap[0], level);
+               level2->volume->color[0] = data[0];
+               level2->volume->color[1] = data[1];
+               level2->volume->color[2] = data[2];
+               level2->volume->color[3] = data[3];
+               grits_volume_set_level(level2->volume, level);
+               GRITS_OBJECT(level2->volume)->hidden = FALSE;
+       } else {
+               GRITS_OBJECT(level2->volume)->hidden = TRUE;
+       }
+}
+
+AWeatherLevel2 *aweather_level2_new(Radar *radar, AWeatherColormap *colormap)
 {
        g_debug("AWeatherLevel2: new - %s", radar->h.radar_name);
+       RSL_sort_radar(radar);
        AWeatherLevel2 *self = g_object_new(AWEATHER_TYPE_LEVEL2, NULL);
-       self->viewer   = viewer;
        self->radar    = radar;
        self->colormap = colormap;
        aweather_level2_set_sweep(self, DZ_INDEX, 0);
+
+       GritsPoint center;
+       Radar_header *h = &radar->h;
+       center.lat  = (double)h->latd + (double)h->latm/60 + (double)h->lats/(60*60);
+       center.lon  = (double)h->lond + (double)h->lonm/60 + (double)h->lons/(60*60);
+       center.elev = h->height;
+       GRITS_OBJECT(self)->center = center;
        return self;
 }
 
-AWeatherLevel2 *aweather_level2_new_from_file(GisViewer *viewer,
-               AWeatherColormap *colormap,
-               const gchar *file, const gchar *site)
+AWeatherLevel2 *aweather_level2_new_from_file(const gchar *file, const gchar *site,
+               AWeatherColormap *colormap)
 {
        g_debug("AWeatherLevel2: new_from_file %s %s", site, file);
 
@@ -275,16 +379,25 @@ AWeatherLevel2 *aweather_level2_new_from_file(GisViewer *viewer,
        if (!radar)
                return NULL;
 
-       return aweather_level2_new(viewer, colormaps, radar);
+       return aweather_level2_new(radar, colormaps);
 }
 
 static void _on_sweep_clicked(GtkRadioButton *button, gpointer _level2)
 {
        AWeatherLevel2 *level2 = _level2;
-       gint type = (gint)g_object_get_data(G_OBJECT(button), "type");
-       gint elev = (gint)g_object_get_data(G_OBJECT(button), "elev");
-       aweather_level2_set_sweep(level2, type, (float)elev/100);
-       //self->colormap = level2->sweep_colors;
+       if (gtk_toggle_button_get_active(GTK_TOGGLE_BUTTON(button))) {
+               gint type = (gint)g_object_get_data(G_OBJECT(button), "type");
+               gint elev = (gint)g_object_get_data(G_OBJECT(button), "elev");
+               aweather_level2_set_sweep(level2, type, (float)elev/100);
+               //self->colormap = level2->sweep_colors;
+       }
+}
+
+static void _on_iso_changed(GtkRange *range, gpointer _level2)
+{
+       AWeatherLevel2 *level2 = _level2;
+       gfloat level = gtk_range_get_value(range);
+       aweather_level2_set_iso(level2, level);
 }
 
 GtkWidget *aweather_level2_get_config(AWeatherLevel2 *level2)
@@ -292,7 +405,7 @@ GtkWidget *aweather_level2_get_config(AWeatherLevel2 *level2)
        Radar *radar = level2->radar;
        g_debug("AWeatherLevel2: get_config - %p, %p", level2, radar);
        /* Clear existing items */
-       gdouble elev;
+       gfloat elev;
        guint rows = 1, cols = 1, cur_cols;
        gchar row_label_str[64], col_label_str[64], button_str[64];
        GtkWidget *row_label, *col_label, *button = NULL, *elev_box = NULL;
@@ -308,6 +421,7 @@ GtkWidget *aweather_level2_get_config(AWeatherLevel2 *level2)
                        0,1, 0,1, GTK_FILL,GTK_FILL, 5,0);
        g_free(date_str);
 
+       /* Add sweeps */
        for (guint vi = 0; vi < radar->h.nvolumes; vi++) {
                Volume *vol = radar->v[vi];
                if (vol == NULL) continue;
@@ -361,22 +475,34 @@ GtkWidget *aweather_level2_get_config(AWeatherLevel2 *level2)
                        g_signal_connect(button, "clicked", G_CALLBACK(_on_sweep_clicked), level2);
                }
        }
+
+       /* Add Iso-surface volume */
+       g_object_get(table, "n-columns", &cols, NULL);
+       row_label = gtk_label_new("<b>Isosurface:</b>");
+       gtk_label_set_use_markup(GTK_LABEL(row_label), TRUE);
+       gtk_misc_set_alignment(GTK_MISC(row_label), 1, 0.5);
+       gtk_table_attach(GTK_TABLE(table), row_label,
+                       0,1, rows,rows+1, GTK_FILL,GTK_FILL, 5,0);
+       GtkWidget *scale = gtk_hscale_new_with_range(ISO_MIN, ISO_MAX, 0.5);
+       gtk_widget_set_size_request(scale, -1, 26);
+       gtk_scale_set_value_pos(GTK_SCALE(scale), GTK_POS_LEFT);
+       gtk_range_set_inverted(GTK_RANGE(scale), TRUE);
+       gtk_range_set_value(GTK_RANGE(scale), ISO_MAX);
+       g_signal_connect(scale, "value-changed", G_CALLBACK(_on_iso_changed), level2);
+       gtk_table_attach(GTK_TABLE(table), scale,
+                       1,cols+1, rows,rows+1, GTK_FILL|GTK_EXPAND,GTK_FILL, 0,0);
+       /* Shove all the buttons to the left, but keep the slider expanded */
+       gtk_table_attach(GTK_TABLE(table), gtk_label_new(""),
+                       cols,cols+1, 0,1, GTK_FILL|GTK_EXPAND,GTK_FILL, 0,0);
        return table;
 }
 
 /****************
  * GObject code *
  ****************/
-G_DEFINE_TYPE(AWeatherLevel2, aweather_level2, GIS_TYPE_CALLBACK);
+G_DEFINE_TYPE(AWeatherLevel2, aweather_level2, GRITS_TYPE_OBJECT);
 static void aweather_level2_init(AWeatherLevel2 *self)
 {
-       GIS_CALLBACK(self)->callback  = _draw_radar;
-       GIS_CALLBACK(self)->user_data = self;
-}
-static void aweather_level2_dispose(GObject *_self)
-{
-       g_debug("AWeatherLevel2: dispose - %p", _self);
-       G_OBJECT_CLASS(aweather_level2_parent_class)->dispose(_self);
 }
 static void aweather_level2_finalize(GObject *_self)
 {
@@ -390,5 +516,5 @@ static void aweather_level2_finalize(GObject *_self)
 static void aweather_level2_class_init(AWeatherLevel2Class *klass)
 {
        G_OBJECT_CLASS(klass)->finalize = aweather_level2_finalize;
-       G_OBJECT_CLASS(klass)->dispose  = aweather_level2_dispose;
+       GRITS_OBJECT_CLASS(klass)->draw   = aweather_level2_draw;
 }