Adding some (commented out) support for generating iso surfaces.

Thanks to Jamie Zawinski and Paul Bourke for the code for this.

diff --git a/src/Makefile.am b/src/Makefile.am
index f1f9329ca7bfd50d72be9df9fdbb95e22e228cbd..620cdcb18f280288915bb1a3f8238f3c0343030d 100644 (file)
--- a/src/Makefile.am
+++ b/src/Makefile.am
@@ -10,6 +10,7 @@ aweather_SOURCES  = main.c \
        aweather-plugin.c aweather-plugin.h \
        data.c            data.h            \
        location.c        location.h        \
+       marching.c        marching.h        \
        plugin-radar.c    plugin-radar.h    \
        plugin-radar-colormap.c             \
        plugin-ridge.c    plugin-ridge.h    \
@@ -30,7 +31,7 @@ MAINTAINERCLEANFILES = Makefile.in
        glib-genmarshal --prefix=aweather_cclosure_marshal --header $< > $@
 
 test: all
-       ./aweather -o -d 7 -s KLSX
+       ./aweather -o -d 7 -s KVNX
 
 gdb: all
        gdb ./aweather

diff --git a/src/aweather-gui.c b/src/aweather-gui.c
index 07fb089b69cdf58b8e35756ae64597f2b151b105..da4c161d066bd6d67831baad3085505d16bbc5fd 100644 (file)
--- a/src/aweather-gui.c
+++ b/src/aweather-gui.c
@@ -23,6 +23,7 @@
 #include <GL/glu.h>
 #include <math.h>
 
+#include "misc.h"
 #include "aweather-gui.h"
 #include "aweather-view.h"
 #include "aweather-plugin.h"
@@ -97,18 +98,22 @@ gboolean on_drawing_button_press(GtkWidget *widget, GdkEventButton *event, AWeat
 }
 gboolean on_drawing_key_press(GtkWidget *widget, GdkEventKey *event, AWeatherGui *gui)
 {
-       g_debug("AWeatherGui: on_drawing_key_press - key=%x, state=%x",
-                       event->keyval, event->state);
+       g_debug("AWeatherGui: on_drawing_key_press - key=%x, state=%x, plus=%x",
+                       event->keyval, event->state, GDK_plus);
        AWeatherView *view = aweather_gui_get_view(gui);
        double x,y,z;
        aweather_view_get_location(view, &x, &y, &z);
        guint kv = event->keyval;
-       if      (kv == GDK_Right || kv == GDK_l) aweather_view_pan(view,  z/10, 0, 0);
-       else if (kv == GDK_Left  || kv == GDK_h) aweather_view_pan(view, -z/10, 0, 0);
-       else if (kv == GDK_Up    || kv == GDK_k) aweather_view_pan(view, 0,  z/10, 0);
+       if      (kv == GDK_Left  || kv == GDK_h) aweather_view_pan(view, -z/10, 0, 0);
        else if (kv == GDK_Down  || kv == GDK_j) aweather_view_pan(view, 0, -z/10, 0);
+       else if (kv == GDK_Up    || kv == GDK_k) aweather_view_pan(view, 0,  z/10, 0);
+       else if (kv == GDK_Right || kv == GDK_l) aweather_view_pan(view,  z/10, 0, 0);
        else if (kv == GDK_minus || kv == GDK_o) aweather_view_zoom(view, 10./9);
        else if (kv == GDK_plus  || kv == GDK_i) aweather_view_zoom(view, 9./10);
+       else if (kv == GDK_H                   ) aweather_view_rotate(view,  0, -10, 0);
+       else if (kv == GDK_J                   ) aweather_view_rotate(view,  10,  0, 0);
+       else if (kv == GDK_K                   ) aweather_view_rotate(view, -10,  0, 0);
+       else if (kv == GDK_L                   ) aweather_view_rotate(view,  0,  10, 0);
        return TRUE;
 }
 
@@ -246,11 +251,10 @@ gboolean on_configure(GtkWidget *da, GdkEventConfigure *event, AWeatherGui *gui)
        /* Perspective */
        glMatrixMode(GL_PROJECTION);
        glLoadIdentity();
-       double rad  = atan((height/2)/500);
-       double deg  = (rad*180)/M_PI;
+       double ang = atan((height/2)/500);
 
-       gluPerspective(deg*2, width/height, -z-20, -z+20);
-       //gluPerspective(deg*2, width/height, 1, 500*1000);
+       //gluPerspective(r2d(ang)*2, width/height, -z-20, -z+20);
+       gluPerspective(r2d(ang)*2, width/height, 1, 500*1000);
 
        aweather_gui_gl_end(gui);
        return FALSE;
@@ -261,15 +265,18 @@ gboolean on_expose(GtkWidget *da, GdkEventExpose *event, AWeatherGui *gui)
        g_debug("AWeatherGui: on_expose - begin");
        aweather_gui_gl_begin(gui);
 
-       double x, y, z;
+       double lx, ly, lz;
+       double rx, ry, rz;
        AWeatherView *view = aweather_gui_get_view(gui);
-       aweather_view_get_location(view, &x, &y, &z);
+       aweather_view_get_location(view, &lx, &ly, &lz);
+       aweather_view_get_rotation(view, &rx, &ry, &rz);
 
        glMatrixMode(GL_MODELVIEW);
        glLoadIdentity();
-       glTranslatef(x, y, z);
-
-       //glRotatef(-45, 1, 0, 0);
+       glTranslatef(lx, ly, lz);
+       glRotatef(rx, 1, 0, 0);
+       glRotatef(ry, 0, 1, 0);
+       glRotatef(rz, 0, 0, 1);
 
        /* Expose plugins */
        for (GList *cur = gui->plugins; cur; cur = cur->next) {
@@ -283,7 +290,7 @@ gboolean on_expose(GtkWidget *da, GdkEventExpose *event, AWeatherGui *gui)
        return FALSE;
 }
 
-void on_location_changed(AWeatherView *view,
+void on_state_changed(AWeatherView *view,
                gdouble x, gdouble y, gdouble z, AWeatherGui *gui)
 {
        /* Reset clipping area and redraw */
@@ -436,7 +443,9 @@ AWeatherGui *aweather_gui_new()
        g_signal_connect(self, "key-press-event",
                        G_CALLBACK(on_gui_key_press), self);
        g_signal_connect(self->view, "location-changed",
-                       G_CALLBACK(on_location_changed), self);
+                       G_CALLBACK(on_state_changed), self);
+       g_signal_connect(self->view, "rotation-changed",
+                       G_CALLBACK(on_state_changed), self);
        g_signal_connect_swapped(self->view, "offline",
                        G_CALLBACK(gtk_toggle_action_set_active),
                        aweather_gui_get_object(self, "offline"));

diff --git a/src/aweather-view.c b/src/aweather-view.c
index fe8bafa462f42384e94ca949b29c80cd7c4991d8..0c765fb32f3e48d50c59d1e7a77073a53eb68c43 100644 (file)
--- a/src/aweather-view.c
+++ b/src/aweather-view.c
@@ -33,6 +33,7 @@ enum {
        SIG_TIME_CHANGED,
        SIG_SITE_CHANGED,
        SIG_LOCATION_CHANGED,
+       SIG_ROTATION_CHANGED,
        SIG_REFRESH,
        SIG_OFFLINE,
        NUM_SIGNALS,
@@ -50,6 +51,9 @@ static void aweather_view_init(AWeatherView *self)
        self->location[0] = 0;
        self->location[1] = 0;
        self->location[2] = -300*1000;
+       self->rotation[0] = 0;
+       self->rotation[1] = 0;
+       self->rotation[2] = 0;
        self->offline = FALSE;
 }
 static void aweather_view_dispose(GObject *gobject)
@@ -143,6 +147,19 @@ static void aweather_view_class_init(AWeatherViewClass *klass)
                        G_TYPE_DOUBLE,
                        G_TYPE_DOUBLE,
                        G_TYPE_DOUBLE);
+       signals[SIG_ROTATION_CHANGED] = g_signal_new(
+                       "rotation-changed",
+                       G_TYPE_FROM_CLASS(gobject_class),
+                       G_SIGNAL_RUN_LAST,
+                       0,
+                       NULL,
+                       NULL,
+                       aweather_cclosure_marshal_VOID__DOUBLE_DOUBLE_DOUBLE,
+                       G_TYPE_NONE,
+                       3,
+                       G_TYPE_DOUBLE,
+                       G_TYPE_DOUBLE,
+                       G_TYPE_DOUBLE);
        signals[SIG_REFRESH] = g_signal_new(
                        "refresh",
                        G_TYPE_FROM_CLASS(gobject_class),
@@ -174,6 +191,13 @@ static void _aweather_view_emit_location_changed(AWeatherView *view)
                        view->location[1],
                        view->location[2]);
 }
+static void _aweather_view_emit_rotation_changed(AWeatherView *view)
+{
+       g_signal_emit(view, signals[SIG_ROTATION_CHANGED], 0, 
+                       view->rotation[0],
+                       view->rotation[1],
+                       view->rotation[2]);
+}
 static void _aweather_view_emit_time_changed(AWeatherView *view)
 {
        g_signal_emit(view, signals[SIG_TIME_CHANGED], 0,
@@ -257,6 +281,35 @@ void aweather_view_zoom(AWeatherView *view, gdouble scale)
        _aweather_view_emit_location_changed(view);
 }
 
+void aweather_view_set_rotation(AWeatherView *view, gdouble x, gdouble y, gdouble z)
+{
+       g_assert(AWEATHER_IS_VIEW(view));
+       g_debug("AWeatherView: set_rotation");
+       view->rotation[0] = x;
+       view->rotation[1] = y;
+       view->rotation[2] = z;
+       _aweather_view_emit_rotation_changed(view);
+}
+
+void aweather_view_get_rotation(AWeatherView *view, gdouble *x, gdouble *y, gdouble *z)
+{
+       g_assert(AWEATHER_IS_VIEW(view));
+       g_debug("AWeatherView: get_rotation");
+       *x = view->rotation[0];
+       *y = view->rotation[1];
+       *z = view->rotation[2];
+}
+
+void aweather_view_rotate(AWeatherView *view, gdouble x, gdouble y, gdouble z)
+{
+       g_assert(AWEATHER_IS_VIEW(view));
+       g_debug("AWeatherView: rotate - x=%.0f, y=%.0f, z=%.0f", x, y, z);
+       view->rotation[0] += x;
+       view->rotation[1] += y;
+       view->rotation[2] += z;
+       _aweather_view_emit_rotation_changed(view);
+}
+
 void aweather_view_refresh(AWeatherView *view)
 {
        g_debug("AWeatherView: refresh");

diff --git a/src/aweather-view.h b/src/aweather-view.h
index 7fc964b8d1f28be78da0a654db0a2990d1d91adc..c67c8d82b3bac9b14b2fe23872ab585a4da9336e 100644 (file)
--- a/src/aweather-view.h
+++ b/src/aweather-view.h
@@ -38,6 +38,7 @@ struct _AWeatherView {
        gchar   *time;
        gchar   *site;
        gdouble  location[3];
+       gdouble  rotation[3];
        gboolean offline;
 };
 
@@ -60,6 +61,10 @@ void aweather_view_get_location(AWeatherView *view, gdouble *x, gdouble *y, gdou
 void aweather_view_pan         (AWeatherView *view, gdouble  x, gdouble  y, gdouble  z);
 void aweather_view_zoom        (AWeatherView *view, gdouble  scale);
 
+void aweather_view_set_rotation(AWeatherView *view, gdouble  x, gdouble  y, gdouble  z);
+void aweather_view_get_rotation(AWeatherView *view, gdouble *x, gdouble *y, gdouble *z);
+void aweather_view_rotate      (AWeatherView *view, gdouble  x, gdouble  y, gdouble  z);
+
 void aweather_view_refresh(AWeatherView *view);
 
 void aweather_view_set_offline(AWeatherView *view, gboolean offline);

diff --git a/src/marching.c b/src/marching.c
new file mode 100644 (file)
index 0000000..7c88311
--- /dev/null
+++ b/src/marching.c
@@ -0,0 +1,645 @@
+/**
+ * Copyright (c) 2002 Jamie Zawinski <jwz@jwz.org>
+ * Copyright (c) 2009 Andy Spencer <andy753421@gmail.com>
+ *
+ * Utility functions to create "marching cubes" meshes from 3d fields.
+ *
+ * Permission to use, copy, modify, distribute, and sell this software and its
+ * documentation for any purpose is hereby granted without fee, provided that
+ * the above copyright notice appear in all copies and that both that
+ * copyright notice and this permission notice appear in supporting
+ * documentation.  No representations are made about the suitability of this
+ * software for any purpose.  It is provided "as is" without express or 
+ * implied warranty.
+ *
+ * Marching cubes implementation by Paul Bourke <pbourke@swin.edu.au>
+ * http://astronomy.swin.edu.au/~pbourke/modelling/polygonise/
+ */
+
+#include <math.h>
+#include <glib.h>
+
+#include <GL/gl.h>
+
+#include "marching.h"
+
+/* Calculate the unit normal at p given two other points p1,p2 on the
+ * surface. The normal points in the direction of p1 crossproduct p2
+ */
+XYZ calc_normal (XYZ p, XYZ p1, XYZ p2)
+{
+       XYZ n, pa, pb;
+       pa.x = p1.x - p.x;
+       pa.y = p1.y - p.y;
+       pa.z = p1.z - p.z;
+       pb.x = p2.x - p.x;
+       pb.y = p2.y - p.y;
+       pb.z = p2.z - p.z;
+       n.x = pa.y * pb.z - pa.z * pb.y;
+       n.y = pa.z * pb.x - pa.x * pb.z;
+       n.z = pa.x * pb.y - pa.y * pb.x;
+       return (n);
+}
+
+/* Call glNormal3f() with a normal of the indicated triangle.
+ */
+void do_normal(GLfloat x1, GLfloat y1, GLfloat z1,
+               GLfloat x2, GLfloat y2, GLfloat z2,
+               GLfloat x3, GLfloat y3, GLfloat z3)
+{
+       XYZ p1, p2, p3, p;
+       p1.x = x1; p1.y = y1; p1.z = z1;
+       p2.x = x2; p2.y = y2; p2.z = z2;
+       p3.x = x3; p3.y = y3; p3.z = z3;
+       p = calc_normal (p1, p2, p3);
+       glNormal3f (p.x, p.y, p.z);
+}
+
+/* Indexing convention:
+ * 
+ *           Vertices:                    Edges:
+ * 
+ *        4  ______________ 5           ______________
+ *         /|            /|           /|     4      /|
+ *        / |         6 / |       7  / |8        5 / |
+ *    7  /_____________/  |        /______________/  | 9
+ *      |   |         |   |        |   |   6     |   |
+ *      | 0 |_________|___| 1      |   |_________|10_|
+ *      |  /          |  /      11 | 3/     0    |  /
+ *      | /           | /          | /           | / 1
+ *    3 |/____________|/ 2         |/____________|/
+ *                                        2
+ */
+
+static const int edgeTable[256] = {
+       0x0  , 0x109, 0x203, 0x30a, 0x406, 0x50f, 0x605, 0x70c,
+       0x80c, 0x905, 0xa0f, 0xb06, 0xc0a, 0xd03, 0xe09, 0xf00,
+       0x190, 0x99 , 0x393, 0x29a, 0x596, 0x49f, 0x795, 0x69c,
+       0x99c, 0x895, 0xb9f, 0xa96, 0xd9a, 0xc93, 0xf99, 0xe90,
+       0x230, 0x339, 0x33 , 0x13a, 0x636, 0x73f, 0x435, 0x53c,
+       0xa3c, 0xb35, 0x83f, 0x936, 0xe3a, 0xf33, 0xc39, 0xd30,
+       0x3a0, 0x2a9, 0x1a3, 0xaa , 0x7a6, 0x6af, 0x5a5, 0x4ac,
+       0xbac, 0xaa5, 0x9af, 0x8a6, 0xfaa, 0xea3, 0xda9, 0xca0,
+       0x460, 0x569, 0x663, 0x76a, 0x66 , 0x16f, 0x265, 0x36c,
+       0xc6c, 0xd65, 0xe6f, 0xf66, 0x86a, 0x963, 0xa69, 0xb60,
+       0x5f0, 0x4f9, 0x7f3, 0x6fa, 0x1f6, 0xff , 0x3f5, 0x2fc,
+       0xdfc, 0xcf5, 0xfff, 0xef6, 0x9fa, 0x8f3, 0xbf9, 0xaf0,
+       0x650, 0x759, 0x453, 0x55a, 0x256, 0x35f, 0x55 , 0x15c,
+       0xe5c, 0xf55, 0xc5f, 0xd56, 0xa5a, 0xb53, 0x859, 0x950,
+       0x7c0, 0x6c9, 0x5c3, 0x4ca, 0x3c6, 0x2cf, 0x1c5, 0xcc ,
+       0xfcc, 0xec5, 0xdcf, 0xcc6, 0xbca, 0xac3, 0x9c9, 0x8c0,
+       0x8c0, 0x9c9, 0xac3, 0xbca, 0xcc6, 0xdcf, 0xec5, 0xfcc,
+       0xcc , 0x1c5, 0x2cf, 0x3c6, 0x4ca, 0x5c3, 0x6c9, 0x7c0,
+       0x950, 0x859, 0xb53, 0xa5a, 0xd56, 0xc5f, 0xf55, 0xe5c,
+       0x15c, 0x55 , 0x35f, 0x256, 0x55a, 0x453, 0x759, 0x650,
+       0xaf0, 0xbf9, 0x8f3, 0x9fa, 0xef6, 0xfff, 0xcf5, 0xdfc,
+       0x2fc, 0x3f5, 0xff , 0x1f6, 0x6fa, 0x7f3, 0x4f9, 0x5f0,
+       0xb60, 0xa69, 0x963, 0x86a, 0xf66, 0xe6f, 0xd65, 0xc6c,
+       0x36c, 0x265, 0x16f, 0x66 , 0x76a, 0x663, 0x569, 0x460,
+       0xca0, 0xda9, 0xea3, 0xfaa, 0x8a6, 0x9af, 0xaa5, 0xbac,
+       0x4ac, 0x5a5, 0x6af, 0x7a6, 0xaa , 0x1a3, 0x2a9, 0x3a0,
+       0xd30, 0xc39, 0xf33, 0xe3a, 0x936, 0x83f, 0xb35, 0xa3c,
+       0x53c, 0x435, 0x73f, 0x636, 0x13a, 0x33 , 0x339, 0x230,
+       0xe90, 0xf99, 0xc93, 0xd9a, 0xa96, 0xb9f, 0x895, 0x99c,
+       0x69c, 0x795, 0x49f, 0x596, 0x29a, 0x393, 0x99 , 0x190,
+       0xf00, 0xe09, 0xd03, 0xc0a, 0xb06, 0xa0f, 0x905, 0x80c,
+       0x70c, 0x605, 0x50f, 0x406, 0x30a, 0x203, 0x109, 0x0
+};
+
+static const int triTable[256][16] = {
+       {-1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 0,  8,  3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 0,  1,  9, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 1,  8,  3,  9,  8,  1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 1,  2, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 0,  8,  3,  1,  2, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 9,  2, 10,  0,  2,  9, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 2,  8,  3,  2, 10,  8, 10,  9,  8, -1, -1, -1, -1, -1, -1, -1},
+       { 3, 11,  2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 0, 11,  2,  8, 11,  0, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 1,  9,  0,  2,  3, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 1, 11,  2,  1,  9, 11,  9,  8, 11, -1, -1, -1, -1, -1, -1, -1},
+       { 3, 10,  1, 11, 10,  3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 0, 10,  1,  0,  8, 10,  8, 11, 10, -1, -1, -1, -1, -1, -1, -1},
+       { 3,  9,  0,  3, 11,  9, 11, 10,  9, -1, -1, -1, -1, -1, -1, -1},
+       { 9,  8, 10, 10,  8, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 4,  7,  8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 4,  3,  0,  7,  3,  4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 0,  1,  9,  8,  4,  7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 4,  1,  9,  4,  7,  1,  7,  3,  1, -1, -1, -1, -1, -1, -1, -1},
+       { 1,  2, 10,  8,  4,  7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 3,  4,  7,  3,  0,  4,  1,  2, 10, -1, -1, -1, -1, -1, -1, -1},
+       { 9,  2, 10,  9,  0,  2,  8,  4,  7, -1, -1, -1, -1, -1, -1, -1},
+       { 2, 10,  9,  2,  9,  7,  2,  7,  3,  7,  9,  4, -1, -1, -1, -1},
+       { 8,  4,  7,  3, 11,  2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       {11,  4,  7, 11,  2,  4,  2,  0,  4, -1, -1, -1, -1, -1, -1, -1},
+       { 9,  0,  1,  8,  4,  7,  2,  3, 11, -1, -1, -1, -1, -1, -1, -1},
+       { 4,  7, 11,  9,  4, 11,  9, 11,  2,  9,  2,  1, -1, -1, -1, -1},
+       { 3, 10,  1,  3, 11, 10,  7,  8,  4, -1, -1, -1, -1, -1, -1, -1},
+       { 1, 11, 10,  1,  4, 11,  1,  0,  4,  7, 11,  4, -1, -1, -1, -1},
+       { 4,  7,  8,  9,  0, 11,  9, 11, 10, 11,  0,  3, -1, -1, -1, -1},
+       { 4,  7, 11,  4, 11,  9,  9, 11, 10, -1, -1, -1, -1, -1, -1, -1},
+       { 9,  5,  4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 9,  5,  4,  0,  8,  3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 0,  5,  4,  1,  5,  0, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 8,  5,  4,  8,  3,  5,  3,  1,  5, -1, -1, -1, -1, -1, -1, -1},
+       { 1,  2, 10,  9,  5,  4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 3,  0,  8,  1,  2, 10,  4,  9,  5, -1, -1, -1, -1, -1, -1, -1},
+       { 5,  2, 10,  5,  4,  2,  4,  0,  2, -1, -1, -1, -1, -1, -1, -1},
+       { 2, 10,  5,  3,  2,  5,  3,  5,  4,  3,  4,  8, -1, -1, -1, -1},
+       { 9,  5,  4,  2,  3, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 0, 11,  2,  0,  8, 11,  4,  9,  5, -1, -1, -1, -1, -1, -1, -1},
+       { 0,  5,  4,  0,  1,  5,  2,  3, 11, -1, -1, -1, -1, -1, -1, -1},
+       { 2,  1,  5,  2,  5,  8,  2,  8, 11,  4,  8,  5, -1, -1, -1, -1},
+       {10,  3, 11, 10,  1,  3,  9,  5,  4, -1, -1, -1, -1, -1, -1, -1},
+       { 4,  9,  5,  0,  8,  1,  8, 10,  1,  8, 11, 10, -1, -1, -1, -1},
+       { 5,  4,  0,  5,  0, 11,  5, 11, 10, 11,  0,  3, -1, -1, -1, -1},
+       { 5,  4,  8,  5,  8, 10, 10,  8, 11, -1, -1, -1, -1, -1, -1, -1},
+       { 9,  7,  8,  5,  7,  9, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 9,  3,  0,  9,  5,  3,  5,  7,  3, -1, -1, -1, -1, -1, -1, -1},
+       { 0,  7,  8,  0,  1,  7,  1,  5,  7, -1, -1, -1, -1, -1, -1, -1},
+       { 1,  5,  3,  3,  5,  7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 9,  7,  8,  9,  5,  7, 10,  1,  2, -1, -1, -1, -1, -1, -1, -1},
+       {10,  1,  2,  9,  5,  0,  5,  3,  0,  5,  7,  3, -1, -1, -1, -1},
+       { 8,  0,  2,  8,  2,  5,  8,  5,  7, 10,  5,  2, -1, -1, -1, -1},
+       { 2, 10,  5,  2,  5,  3,  3,  5,  7, -1, -1, -1, -1, -1, -1, -1},
+       { 7,  9,  5,  7,  8,  9,  3, 11,  2, -1, -1, -1, -1, -1, -1, -1},
+       { 9,  5,  7,  9,  7,  2,  9,  2,  0,  2,  7, 11, -1, -1, -1, -1},
+       { 2,  3, 11,  0,  1,  8,  1,  7,  8,  1,  5,  7, -1, -1, -1, -1},
+       {11,  2,  1, 11,  1,  7,  7,  1,  5, -1, -1, -1, -1, -1, -1, -1},
+       { 9,  5,  8,  8,  5,  7, 10,  1,  3, 10,  3, 11, -1, -1, -1, -1},
+       { 5,  7,  0,  5,  0,  9,  7, 11,  0,  1,  0, 10, 11, 10,  0, -1},
+       {11, 10,  0, 11,  0,  3, 10,  5,  0,  8,  0,  7,  5,  7,  0, -1},
+       {11, 10,  5,  7, 11,  5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       {10,  6,  5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 0,  8,  3,  5, 10,  6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 9,  0,  1,  5, 10,  6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 1,  8,  3,  1,  9,  8,  5, 10,  6, -1, -1, -1, -1, -1, -1, -1},
+       { 1,  6,  5,  2,  6,  1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 1,  6,  5,  1,  2,  6,  3,  0,  8, -1, -1, -1, -1, -1, -1, -1},
+       { 9,  6,  5,  9,  0,  6,  0,  2,  6, -1, -1, -1, -1, -1, -1, -1},
+       { 5,  9,  8,  5,  8,  2,  5,  2,  6,  3,  2,  8, -1, -1, -1, -1},
+       { 2,  3, 11, 10,  6,  5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       {11,  0,  8, 11,  2,  0, 10,  6,  5, -1, -1, -1, -1, -1, -1, -1},
+       { 0,  1,  9,  2,  3, 11,  5, 10,  6, -1, -1, -1, -1, -1, -1, -1},
+       { 5, 10,  6,  1,  9,  2,  9, 11,  2,  9,  8, 11, -1, -1, -1, -1},
+       { 6,  3, 11,  6,  5,  3,  5,  1,  3, -1, -1, -1, -1, -1, -1, -1},
+       { 0,  8, 11,  0, 11,  5,  0,  5,  1,  5, 11,  6, -1, -1, -1, -1},
+       { 3, 11,  6,  0,  3,  6,  0,  6,  5,  0,  5,  9, -1, -1, -1, -1},
+       { 6,  5,  9,  6,  9, 11, 11,  9,  8, -1, -1, -1, -1, -1, -1, -1},
+       { 5, 10,  6,  4,  7,  8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 4,  3,  0,  4,  7,  3,  6,  5, 10, -1, -1, -1, -1, -1, -1, -1},
+       { 1,  9,  0,  5, 10,  6,  8,  4,  7, -1, -1, -1, -1, -1, -1, -1},
+       {10,  6,  5,  1,  9,  7,  1,  7,  3,  7,  9,  4, -1, -1, -1, -1},
+       { 6,  1,  2,  6,  5,  1,  4,  7,  8, -1, -1, -1, -1, -1, -1, -1},
+       { 1,  2,  5,  5,  2,  6,  3,  0,  4,  3,  4,  7, -1, -1, -1, -1},
+       { 8,  4,  7,  9,  0,  5,  0,  6,  5,  0,  2,  6, -1, -1, -1, -1},
+       { 7,  3,  9,  7,  9,  4,  3,  2,  9,  5,  9,  6,  2,  6,  9, -1},
+       { 3, 11,  2,  7,  8,  4, 10,  6,  5, -1, -1, -1, -1, -1, -1, -1},
+       { 5, 10,  6,  4,  7,  2,  4,  2,  0,  2,  7, 11, -1, -1, -1, -1},
+       { 0,  1,  9,  4,  7,  8,  2,  3, 11,  5, 10,  6, -1, -1, -1, -1},
+       { 9,  2,  1,  9, 11,  2,  9,  4, 11,  7, 11,  4,  5, 10,  6, -1},
+       { 8,  4,  7,  3, 11,  5,  3,  5,  1,  5, 11,  6, -1, -1, -1, -1},
+       { 5,  1, 11,  5, 11,  6,  1,  0, 11,  7, 11,  4,  0,  4, 11, -1},
+       { 0,  5,  9,  0,  6,  5,  0,  3,  6, 11,  6,  3,  8,  4,  7, -1},
+       { 6,  5,  9,  6,  9, 11,  4,  7,  9,  7, 11,  9, -1, -1, -1, -1},
+       {10,  4,  9,  6,  4, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 4, 10,  6,  4,  9, 10,  0,  8,  3, -1, -1, -1, -1, -1, -1, -1},
+       {10,  0,  1, 10,  6,  0,  6,  4,  0, -1, -1, -1, -1, -1, -1, -1},
+       { 8,  3,  1,  8,  1,  6,  8,  6,  4,  6,  1, 10, -1, -1, -1, -1},
+       { 1,  4,  9,  1,  2,  4,  2,  6,  4, -1, -1, -1, -1, -1, -1, -1},
+       { 3,  0,  8,  1,  2,  9,  2,  4,  9,  2,  6,  4, -1, -1, -1, -1},
+       { 0,  2,  4,  4,  2,  6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 8,  3,  2,  8,  2,  4,  4,  2,  6, -1, -1, -1, -1, -1, -1, -1},
+       {10,  4,  9, 10,  6,  4, 11,  2,  3, -1, -1, -1, -1, -1, -1, -1},
+       { 0,  8,  2,  2,  8, 11,  4,  9, 10,  4, 10,  6, -1, -1, -1, -1},
+       { 3, 11,  2,  0,  1,  6,  0,  6,  4,  6,  1, 10, -1, -1, -1, -1},
+       { 6,  4,  1,  6,  1, 10,  4,  8,  1,  2,  1, 11,  8, 11,  1, -1},
+       { 9,  6,  4,  9,  3,  6,  9,  1,  3, 11,  6,  3, -1, -1, -1, -1},
+       { 8, 11,  1,  8,  1,  0, 11,  6,  1,  9,  1,  4,  6,  4,  1, -1},
+       { 3, 11,  6,  3,  6,  0,  0,  6,  4, -1, -1, -1, -1, -1, -1, -1},
+       { 6,  4,  8, 11,  6,  8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 7, 10,  6,  7,  8, 10,  8,  9, 10, -1, -1, -1, -1, -1, -1, -1},
+       { 0,  7,  3,  0, 10,  7,  0,  9, 10,  6,  7, 10, -1, -1, -1, -1},
+       {10,  6,  7,  1, 10,  7,  1,  7,  8,  1,  8,  0, -1, -1, -1, -1},
+       {10,  6,  7, 10,  7,  1,  1,  7,  3, -1, -1, -1, -1, -1, -1, -1},
+       { 1,  2,  6,  1,  6,  8,  1,  8,  9,  8,  6,  7, -1, -1, -1, -1},
+       { 2,  6,  9,  2,  9,  1,  6,  7,  9,  0,  9,  3,  7,  3,  9, -1},
+       { 7,  8,  0,  7,  0,  6,  6,  0,  2, -1, -1, -1, -1, -1, -1, -1},
+       { 7,  3,  2,  6,  7,  2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 2,  3, 11, 10,  6,  8, 10,  8,  9,  8,  6,  7, -1, -1, -1, -1},
+       { 2,  0,  7,  2,  7, 11,  0,  9,  7,  6,  7, 10,  9, 10,  7, -1},
+       { 1,  8,  0,  1,  7,  8,  1, 10,  7,  6,  7, 10,  2,  3, 11, -1},
+       {11,  2,  1, 11,  1,  7, 10,  6,  1,  6,  7,  1, -1, -1, -1, -1},
+       { 8,  9,  6,  8,  6,  7,  9,  1,  6, 11,  6,  3,  1,  3,  6, -1},
+       { 0,  9,  1, 11,  6,  7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 7,  8,  0,  7,  0,  6,  3, 11,  0, 11,  6,  0, -1, -1, -1, -1},
+       { 7, 11,  6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 7,  6, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 3,  0,  8, 11,  7,  6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 0,  1,  9, 11,  7,  6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 8,  1,  9,  8,  3,  1, 11,  7,  6, -1, -1, -1, -1, -1, -1, -1},
+       {10,  1,  2,  6, 11,  7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 1,  2, 10,  3,  0,  8,  6, 11,  7, -1, -1, -1, -1, -1, -1, -1},
+       { 2,  9,  0,  2, 10,  9,  6, 11,  7, -1, -1, -1, -1, -1, -1, -1},
+       { 6, 11,  7,  2, 10,  3, 10,  8,  3, 10,  9,  8, -1, -1, -1, -1},
+       { 7,  2,  3,  6,  2,  7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 7,  0,  8,  7,  6,  0,  6,  2,  0, -1, -1, -1, -1, -1, -1, -1},
+       { 2,  7,  6,  2,  3,  7,  0,  1,  9, -1, -1, -1, -1, -1, -1, -1},
+       { 1,  6,  2,  1,  8,  6,  1,  9,  8,  8,  7,  6, -1, -1, -1, -1},
+       {10,  7,  6, 10,  1,  7,  1,  3,  7, -1, -1, -1, -1, -1, -1, -1},
+       {10,  7,  6,  1,  7, 10,  1,  8,  7,  1,  0,  8, -1, -1, -1, -1},
+       { 0,  3,  7,  0,  7, 10,  0, 10,  9,  6, 10,  7, -1, -1, -1, -1},
+       { 7,  6, 10,  7, 10,  8,  8, 10,  9, -1, -1, -1, -1, -1, -1, -1},
+       { 6,  8,  4, 11,  8,  6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 3,  6, 11,  3,  0,  6,  0,  4,  6, -1, -1, -1, -1, -1, -1, -1},
+       { 8,  6, 11,  8,  4,  6,  9,  0,  1, -1, -1, -1, -1, -1, -1, -1},
+       { 9,  4,  6,  9,  6,  3,  9,  3,  1, 11,  3,  6, -1, -1, -1, -1},
+       { 6,  8,  4,  6, 11,  8,  2, 10,  1, -1, -1, -1, -1, -1, -1, -1},
+       { 1,  2, 10,  3,  0, 11,  0,  6, 11,  0,  4,  6, -1, -1, -1, -1},
+       { 4, 11,  8,  4,  6, 11,  0,  2,  9,  2, 10,  9, -1, -1, -1, -1},
+       {10,  9,  3, 10,  3,  2,  9,  4,  3, 11,  3,  6,  4,  6,  3, -1},
+       { 8,  2,  3,  8,  4,  2,  4,  6,  2, -1, -1, -1, -1, -1, -1, -1},
+       { 0,  4,  2,  4,  6,  2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 1,  9,  0,  2,  3,  4,  2,  4,  6,  4,  3,  8, -1, -1, -1, -1},
+       { 1,  9,  4,  1,  4,  2,  2,  4,  6, -1, -1, -1, -1, -1, -1, -1},
+       { 8,  1,  3,  8,  6,  1,  8,  4,  6,  6, 10,  1, -1, -1, -1, -1},
+       {10,  1,  0, 10,  0,  6,  6,  0,  4, -1, -1, -1, -1, -1, -1, -1},
+       { 4,  6,  3,  4,  3,  8,  6, 10,  3,  0,  3,  9, 10,  9,  3, -1},
+       {10,  9,  4,  6, 10,  4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 4,  9,  5,  7,  6, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 0,  8,  3,  4,  9,  5, 11,  7,  6, -1, -1, -1, -1, -1, -1, -1},
+       { 5,  0,  1,  5,  4,  0,  7,  6, 11, -1, -1, -1, -1, -1, -1, -1},
+       {11,  7,  6,  8,  3,  4,  3,  5,  4,  3,  1,  5, -1, -1, -1, -1},
+       { 9,  5,  4, 10,  1,  2,  7,  6, 11, -1, -1, -1, -1, -1, -1, -1},
+       { 6, 11,  7,  1,  2, 10,  0,  8,  3,  4,  9,  5, -1, -1, -1, -1},
+       { 7,  6, 11,  5,  4, 10,  4,  2, 10,  4,  0,  2, -1, -1, -1, -1},
+       { 3,  4,  8,  3,  5,  4,  3,  2,  5, 10,  5,  2, 11,  7,  6, -1},
+       { 7,  2,  3,  7,  6,  2,  5,  4,  9, -1, -1, -1, -1, -1, -1, -1},
+       { 9,  5,  4,  0,  8,  6,  0,  6,  2,  6,  8,  7, -1, -1, -1, -1},
+       { 3,  6,  2,  3,  7,  6,  1,  5,  0,  5,  4,  0, -1, -1, -1, -1},
+       { 6,  2,  8,  6,  8,  7,  2,  1,  8,  4,  8,  5,  1,  5,  8, -1},
+       { 9,  5,  4, 10,  1,  6,  1,  7,  6,  1,  3,  7, -1, -1, -1, -1},
+       { 1,  6, 10,  1,  7,  6,  1,  0,  7,  8,  7,  0,  9,  5,  4, -1},
+       { 4,  0, 10,  4, 10,  5,  0,  3, 10,  6, 10,  7,  3,  7, 10, -1},
+       { 7,  6, 10,  7, 10,  8,  5,  4, 10,  4,  8, 10, -1, -1, -1, -1},
+       { 6,  9,  5,  6, 11,  9, 11,  8,  9, -1, -1, -1, -1, -1, -1, -1},
+       { 3,  6, 11,  0,  6,  3,  0,  5,  6,  0,  9,  5, -1, -1, -1, -1},
+       { 0, 11,  8,  0,  5, 11,  0,  1,  5,  5,  6, 11, -1, -1, -1, -1},
+       { 6, 11,  3,  6,  3,  5,  5,  3,  1, -1, -1, -1, -1, -1, -1, -1},
+       { 1,  2, 10,  9,  5, 11,  9, 11,  8, 11,  5,  6, -1, -1, -1, -1},
+       { 0, 11,  3,  0,  6, 11,  0,  9,  6,  5,  6,  9,  1,  2, 10, -1},
+       {11,  8,  5, 11,  5,  6,  8,  0,  5, 10,  5,  2,  0,  2,  5, -1},
+       { 6, 11,  3,  6,  3,  5,  2, 10,  3, 10,  5,  3, -1, -1, -1, -1},
+       { 5,  8,  9,  5,  2,  8,  5,  6,  2,  3,  8,  2, -1, -1, -1, -1},
+       { 9,  5,  6,  9,  6,  0,  0,  6,  2, -1, -1, -1, -1, -1, -1, -1},
+       { 1,  5,  8,  1,  8,  0,  5,  6,  8,  3,  8,  2,  6,  2,  8, -1},
+       { 1,  5,  6,  2,  1,  6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 1,  3,  6,  1,  6, 10,  3,  8,  6,  5,  6,  9,  8,  9,  6, -1},
+       {10,  1,  0, 10,  0,  6,  9,  5,  0,  5,  6,  0, -1, -1, -1, -1},
+       { 0,  3,  8,  5,  6, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       {10,  5,  6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       {11,  5, 10,  7,  5, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       {11,  5, 10, 11,  7,  5,  8,  3,  0, -1, -1, -1, -1, -1, -1, -1},
+       { 5, 11,  7,  5, 10, 11,  1,  9,  0, -1, -1, -1, -1, -1, -1, -1},
+       {10,  7,  5, 10, 11,  7,  9,  8,  1,  8,  3,  1, -1, -1, -1, -1},
+       {11,  1,  2, 11,  7,  1,  7,  5,  1, -1, -1, -1, -1, -1, -1, -1},
+       { 0,  8,  3,  1,  2,  7,  1,  7,  5,  7,  2, 11, -1, -1, -1, -1},
+       { 9,  7,  5,  9,  2,  7,  9,  0,  2,  2, 11,  7, -1, -1, -1, -1},
+       { 7,  5,  2,  7,  2, 11,  5,  9,  2,  3,  2,  8,  9,  8,  2, -1},
+       { 2,  5, 10,  2,  3,  5,  3,  7,  5, -1, -1, -1, -1, -1, -1, -1},
+       { 8,  2,  0,  8,  5,  2,  8,  7,  5, 10,  2,  5, -1, -1, -1, -1},
+       { 9,  0,  1,  5, 10,  3,  5,  3,  7,  3, 10,  2, -1, -1, -1, -1},
+       { 9,  8,  2,  9,  2,  1,  8,  7,  2, 10,  2,  5,  7,  5,  2, -1},
+       { 1,  3,  5,  3,  7,  5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 0,  8,  7,  0,  7,  1,  1,  7,  5, -1, -1, -1, -1, -1, -1, -1},
+       { 9,  0,  3,  9,  3,  5,  5,  3,  7, -1, -1, -1, -1, -1, -1, -1},
+       { 9,  8,  7,  5,  9,  7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 5,  8,  4,  5, 10,  8, 10, 11,  8, -1, -1, -1, -1, -1, -1, -1},
+       { 5,  0,  4,  5, 11,  0,  5, 10, 11, 11,  3,  0, -1, -1, -1, -1},
+       { 0,  1,  9,  8,  4, 10,  8, 10, 11, 10,  4,  5, -1, -1, -1, -1},
+       {10, 11,  4, 10,  4,  5, 11,  3,  4,  9,  4,  1,  3,  1,  4, -1},
+       { 2,  5,  1,  2,  8,  5,  2, 11,  8,  4,  5,  8, -1, -1, -1, -1},
+       { 0,  4, 11,  0, 11,  3,  4,  5, 11,  2, 11,  1,  5,  1, 11, -1},
+       { 0,  2,  5,  0,  5,  9,  2, 11,  5,  4,  5,  8, 11,  8,  5, -1},
+       { 9,  4,  5,  2, 11,  3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 2,  5, 10,  3,  5,  2,  3,  4,  5,  3,  8,  4, -1, -1, -1, -1},
+       { 5, 10,  2,  5,  2,  4,  4,  2,  0, -1, -1, -1, -1, -1, -1, -1},
+       { 3, 10,  2,  3,  5, 10,  3,  8,  5,  4,  5,  8,  0,  1,  9, -1},
+       { 5, 10,  2,  5,  2,  4,  1,  9,  2,  9,  4,  2, -1, -1, -1, -1},
+       { 8,  4,  5,  8,  5,  3,  3,  5,  1, -1, -1, -1, -1, -1, -1, -1},
+       { 0,  4,  5,  1,  0,  5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 8,  4,  5,  8,  5,  3,  9,  0,  5,  0,  3,  5, -1, -1, -1, -1},
+       { 9,  4,  5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 4, 11,  7,  4,  9, 11,  9, 10, 11, -1, -1, -1, -1, -1, -1, -1},
+       { 0,  8,  3,  4,  9,  7,  9, 11,  7,  9, 10, 11, -1, -1, -1, -1},
+       { 1, 10, 11,  1, 11,  4,  1,  4,  0,  7,  4, 11, -1, -1, -1, -1},
+       { 3,  1,  4,  3,  4,  8,  1, 10,  4,  7,  4, 11, 10, 11,  4, -1},
+       { 4, 11,  7,  9, 11,  4,  9,  2, 11,  9,  1,  2, -1, -1, -1, -1},
+       { 9,  7,  4,  9, 11,  7,  9,  1, 11,  2, 11,  1,  0,  8,  3, -1},
+       {11,  7,  4, 11,  4,  2,  2,  4,  0, -1, -1, -1, -1, -1, -1, -1},
+       {11,  7,  4, 11,  4,  2,  8,  3,  4,  3,  2,  4, -1, -1, -1, -1},
+       { 2,  9, 10,  2,  7,  9,  2,  3,  7,  7,  4,  9, -1, -1, -1, -1},
+       { 9, 10,  7,  9,  7,  4, 10,  2,  7,  8,  7,  0,  2,  0,  7, -1},
+       { 3,  7, 10,  3, 10,  2,  7,  4, 10,  1, 10,  0,  4,  0, 10, -1},
+       { 1, 10,  2,  8,  7,  4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 4,  9,  1,  4,  1,  7,  7,  1,  3, -1, -1, -1, -1, -1, -1, -1},
+       { 4,  9,  1,  4,  1,  7,  0,  8,  1,  8,  7,  1, -1, -1, -1, -1},
+       { 4,  0,  3,  7,  4,  3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 4,  8,  7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 9, 10,  8, 10, 11,  8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 3,  0,  9,  3,  9, 11, 11,  9, 10, -1, -1, -1, -1, -1, -1, -1},
+       { 0,  1, 10,  0, 10,  8,  8, 10, 11, -1, -1, -1, -1, -1, -1, -1},
+       { 3,  1, 10, 11,  3, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 1,  2, 11,  1, 11,  9,  9, 11,  8, -1, -1, -1, -1, -1, -1, -1},
+       { 3,  0,  9,  3,  9, 11,  1,  2,  9,  2, 11,  9, -1, -1, -1, -1},
+       { 0,  2, 11,  8,  0, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 3,  2, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 2,  3,  8,  2,  8, 10, 10,  8,  9, -1, -1, -1, -1, -1, -1, -1},
+       { 9, 10,  2,  0,  9,  2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 2,  3,  8,  2,  8, 10,  0,  1,  8,  1, 10,  8, -1, -1, -1, -1},
+       { 1, 10,  2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 1,  3,  8,  9,  1,  8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 0,  9,  1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       { 0,  3,  8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+       {-1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1}
+};
+
+
+
+/* Linearly interpolate the position where an isosurface cuts
+ * an edge between two vertices, each with their own scalar value
+ */
+static XYZ interp_vertex (double isolevel, XYZ p1, XYZ p2, double valp1, double valp2)
+{
+       double mu;
+       XYZ p;
+
+       if (ABS(isolevel-valp1) < 0.00001)
+               return(p1);
+       if (ABS(isolevel-valp2) < 0.00001)
+               return(p2);
+       if (ABS(valp1-valp2) < 0.00001)
+               return(p1);
+       mu = (isolevel - valp1) / (valp2 - valp1);
+       p.x = p1.x + mu * (p2.x - p1.x);
+       p.y = p1.y + mu * (p2.y - p1.y);
+       p.z = p1.z + mu * (p2.z - p1.z);
+
+       return(p);
+}
+
+
+/* Given a grid cell and an isolevel, calculate the triangular
+ * facets required to represent the isosurface through the cell.
+ * Return the number of triangular facets.
+ * `triangles' will be loaded up with the vertices at most 5 triangular facets.
+ * 0 will be returned if the grid cell is either totally above
+ * of totally below the isolevel.
+ *
+ * By Paul Bourke <pbourke@swin.edu.au>
+ */
+int march_one_cube (GRIDCELL grid, double isolevel, TRIANGLE *triangles)
+{
+       int i, ntriang;
+       int cubeindex;
+       XYZ vertlist[12];
+
+       /*
+        * Determine the index into the edge table which
+        * tells us which vertices are inside of the surface
+        */
+       cubeindex = 0;
+       if (grid.val[0] < isolevel) cubeindex |= 1;
+       if (grid.val[1] < isolevel) cubeindex |= 2;
+       if (grid.val[2] < isolevel) cubeindex |= 4;
+       if (grid.val[3] < isolevel) cubeindex |= 8;
+       if (grid.val[4] < isolevel) cubeindex |= 16;
+       if (grid.val[5] < isolevel) cubeindex |= 32;
+       if (grid.val[6] < isolevel) cubeindex |= 64;
+       if (grid.val[7] < isolevel) cubeindex |= 128;
+
+       /* Cube is entirely in/out of the surface */
+       if (edgeTable[cubeindex] == 0)
+               return(0);
+
+       /* Find the vertices where the surface intersects the cube */
+       if (edgeTable[cubeindex] & 1)
+               vertlist[0] =
+                       interp_vertex (isolevel,grid.p[0],grid.p[1],grid.val[0],grid.val[1]);
+       if (edgeTable[cubeindex] & 2)
+               vertlist[1] =
+                       interp_vertex (isolevel,grid.p[1],grid.p[2],grid.val[1],grid.val[2]);
+       if (edgeTable[cubeindex] & 4)
+               vertlist[2] =
+                       interp_vertex (isolevel,grid.p[2],grid.p[3],grid.val[2],grid.val[3]);
+       if (edgeTable[cubeindex] & 8)
+               vertlist[3] =
+                       interp_vertex (isolevel,grid.p[3],grid.p[0],grid.val[3],grid.val[0]);
+       if (edgeTable[cubeindex] & 16)
+               vertlist[4] =
+                       interp_vertex (isolevel,grid.p[4],grid.p[5],grid.val[4],grid.val[5]);
+       if (edgeTable[cubeindex] & 32)
+               vertlist[5] =
+                       interp_vertex (isolevel,grid.p[5],grid.p[6],grid.val[5],grid.val[6]);
+       if (edgeTable[cubeindex] & 64)
+               vertlist[6] =
+                       interp_vertex (isolevel,grid.p[6],grid.p[7],grid.val[6],grid.val[7]);
+       if (edgeTable[cubeindex] & 128)
+               vertlist[7] =
+                       interp_vertex (isolevel,grid.p[7],grid.p[4],grid.val[7],grid.val[4]);
+       if (edgeTable[cubeindex] & 256)
+               vertlist[8] =
+                       interp_vertex (isolevel,grid.p[0],grid.p[4],grid.val[0],grid.val[4]);
+       if (edgeTable[cubeindex] & 512)
+               vertlist[9] =
+                       interp_vertex (isolevel,grid.p[1],grid.p[5],grid.val[1],grid.val[5]);
+       if (edgeTable[cubeindex] & 1024)
+               vertlist[10] =
+                       interp_vertex (isolevel,grid.p[2],grid.p[6],grid.val[2],grid.val[6]);
+       if (edgeTable[cubeindex] & 2048)
+               vertlist[11] =
+                       interp_vertex (isolevel,grid.p[3],grid.p[7],grid.val[3],grid.val[7]);
+
+       /* Create the triangle */
+       ntriang = 0;
+       for (i=0; triTable[cubeindex][i] != -1; i+=3)
+       {
+               triangles[ntriang].p[0] = vertlist[triTable[cubeindex][i  ]];
+               triangles[ntriang].p[1] = vertlist[triTable[cubeindex][i+1]];
+               triangles[ntriang].p[2] = vertlist[triTable[cubeindex][i+2]];
+               ntriang++;
+       }
+
+       return(ntriang);
+}
+
+
+/* Walking the grid.  By jwz.
+ */
+
+
+/* Computes the normal of the scalar field at the given point,
+ * for vertex normals (as opposed to face normals.)
+ */
+void do_function_normal (double x, double y, double z,
+               double (*compute_fn) (double x, double y, double z,
+                       void *closure),
+               void *c)
+{
+       XYZ n;
+       double off = 0.5;
+       n.x = compute_fn (x-off, y, z, c) - compute_fn (x+off, y, z, c);
+       n.y = compute_fn (x, y-off, z, c) - compute_fn (x, y+off, z, c);
+       n.z = compute_fn (x, y, z-off, c) - compute_fn (x, y, z+off, c);
+       /* normalize (&n); */
+       glNormal3f (n.x, n.y, n.z);
+}
+
+
+/* Given a function capable of generating a value at any XYZ position,
+ * creates OpenGL faces for the solids defined.
+ *
+ * init_fn is called at the beginning for initial, and returns an object.
+ * free_fn is called at the end.
+ *
+ * compute_fn is called for each XYZ in the specified grid, and returns
+ * the double value of that coordinate.  If smoothing is on, then
+ * compute_fn will also be called twice more for each emitted vertex,
+ * in order to calculate vertex normals (so don't assume it will only
+ * be called with values falling on the grid boundaries.)
+ *
+ * Points are inside an object if the are less than `isolevel', and
+ * outside otherwise.
+ */
+void
+marching_cubes (int grid_size,     /* density of the mesh */
+               double isolevel,   /* cutoff point for "in" versus "out" */
+               int wireframe_p,   /* wireframe, or solid */
+               int smooth_p,      /* smooth, or faceted */
+
+               void * (*init_fn)    (double grid_size, void *closure1),
+               double (*compute_fn) (double x, double y, double z,
+                       void *closure2),
+               void   (*free_fn)    (void *closure2),
+               void *closure1,
+
+               unsigned long *polygon_count)
+{
+       int planesize = grid_size * grid_size;
+       int x, y, z;
+       void *closure2 = 0;
+       unsigned long polys = 0;
+       double *layers;
+
+       layers = (double *) g_malloc0(sizeof (*layers) * planesize * 2);
+       if (!layers)
+       {
+               g_error("out of memory for %dx%dx%d grid\n",
+                               grid_size, grid_size, 2);
+       }
+
+       if (init_fn)
+               closure2 = init_fn (grid_size, closure1);
+
+       glFrontFace(GL_CCW);
+       if (!wireframe_p)
+               glBegin (GL_TRIANGLES);
+
+       for (z = 0; z < grid_size; z++)
+       {
+               double *layer0 = (z & 1 ? layers+planesize : layers);
+               double *layer1 = (z & 1 ? layers           : layers+planesize);
+               double *row;
+
+               /* Fill in the XY grid on the currently-bottommost layer. */
+               row = layer1;
+               for (y = 0; y < grid_size; y++, row += grid_size)
+               {
+                       double *cell = row;
+                       for (x = 0; x < grid_size; x++, cell++)
+                               *cell = compute_fn (x, y, z, closure2);
+               }
+
+               /* Now we've completed one layer (an XY slice of Z.)  Now we can
+                  generate the polygons that fill the space between this layer
+                  and the previous one (unless this is the first layer.)
+                */
+               if (z == 0) continue;
+
+               for (y = 1; y < grid_size; y += 1)
+                       for (x = 1; x < grid_size; x += 1)
+                       {
+                               TRIANGLE tri[6];
+                               int i, ntri;
+                               GRIDCELL cell;
+
+                               /* This is kinda hokey, there ought to be a more efficient
+                                  way to do this... */
+                               cell.p[0].x = x-1; cell.p[0].y = y-1; cell.p[0].z = z-1;
+                               cell.p[1].x = x  ; cell.p[1].y = y-1; cell.p[1].z = z-1;
+                               cell.p[2].x = x  ; cell.p[2].y = y  ; cell.p[2].z = z-1;
+                               cell.p[3].x = x-1; cell.p[3].y = y  ; cell.p[3].z = z-1;
+                               cell.p[4].x = x-1; cell.p[4].y = y-1; cell.p[4].z = z  ;
+                               cell.p[5].x = x  ; cell.p[5].y = y-1; cell.p[5].z = z  ;
+                               cell.p[6].x = x  ; cell.p[6].y = y  ; cell.p[6].z = z  ;
+                               cell.p[7].x = x-1; cell.p[7].y = y  ; cell.p[7].z = z  ;
+
+# define GRID(X,Y,WHICH) ((WHICH) \
+               ? layer1[((Y)*grid_size) + ((X))] \
+               : layer0[((Y)*grid_size) + ((X))])
+
+                               cell.val[0] = GRID (x-1, y-1, 0);
+                               cell.val[1] = GRID (x  , y-1, 0);
+                               cell.val[2] = GRID (x  , y  , 0);
+                               cell.val[3] = GRID (x-1, y  , 0);
+                               cell.val[4] = GRID (x-1, y-1, 1);
+                               cell.val[5] = GRID (x  , y-1, 1);
+                               cell.val[6] = GRID (x  , y  , 1);
+                               cell.val[7] = GRID (x-1, y  , 1);
+# undef GRID
+
+                               /* Now generate the triangles for this cubic segment,
+                                  and emit the GL faces.
+                                */
+                               ntri = march_one_cube (cell, isolevel, tri);
+                               polys += ntri;
+                               for (i = 0; i < ntri; i++)
+                               {
+                                       if (wireframe_p) glBegin (GL_LINE_LOOP);
+
+                                       /* If we're smoothing, we need to call the field function
+                                          again for each vertex (via function_normal().)  If we're
+                                          not smoothing, then we can just compute the normal from
+                                          this triangle.
+                                        */
+                                       if (!smooth_p)
+                                               do_normal (tri[i].p[0].x, tri[i].p[0].y, tri[i].p[0].z,
+                                                               tri[i].p[1].x, tri[i].p[1].y, tri[i].p[1].z,
+                                                               tri[i].p[2].x, tri[i].p[2].y, tri[i].p[2].z);
+
+# define VERT(X,Y,Z) \
+                                       if (smooth_p) \
+                                       do_function_normal ((X), (Y), (Z), compute_fn, closure2); \
+                                       glVertex3f ((X), (Y), (Z))
+
+                                       VERT (tri[i].p[0].x, tri[i].p[0].y, tri[i].p[0].z);
+                                       VERT (tri[i].p[1].x, tri[i].p[1].y, tri[i].p[1].z);
+                                       VERT (tri[i].p[2].x, tri[i].p[2].y, tri[i].p[2].z);
+# undef VERT
+                                       if (wireframe_p) glEnd ();
+                               }
+                       }
+       }
+
+       if (!wireframe_p)
+               glEnd ();
+
+       g_free(layers);
+
+       if (free_fn)
+               free_fn (closure2);
+
+       if (polygon_count)
+               *polygon_count = polys;
+}

diff --git a/src/marching.h b/src/marching.h
new file mode 100644 (file)
index 0000000..b1186f5
--- /dev/null
+++ b/src/marching.h
@@ -0,0 +1,80 @@
+/**
+ * Copyright (c) 2002 Jamie Zawinski <jwz@jwz.org>
+ * Copyright (c) 2009 Andy Spencer <andy753421@gmail.com>
+ *
+ * Utility functions to create "marching cubes" meshes from 3d fields.
+ *
+ * Permission to use, copy, modify, distribute, and sell this software and its
+ * documentation for any purpose is hereby granted without fee, provided that
+ * the above copyright notice appear in all copies and that both that
+ * copyright notice and this permission notice appear in supporting
+ * documentation.  No representations are made about the suitability of this
+ * software for any purpose.  It is provided "as is" without express or 
+ * implied warranty.
+ */
+
+#ifndef __MARCHING_H__
+#define __MARCHING_H__
+
+typedef struct {
+       double x;
+       double y;
+       double z;
+} XYZ;
+
+typedef struct {
+   XYZ p[3];
+} TRIANGLE;
+
+typedef struct {
+   XYZ p[8];
+   double val[8];
+} GRIDCELL;
+
+int march_one_cube(GRIDCELL grid, double isolevel, TRIANGLE *triangles);
+
+void do_function_normal (double x, double y, double z,
+               double (*compute_fn) (double x, double y, double z,
+                       void *closure),
+               void *c);
+
+XYZ calc_normal (XYZ p, XYZ p1, XYZ p2);
+
+void do_normal(float x1, float y1, float z1,
+               float x2, float y2, float z2,
+               float x3, float y3, float z3);
+
+
+/* Given a function capable of generating a value at any XYZ position,
+   creates OpenGL faces for the solids defined.
+
+   init_fn is called at the beginning for initial, and returns an object.
+   free_fn is called at the end.
+
+   compute_fn is called for each XYZ in the specified grid, and returns
+   the double value of that coordinate.  If smoothing is on, then
+   compute_fn will also be called twice more for each emitted vertex,
+   in order to calculate vertex normals (so don't assume it will only
+   be called with values falling on the grid boundaries.)
+
+   Points are inside an object if the are less than `isolevel', and
+   outside otherwise.
+
+   If polygon_count is specified, the number of faces generated will be
+   returned there.
+*/
+extern void
+marching_cubes (int grid_size,     /* density of the mesh */
+                double isolevel,   /* cutoff point for "in" versus "out" */
+                int wireframe_p,   /* wireframe, or solid */
+                int smooth_p,      /* smooth, or faceted */
+
+                void * (*init_fn)    (double grid_size, void *closure1),
+                double (*compute_fn) (double x, double y, double z,
+                                      void *closure2),
+                void   (*free_fn)    (void *closure2),
+                void *closure1,
+
+                unsigned long *polygon_count);
+
+#endif /* __MARCHING_H__ */

diff --git a/src/misc.h b/src/misc.h
new file mode 100644 (file)
index 0000000..ab17c74
--- /dev/null
+++ b/src/misc.h
@@ -0,0 +1,31 @@
+#define d2r(deg) (((deg)*M_PI)/180.0)
+#define r2d(rad) (((rad)*180.0)/M_PI)
+
+#define RSL_FOREACH_VOL(radar, volume, count, index) \
+       guint count = 0; \
+       for (guint index = 0; index < radar->h.nvolumes; index++) { \
+               Volume *volume = radar->v[index]; \
+               if (volume == NULL) continue; \
+               count++;
+
+#define RSL_FOREACH_SWEEP(volume, sweep, count, index) \
+       guint count = 0; \
+       for (guint index = 0; index < volume->h.nsweeps; index++) { \
+               Sweep *sweep = volume->sweep[index]; \
+               if (sweep == NULL || sweep->h.elev == 0) continue; \
+               count++;
+
+#define RSL_FOREACH_RAY(sweep, ray, count, index) \
+       guint count = 0; \
+       for (guint index = 0; index < sweep->h.nrays; index++) { \
+               Ray *ray = sweep->ray[index]; \
+               if (ray == NULL) continue; \
+               count++;
+
+#define RSL_FOREACH_BIN(ray, bin, count, index) \
+       guint count = 0; \
+       for (guint index = 0; index < ray->h.nbins; index++) { \
+               Range bin = ray->range[index]; \
+               count++;
+
+#define RSL_FOREACH_END }

diff --git a/src/plugin-radar.c b/src/plugin-radar.c
index 18e92796c0e363cd4fbd579a50e2422011edb3d9..23b6d3ab774fdf1976e6c22bba71503ce413e4a8 100644 (file)
--- a/src/plugin-radar.c
+++ b/src/plugin-radar.c
@@ -23,9 +23,11 @@
 #include <math.h>
 #include <rsl.h>
 
+#include "misc.h"
 #include "aweather-gui.h"
 #include "plugin-radar.h"
 #include "data.h"
+#include "marching.h"
 
 static char *nexrad_base = "http://mesonet.agron.iastate.edu/data/";
 
@@ -49,8 +51,10 @@ static void aweather_radar_init(AWeatherRadar *radar)
 {
        g_debug("AWeatherRadar: class_init");
        /* Set defaults */
-       radar->gui  = NULL;
-       radar->soup = NULL;
+       radar->gui           = NULL;
+       radar->soup          = NULL;
+       radar->cur_triangles = NULL;
+       radar->cur_num_triangles = 0;
 }
 static void aweather_radar_dispose(GObject *gobject)
 {
@@ -83,18 +87,17 @@ static void bscan_sweep(AWeatherRadar *self, Sweep *sweep, colormap_t *colormap,
                guint8 **data, int *width, int *height)
 {
        /* Calculate max number of bins */
-       int i, max_bins = 0;
-       for (i = 0; i < sweep->h.nrays; i++)
+       int max_bins = 0;
+       for (int i = 0; i < sweep->h.nrays; i++)
                max_bins = MAX(max_bins, sweep->ray[i]->h.nbins);
 
        /* Allocate buffer using max number of bins for each ray */
        guint8 *buf = g_malloc0(sweep->h.nrays * max_bins * 4);
 
        /* Fill the data */
-       int ri, bi;
-       for (ri = 0; ri < sweep->h.nrays; ri++) {
+       for (int ri = 0; ri < sweep->h.nrays; ri++) {
                Ray *ray  = sweep->ray[ri];
-               for (bi = 0; bi < ray->h.nbins; bi++) {
+               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]);
@@ -102,7 +105,7 @@ static void bscan_sweep(AWeatherRadar *self, Sweep *sweep, colormap_t *colormap,
                        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];
+                       buf[buf_i+3] = colormap->data[val][3]; // TESTING
                        if (val == BADVAL     || val == RFVAL      || val == APFLAG ||
                            val == NOTFOUND_H || val == NOTFOUND_V || val == NOECHO) {
                                buf[buf_i+3] = 0x00; // transparent
@@ -187,7 +190,7 @@ static void load_radar_gui(AWeatherRadar *self, Radar *radar)
                                        g_snprintf(col_label_str, 64, "<b>%.2f°</b>", elev);
                                        col_label = gtk_label_new(col_label_str);
                                        gtk_label_set_use_markup(GTK_LABEL(col_label), TRUE);
-                                       gtk_widget_set_size_request(col_label, 40, -1);
+                                       gtk_widget_set_size_request(col_label, 70, -1);
                                        gtk_table_attach(GTK_TABLE(table), col_label,
                                                        cols-1,cols, 0,1, GTK_FILL,GTK_FILL, 0,0);
                                }
@@ -217,6 +220,37 @@ static void load_radar_gui(AWeatherRadar *self, Radar *radar)
        gtk_widget_show_all(table);
 }
 
+static void _aweather_radar_grid_set(GRIDCELL *grid, int gi, Ray *ray, int bi)
+{
+       Range range = ray->range[bi];
+
+       double angle = d2r(ray->h.azimuth);
+       double tilt  = d2r(ray->h.elev);
+
+       double lx    = sin(angle);
+       double ly    = cos(angle);
+       double lz    = sin(tilt);
+
+       double dist   = bi*ray->h.gate_size + ray->h.range_bin1;
+               
+       grid->p[gi].x = lx*dist;
+       grid->p[gi].y = ly*dist;
+       grid->p[gi].z = lz*dist;
+
+       guint8 val = (guint8)ray->h.f(ray->range[bi]);
+       if (val == BADVAL     || val == RFVAL      || val == APFLAG ||
+           val == NOTFOUND_H || val == NOTFOUND_V || val == NOECHO ||
+           val > 80)
+               val = 0;
+       grid->val[gi] = (float)val;
+       //g_debug("(%.2f,%.2f,%.2f) - (%.0f,%.0f,%.0f) = %d", 
+       //      angle, tilt, dist,
+       //      grid->p[gi].x,
+       //      grid->p[gi].y,
+       //      grid->p[gi].z,
+       //      val);
+}
+
 /* Load a radar from a decompressed file */
 static void load_radar(AWeatherRadar *self, gchar *radar_file)
 {
@@ -233,6 +267,81 @@ static void load_radar(AWeatherRadar *self, gchar *radar_file)
        }
        g_free(site);
 
+#ifdef MARCHING
+       /* Load the surface */
+       if (self->cur_triangles) {
+               g_free(self->cur_triangles);
+               self->cur_triangles = NULL;
+       }
+       self->cur_num_triangles = 0;
+       int x = 1;
+       for (guint vi = 0; vi < radar->h.nvolumes; vi++) {
+               if (radar->v[vi] == NULL) continue;
+
+               for (guint si = 0; si+1 < radar->v[vi]->h.nsweeps; si++) {
+                       Sweep *sweep0 = radar->v[vi]->sweep[si+0];
+                       Sweep *sweep1 = radar->v[vi]->sweep[si+1];
+
+                       //g_debug("_aweather_radar_expose: sweep[%3d-%3d] -- nrays = %d, %d",
+                       //      si, si+1,sweep0->h.nrays, sweep1->h.nrays);
+
+                       /* Skip super->regular resolution switch for now */
+                       if (sweep0 == NULL || sweep0->h.elev == 0 ||
+                           sweep1 == NULL || sweep1->h.elev == 0 ||
+                           sweep0->h.nrays != sweep1->h.nrays)
+                               continue;
+
+                       /* We repack the arrays so that raysX[0] is always north, etc */
+                       Ray **rays0 = g_malloc0(sizeof(Ray*)*sweep0->h.nrays);
+                       Ray **rays1 = g_malloc0(sizeof(Ray*)*sweep1->h.nrays);
+
+                       for (guint ri = 0; ri < sweep0->h.nrays; ri++)
+                               rays0[(guint)(sweep0->ray[ri]->h.azimuth * sweep0->h.nrays / 360)] =
+                                       sweep0->ray[ri];
+                       for (guint ri = 0; ri < sweep1->h.nrays; ri++)
+                               rays1[(guint)(sweep1->ray[ri]->h.azimuth * sweep1->h.nrays / 360)] =
+                                       sweep1->ray[ri];
+
+                       for (guint ri = 0; ri+x < sweep0->h.nrays; ri+=x) {
+                               //g_debug("_aweather_radar_expose: ray[%3d-%3d] -- nbins = %d, %d, %d, %d",
+                               //      ri, ri+x,
+                               //      rays0[ri  ]->h.nbins, 
+                               //      rays0[ri+1]->h.nbins, 
+                               //      rays1[ri  ]->h.nbins, 
+                               //      rays1[ri+1]->h.nbins);
+
+                               for (guint bi = 0; bi+x < rays1[ri]->h.nbins; bi+=x) {
+                                       GRIDCELL grid = {};
+                                       _aweather_radar_grid_set(&grid, 7, rays0[(ri  )%sweep0->h.nrays], bi+x);
+                                       _aweather_radar_grid_set(&grid, 6, rays0[(ri+x)%sweep0->h.nrays], bi+x);
+                                       _aweather_radar_grid_set(&grid, 5, rays0[(ri+x)%sweep0->h.nrays], bi  );
+                                       _aweather_radar_grid_set(&grid, 4, rays0[(ri  )%sweep0->h.nrays], bi  );
+                                       _aweather_radar_grid_set(&grid, 3, rays1[(ri  )%sweep0->h.nrays], bi+x);
+                                       _aweather_radar_grid_set(&grid, 2, rays1[(ri+x)%sweep0->h.nrays], bi+x);
+                                       _aweather_radar_grid_set(&grid, 1, rays1[(ri+x)%sweep0->h.nrays], bi  );
+                                       _aweather_radar_grid_set(&grid, 0, rays1[(ri  )%sweep0->h.nrays], bi  );
+                                       
+                                       TRIANGLE tris[10];
+                                       int n = march_one_cube(grid, 40, tris);
+
+                                       self->cur_triangles = g_realloc(self->cur_triangles,
+                                               (self->cur_num_triangles+n)*sizeof(TRIANGLE));
+                                       for (int i = 0; i < n; i++) {
+                                               //g_debug("triangle: ");
+                                               //g_debug("\t(%f,%f,%f)", tris[i].p[0].x, tris[i].p[0].y, tris[i].p[0].z);
+                                               //g_debug("\t(%f,%f,%f)", tris[i].p[1].x, tris[i].p[1].y, tris[i].p[1].z);
+                                               //g_debug("\t(%f,%f,%f)", tris[i].p[2].x, tris[i].p[2].y, tris[i].p[2].z);
+                                               self->cur_triangles[self->cur_num_triangles+i] = tris[i];
+                                       }
+                                       self->cur_num_triangles += n;
+                                       //g_debug(" ");
+                               }
+                       }
+               }
+               break; // Exit after first volume (reflectivity)
+       }
+#endif
+
        /* Load the first sweep by default */
        if (radar->h.nvolumes < 1 || radar->v[0]->h.nsweeps < 1) {
                g_warning("No sweeps found\n");
@@ -495,8 +604,39 @@ static void _aweather_radar_expose(AWeatherPlugin *_self)
                return;
        Sweep *sweep = self->cur_sweep;
 
-       /* Draw the rays */
+#ifdef MARCHING
+       /* Draw the surface */
+       glMatrixMode(GL_MODELVIEW);
+       glPushMatrix();
+       glDisable(GL_TEXTURE_2D);
+       float light_ambient[]  = {0.1f, 0.1f, 0.0f};
+       float light_diffuse[]  = {0.9f, 0.9f, 0.9f};
+       float light_position[] = {-300000.0f, 500000.0f, 400000.0f, 1.0f};
+       glLightfv(GL_LIGHT0, GL_AMBIENT,  light_ambient);
+       glLightfv(GL_LIGHT0, GL_DIFFUSE,  light_diffuse);
+       glLightfv(GL_LIGHT0, GL_POSITION, light_position);
+       glEnable(GL_LIGHT0);
+       glEnable(GL_LIGHTING);
+       glEnable(GL_COLOR_MATERIAL);
+       glColor4f(1,1,1,0.75);
+       g_debug("ntri=%d", self->cur_num_triangles);
+       glBegin(GL_TRIANGLES);
+       for (int i = 0; i < self->cur_num_triangles; i++) {
+               TRIANGLE t = self->cur_triangles[i];
+               do_normal(t.p[0].x, t.p[0].y, t.p[0].z,
+                         t.p[1].x, t.p[1].y, t.p[1].z,
+                         t.p[2].x, t.p[2].y, t.p[2].z);
+               glVertex3f(t.p[0].x, t.p[0].y, t.p[0].z);
+               glVertex3f(t.p[1].x, t.p[1].y, t.p[1].z);
+               glVertex3f(t.p[2].x, t.p[2].y, t.p[2].z);
+       }
+       glEnd();
+       glPopMatrix();
+#endif
 
+       /* Draw the rays */
+       glDisable(GL_LIGHTING);
+       glDisable(GL_COLOR_MATERIAL);
        glMatrixMode(GL_MODELVIEW);
        glPushMatrix();
        glBindTexture(GL_TEXTURE_2D, self->cur_sweep_tex);
@@ -509,11 +649,11 @@ static void _aweather_radar_expose(AWeatherPlugin *_self)
                double angle = 0;
                if (ri < sweep->h.nrays) {
                        ray = sweep->ray[ri];
-                       angle = ((ray->h.azimuth - ((double)ray->h.beam_width/2.))*M_PI)/180.0; 
+                       angle = d2r(ray->h.azimuth - ((double)ray->h.beam_width/2.));
                } else {
                        /* Do the right side of the last sweep */
                        ray = sweep->ray[ri-1];
-                       angle = ((ray->h.azimuth + ((double)ray->h.beam_width/2.))*M_PI)/180.0; 
+                       angle = d2r(ray->h.azimuth + ((double)ray->h.beam_width/2.));
                }
 
                double lx = sin(angle);
@@ -528,8 +668,10 @@ static void _aweather_radar_expose(AWeatherPlugin *_self)
                glVertex3f(lx*near_dist, ly*near_dist, 2.0);
 
                // far  left
+               // todo: correct range-height function
+               double height = sin(d2r(ray->h.elev)) * far_dist;
                glTexCoord2f(1.0, (double)ri/sweep->h.nrays-0.01);
-               glVertex3f(lx*far_dist,  ly*far_dist,  2.0);
+               glVertex3f(lx*far_dist,  ly*far_dist, height);
        }
        //g_print("ri=%d, nr=%d, bw=%f\n", _ri, sweep->h.nrays, sweep->h.beam_width);
        glEnd();

diff --git a/src/plugin-radar.h b/src/plugin-radar.h
index ca8dbea6ce514ac93064a56256d8786b70b81f47..2db7980c7fc4edd4f23bd8ee00e67ac3c3c5f016 100644 (file)
--- a/src/plugin-radar.h
+++ b/src/plugin-radar.h
@@ -22,6 +22,8 @@
 #include <libsoup/soup.h>
 #include <rsl.h>
 
+#include "marching.h"
+
 /* TODO: convert */
 typedef struct {
        char *name;
@@ -54,6 +56,8 @@ struct _AWeatherRadar {
        Sweep       *cur_sweep;
        colormap_t  *cur_colormap;
        guint        cur_sweep_tex;
+       TRIANGLE    *cur_triangles;
+       guint        cur_num_triangles;
 };
 
 struct _AWeatherRadarClass {