[grits] Add volume testing code

[Andy Spencer <andys src gnome org>]

commit 8f40fd04fdb9316e53cf0b765003cffba04af481
Author: Andy Spencer <andy753421 gmail com>
Date:   Mon Jan 24 03:36:33 2011 +0000

    Add volume testing code
    
    Render some iso-balls with optional radar rendering
    
    Based on GritsTester

 examples/.gitignore      |    1 +
 examples/volume/mkfile   |   25 ++++++
 examples/volume/volume.c |  193 ++++++++++++++++++++++++++++++++++++++++++++++
 3 files changed, 219 insertions(+), 0 deletions(-)
---
diff --git a/examples/.gitignore b/examples/.gitignore
index 0f06fc9..a8ee815 100644
--- a/examples/.gitignore
+++ b/examples/.gitignore
@@ -9,3 +9,4 @@ sort/sort
 tess/convex
 tess/split
 tex/tex
+volume/volume
diff --git a/examples/volume/mkfile b/examples/volume/mkfile
new file mode 100644
index 0000000..d89b78b
--- /dev/null
+++ b/examples/volume/mkfile
@@ -0,0 +1,25 @@
+MKSHELL=/usr/lib/plan9/bin/rc
+
+PKG_CONFIG_PATH=../../src/
+LD_LIBRARY_PATH=../../src/.libs/
+PKGS=grits
+
+CFLAGS=-Wall -Wno-unused -Werror -g -p -pg --std=gnu99 -I../
+PROGS=volume
+default:V: volume-run
+
+volume_libs=-lrsl
+volume: volume.o ../tester.o
+
+memcheck: all
+	G_SLICE=always-malloc                               \
+	G_DEBUG=gc-friendly,resident-modules                \
+	valgrind '--track-origins=yes'                      \
+	         '--leak-check=full'                        \
+	         '--leak-resolution=high'                   \
+	         '--num-callers=50'                         \
+	         '--suppressions=../../src/gtk.suppression' \
+	         ./volume                                   \
+	>[2] valgrind.out
+
+<$HOME/lib/mkcommon
diff --git a/examples/volume/volume.c b/examples/volume/volume.c
new file mode 100644
index 0000000..a3f421f
--- /dev/null
+++ b/examples/volume/volume.c
@@ -0,0 +1,193 @@
+/*
+ * Copyright (C) 2009-2010 Andy Spencer <andy753421 gmail com>
+ *
+ * This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program.  If not, see <http://www.gnu.org/licenses/>.
+ */
+
+#include <math.h>
+#include <glib.h>
+#include <gdk/gdkgl.h>
+#include <gdk/gdkkeysyms.h>
+
+#include <objects/grits-volume.h>
+#include <objects/grits-callback.h>
+#include <objects/marching.h>
+
+#include <GL/gl.h>
+#include <rsl.h>
+#include <tester.h>
+
+/******************
+ * iso ball setup *
+ ******************/
+static double dist  = 0.75;
+
+static gdouble distp(VolPoint *a,
+		gdouble bx, gdouble by, gdouble bz)
+{
+	return 1/((a->c.x-bx)*(a->c.x-bx) +
+	          (a->c.y-by)*(a->c.y-by) +
+	          (a->c.z-bz)*(a->c.z-bz)) * 0.10;
+	//return 1-MIN(1,sqrt((a->c.x-bx)*(a->c.x-bx) +
+	//                    (a->c.y-by)*(a->c.y-by) +
+	//                    (a->c.z-bz)*(a->c.z-bz)));
+}
+
+static VolGrid *load_balls(float dist, int xs, int ys, int zs)
+{
+	VolGrid *grid = vol_grid_new(xs, ys, zs);
+	for (int x = 0; x < xs; x++)
+	for (int y = 0; y < ys; y++)
+	for (int z = 0; z < zs; z++) {
+		VolPoint *point = vol_grid_get(grid, x, y, z);
+		point->c.x = ((double)x/(xs-1)*2-1);
+		point->c.y = ((double)y/(ys-1)*2-1);
+		point->c.z = ((double)z/(zs-1)*2-1);
+		point->value =
+			distp(point, -dist, 0, 0) +
+			distp(point,  dist, 0, 0);
+		point->value *= 100;
+	}
+	return grid;
+}
+
+/***************
+ * radar setup *
+ ***************/
+/* 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;
+}
+
+static VolGrid *load_radar(gchar *file, gchar *site)
+{
+	/* Load radar file */
+	RSL_read_these_sweeps("all", NULL);
+	Radar  *rad = RSL_wsr88d_to_radar(file, site);
+	Volume *vol = RSL_get_volume(rad, DZ_INDEX);
+	RSL_sort_rays_in_volume(vol);
+
+	/* Count dimensions */
+	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, 100);
+
+	/* Convert to VolGrid */
+	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 azimuth */
+		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);
+	} } }
+
+	/* Convert to spherical coords */
+	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;
+}
+
+/**********
+ * Common *
+ **********/
+static gboolean key_press(GritsTester *tester, GdkEventKey *event, GritsVolume *volume)
+{
+	     if (event->keyval == GDK_v) grits_volume_set_level(volume, volume->level-0.5);
+	else if (event->keyval == GDK_V) grits_volume_set_level(volume, volume->level+0.5);
+	else if (event->keyval == GDK_d) dist  += 0.5;
+	else if (event->keyval == GDK_D) dist  -= 0.5;
+	return FALSE;
+}
+
+/********
+ * Main *
+ ********/
+int main(int argc, char **argv)
+{
+	gtk_init(&argc, &argv);
+	GtkWidget   *window = gtk_window_new(GTK_WINDOW_TOPLEVEL);
+	GritsTester *tester = grits_tester_new();
+	gtk_container_add(GTK_CONTAINER(window), GTK_WIDGET(tester));
+	gtk_widget_show_all(window);
+
+	/* Grits Volume */
+	VolGrid *balls_grid = load_balls(dist, 50, 50, 50);
+	GritsVolume *balls = grits_volume_new(balls_grid);
+	balls->proj = GRITS_VOLUME_CARTESIAN;
+	balls->disp = GRITS_VOLUME_SURFACE;
+	balls->color[0] = (0.8)*0xff;
+	balls->color[1] = (0.6)*0xff;
+	balls->color[2] = (0.2)*0xff;
+	balls->color[3] = (1.0)*0xff;
+	grits_volume_set_level(balls, 50);
+	grits_tester_add(tester, GRITS_OBJECT(balls));
+	g_signal_connect(tester, "key-press-event", G_CALLBACK(key_press), balls);
+
+	/* Grits Volume */
+	//char *file = "/home/andy/.cache/grits/nexrad/level2/KGWX/KGWX_20101130_0459.raw";
+	//char *file = "/home/andy/.cache/grits/nexrad/level2/KTLX/KTLX_19990503_2351.raw";
+	//VolGrid *radar_grid = load_radar(file, "KTLX");
+	//GritsVolume *radar = grits_volume_new(radar_grid);
+	//radar->proj = GRITS_VOLUME_SPHERICAL;
+	//radar->disp = GRITS_VOLUME_SURFACE;
+	//radar->color[0] = (0.8)*0xff;
+	//radar->color[1] = (0.6)*0xff;
+	//radar->color[2] = (0.2)*0xff;
+	//radar->color[3] = (1.0)*0xff;
+	//grits_volume_set_level(radar, 50);
+	//grits_tester_add(tester, GRITS_OBJECT(radar));
+	//g_signal_connect(tester, "key-press-event", G_CALLBACK(key_press), radar);
+
+	/* Go */
+	gtk_main();
+	return 0;
+}