Add Blue Marble Next Gen plugin and tile rendering code
authorAndy Spencer <andy753421@gmail.com>
Mon, 9 Nov 2009 12:37:58 +0000 (12:37 +0000)
committerAndy Spencer <andy753421@gmail.com>
Mon, 9 Nov 2009 12:37:58 +0000 (12:37 +0000)
Tile Splitting:
  (These functions were skipped in the previous commit for GisTile)

  To update a set of tiles, the plugin responsible for the tiles needs
  to periodically call the update function which will split any tiles
  that are below the viewers desired resolution.

  Unlike the WmsInfo method, the tiles are not automatically updated
  when a particular tile is requested. (In fact, requesting a particular
  tile is rarely done).

Tile Rendering:
  Rendering with GisTile is significantly different than with WmsInfo
  and fans. This difference is easily explained with some pseudocode:

  Old way:
    for triangle in mesh:
        for tile in find_tiles(triangle):
            render(triangle, tile.texture)

  New way:
    for tile in tiles:
        for triangle in find_triangles(tile):
            render(triangle, tile.texture)

  Both find_tiles and find_triangles are O(log n) operations in the
  worst case, but I think using find_triangles should result in faster
  code because find_triangles can quickly return a large group of tiles,
  while the find_tile function typically reruns a single tile. There's
  additional discussion of this in the source code.

  From an implementation standpoint, the new way makes it easier to add
  plugins. With the old way, the find_tiles function was invoked by a
  function pointer attached to the sphere, which makes rendering
  multiples layers of tiles difficult. In the new, each plugin just
  calls the render_tiles function which renders a set of tiles
  (recursively) on top of the sphere.

src/gis-opengl.c
src/gis-opengl.h
src/gis-tile.c
src/gis-tile.h
src/gis_test.c
src/plugins/bmng.c
src/plugins/bmng.h
src/roam.c
src/roam.h

index 03a8f2f..4d3b7e4 100644 (file)
@@ -271,6 +271,70 @@ void gis_opengl_center_position(GisOpenGL *self, gdouble lat, gdouble lon, gdoub
        glTranslatef(0, 0, elev2rad(elev));
 }
 
+void gis_opengl_render_tile(GisOpenGL *self, GisTile *tile)
+{
+       if (!tile || !tile->data)
+               return;
+       GList *triangles = roam_sphere_get_intersect(self->sphere,
+                       tile->edge.n, tile->edge.s, tile->edge.e, tile->edge.w);
+       if (!triangles)
+               g_warning("GisOpenGL: render_tiles - No triangles to draw: edges=%f,%f,%f,%f",
+                       tile->edge.n, tile->edge.s, tile->edge.e, tile->edge.w);
+       //g_message("rendering %4d triangles for tile edges=%7.2f,%7.2f,%7.2f,%7.2f",
+       //              g_list_length(triangles), tile->edge.n, tile->edge.s, tile->edge.e, tile->edge.w);
+       for (GList *cur = triangles; cur; cur = cur->next) {
+               RoamTriangle *tri = cur->data;
+
+               gdouble lat[3] = {tri->p.r->lat, tri->p.m->lat, tri->p.l->lat};
+               gdouble lon[3] = {tri->p.r->lon, tri->p.m->lon, tri->p.l->lon};
+
+               if (lon[0] < -90 || lon[1] < -90 || lon[2] < -90) {
+                       if (lon[0] > 90) lon[0] -= 360;
+                       if (lon[1] > 90) lon[1] -= 360;
+                       if (lon[2] > 90) lon[2] -= 360;
+               }
+
+               gdouble n = tile->edge.n;
+               gdouble s = tile->edge.s;
+               gdouble e = tile->edge.e;
+               gdouble w = tile->edge.w;
+
+               gdouble londist = e - w;
+               gdouble latdist = n - s;
+
+               gdouble xy[][3] = {
+                       {(lon[0]-w)/londist, 1-(lat[0]-s)/latdist},
+                       {(lon[1]-w)/londist, 1-(lat[1]-s)/latdist},
+                       {(lon[2]-w)/londist, 1-(lat[2]-s)/latdist},
+               };
+
+               glEnable(GL_TEXTURE_2D);
+               glBindTexture(GL_TEXTURE_2D, *(guint*)tile->data);
+               glBegin(GL_TRIANGLES);
+               glNormal3dv(tri->p.r->norm); glTexCoord2dv(xy[0]); glVertex3dv((double*)tri->p.r);
+               glNormal3dv(tri->p.m->norm); glTexCoord2dv(xy[1]); glVertex3dv((double*)tri->p.m);
+               glNormal3dv(tri->p.l->norm); glTexCoord2dv(xy[2]); glVertex3dv((double*)tri->p.l);
+               glEnd();
+       }
+       g_list_free(triangles);
+}
+
+void gis_opengl_render_tiles(GisOpenGL *opengl, GisTile *tile)
+{
+       /* Only render children if possible */
+       gboolean has_children = TRUE;
+       GisTile *child;
+       gis_tile_foreach(tile, child)
+               if (!child || !child->data)
+                       has_children = FALSE;
+       if (has_children)
+               /* Only render children */
+               gis_tile_foreach(tile, child)
+                       gis_opengl_render_tiles(opengl, child);
+       else
+               /* No children, render this tile */
+               gis_opengl_render_tile(opengl, tile);
+}
 void gis_opengl_redraw(GisOpenGL *self)
 {
        g_debug("GisOpenGL: redraw");
index 49562d1..19feefc 100644 (file)
@@ -36,6 +36,7 @@ typedef struct _GisOpenGLClass GisOpenGLClass;
 #include "gis-view.h"
 #include "gis-world.h"
 #include "gis-plugin.h"
+#include "gis-tile.h"
 #include "roam.h"
 
 struct _GisOpenGL {
@@ -66,6 +67,10 @@ GisOpenGL *gis_opengl_new(GisWorld *world, GisView *view, GisPlugins *plugins);
 void gis_opengl_center_position(GisOpenGL *opengl,
                gdouble lat, gdouble lon, gdouble elev);
 
+void gis_opengl_render_tile(GisOpenGL *opengl, GisTile *tile);
+
+void gis_opengl_render_tiles(GisOpenGL *opengl, GisTile *root);
+
 void gis_opengl_redraw(GisOpenGL *opengl);
 
 void gis_opengl_begin(GisOpenGL *opengl);
index ab45616..80772e6 100644 (file)
@@ -35,6 +35,7 @@ GisTile *gis_tile_new(GisTile *parent,
        self->edge.s = s;
        self->edge.e = e;
        self->edge.w = w;
+       self->atime  = time(NULL);
        return self;
 }
 
@@ -54,6 +55,102 @@ gchar *gis_tile_get_path(GisTile *child)
        return g_string_free(path, FALSE);
 }
 
+gdouble _gis_tile_get_min_dist(GisTile *self,
+               gdouble lat, gdouble lon, gdouble elev)
+{
+       gdouble tlat  = lat > self->edge.n ? self->edge.n :
+                        lat < self->edge.s ? self->edge.s : lat;
+       gdouble tlon  = lon > self->edge.e ? self->edge.e :
+                       lon < self->edge.w ? self->edge.w : lon;
+       gdouble telev = 0; // TODO: elevation at rlat,rlon
+       //if (lat == tlat && lon == tlon)
+       //      return elev; /* Shortcut? */
+       gdouble a[3], b[3];
+       lle2xyz( lat,  lon,  elev, a+0, a+1, a+2);
+       lle2xyz(tlat, tlon, telev, b+0, b+1, b+2);
+       return distd(a, b);
+}
+
+gboolean _gis_tile_needs_split(GisTile *self,
+               gdouble max_res, gint width, gint height,
+               gdouble lat, gdouble lon, gdouble elev)
+{
+       gdouble lat_point = self->edge.n < 0 ? self->edge.n :
+                           self->edge.s > 0 ? self->edge.s : 0;
+       gdouble min_dist  = _gis_tile_get_min_dist(self, lat, lon, elev);
+       gdouble view_res  = MPPX(min_dist);
+       gdouble lon_dist  = self->edge.e - self->edge.w;
+       gdouble tile_res  = ll2m(lon_dist, lat_point)/width;
+
+       /* This isn't really right, but it helps with memory since we don't (yet?) test if the tile
+        * would be drawn */
+       gdouble scale = elev / min_dist;
+       view_res /= scale;
+       //g_message("tile=(%7.2f %7.2f %7.2f %7.2f) "
+       //          "eye=(%9.1f %9.1f %9.1f) "
+       //          "elev=%9.1f / dist=%9.1f = %f",
+       //              self->edge.n, self->edge.s, self->edge.e, self->edge.w,
+       //              lat, lon, elev,
+       //              elev, min_dist, scale);
+
+       if (tile_res < max_res)
+               return FALSE;
+       return view_res < tile_res;
+}
+
+void gis_tile_update(GisTile *self,
+               gdouble res, gint width, gint height,
+               gdouble lat, gdouble lon, gdouble elev,
+               GisTileLoadFunc load_func, gpointer user_data)
+{
+       self->atime = time(NULL);
+       //g_debug("GisTile: update - %p->atime = %u", self, (guint)self->atime);
+       gdouble lat_dist = self->edge.n - self->edge.s;
+       gdouble lon_dist = self->edge.e - self->edge.w;
+       if (_gis_tile_needs_split(self, res, width, height, lat, lon, elev)) {
+               gdouble lat_step = lat_dist / G_N_ELEMENTS(self->children);
+               gdouble lon_step = lon_dist / G_N_ELEMENTS(self->children[0]);
+               int x, y;
+               gis_tile_foreach_index(self, x, y) {
+                       if (!self->children[x][y]) {
+                               self->children[x][y] = gis_tile_new(self,
+                                               self->edge.n-(lat_step*(x+0)),
+                                               self->edge.n-(lat_step*(x+1)),
+                                               self->edge.w+(lon_step*(y+1)),
+                                               self->edge.w+(lon_step*(y+0)));
+                               load_func(self->children[x][y], user_data);
+                       }
+                       gis_tile_update(self->children[x][y],
+                                       res, width, height,
+                                       lat, lon, elev,
+                                       load_func, user_data);
+               }
+       }
+}
+
+GisTile *gis_tile_gc(GisTile *self, time_t atime,
+               GisTileFreeFunc free_func, gpointer user_data)
+{
+       if (!self)
+               return NULL;
+       gboolean has_children = FALSE;
+       int x, y;
+       gis_tile_foreach_index(self, x, y) {
+               self->children[x][y] = gis_tile_gc(
+                               self->children[x][y], atime,
+                               free_func, user_data);
+               if (self->children[x][y])
+                       has_children = TRUE;
+       }
+       //g_debug("GisTile: gc - %p->atime=%u < atime=%u",
+       //              self, (guint)self->atime, (guint)atime);
+       if (!has_children && self->atime < atime && self->data) {
+               free_func(self, user_data);
+               g_free(self);
+               return NULL;
+       }
+       return self;
+}
 
 void gis_tile_free(GisTile *self, GisTileFreeFunc free_func, gpointer user_data)
 {
index 12d397d..cda4442 100644 (file)
@@ -47,6 +47,9 @@ struct _GisTile {
        /* Pointers to parent/child nodes */
        GisTile *parent;
        GisTile *children[2][2];
+
+       /* Last access time (for garbage collection) */
+       time_t atime;
 };
 
 /* Path to string table, keep in sync with tile->children */ 
@@ -59,6 +62,17 @@ GisTile *gis_tile_new(GisTile *parent,
 /* Return a string representation of the tile's path */
 gchar *gis_tile_get_path(GisTile *child);
 
+/* Update a root tile */
+/* Based on eye distance (lat,lon,elev) */
+void gis_tile_update(GisTile *root,
+               gdouble res, gint width, gint height,
+               gdouble lat, gdouble lon, gdouble elev,
+               GisTileLoadFunc load_func, gpointer user_data);
+
+/* Delete nodes that haven't been accessed since atime */
+GisTile *gis_tile_gc(GisTile *root, time_t atime,
+               GisTileFreeFunc free_func, gpointer user_data);
+
 /* Free a tile and all it's children */
 void gis_tile_free(GisTile *root,
                GisTileFreeFunc free_func, gpointer user_data);
index 9d9fe9a..35f2901 100644 (file)
@@ -51,18 +51,16 @@ int main(int argc, char **argv)
        GisView    *view    = gis_view_new();
        GisOpenGL  *opengl  = gis_opengl_new(world, view, plugins);
 
-       //gis_plugins_load(plugins, "radar",   world, view, opengl, prefs);
-       //gis_plugins_load(plugins, "ridge",   world, view, opengl, prefs);
-       gis_plugins_load(plugins, "bmng", world, view, opengl, prefs);
-       gis_plugins_load(plugins, "srtm", world, view, opengl, prefs);
-
        GtkWidget  *window  = gtk_window_new(GTK_WINDOW_TOPLEVEL);
        g_signal_connect(window,  "destroy",         G_CALLBACK(gtk_main_quit), NULL);
        g_signal_connect(window,  "key-press-event", G_CALLBACK(on_key_press),  NULL);
        gtk_container_add(GTK_CONTAINER(window), GTK_WIDGET(opengl));
        gtk_widget_show_all(window);
+       gis_plugins_load(plugins, "bmng", world, view, opengl, prefs);
+       //gis_plugins_load(plugins, "srtm", world, view, opengl, prefs);
 
        gis_view_set_site(view, "KLSX");
+
        gtk_main();
 
        g_object_unref(prefs);
index bc61d06..a1499c6 100644 (file)
  * along with this program.  If not, see <http://www.gnu.org/licenses/>.
  */
 
-#include <gtk/gtkgl.h>
+#include <time.h>
 #include <GL/gl.h>
 
 #include <gis.h>
 
 #include "bmng.h"
 
-/***********
- * Helpers *
- ***********/
-static gboolean rotate(gpointer _self)
+#define MAX_RESOLUTION 500
+#define TILE_WIDTH     1024
+#define TILE_HEIGHT    512
+
+static void _load_tile(GisTile *tile, gpointer _self)
+{
+       GisPluginBmng *self = _self;
+       g_debug("GisPluginBmng: _load_tile start");
+
+       char *path = gis_wms_make_local(self->wms, tile);
+       GdkPixbuf *pixbuf = gdk_pixbuf_new_from_file(path, NULL);
+       g_free(path);
+
+       guchar   *pixels = gdk_pixbuf_get_pixels(pixbuf);
+       gboolean  alpha  = gdk_pixbuf_get_has_alpha(pixbuf);
+       gint      width  = gdk_pixbuf_get_width(pixbuf);
+       gint      height = gdk_pixbuf_get_height(pixbuf);
+
+       guint *tex = g_new0(guint, 1);
+       gis_opengl_begin(self->opengl);
+       glGenTextures(1, tex);
+       glBindTexture(GL_TEXTURE_2D, *tex);
+
+       glPixelStorei(GL_UNPACK_ALIGNMENT, 1);
+       glPixelStorei(GL_PACK_ALIGNMENT, 1);
+       glTexImage2D(GL_TEXTURE_2D, 0, 4, width, height, 0,
+                       (alpha ? GL_RGBA : GL_RGB), GL_UNSIGNED_BYTE, pixels);
+       glTexParameterf(GL_TEXTURE_2D, GL_TEXTURE_MAG_FILTER, GL_NEAREST);
+       glTexParameterf(GL_TEXTURE_2D, GL_TEXTURE_MIN_FILTER, GL_NEAREST);
+       glTexParameteri(GL_TEXTURE_2D, GL_TEXTURE_WRAP_S, GL_CLAMP_TO_BORDER);
+       glTexParameteri(GL_TEXTURE_2D, GL_TEXTURE_WRAP_T, GL_CLAMP_TO_BORDER);
+       glFlush();
+       gis_opengl_end(self->opengl);
+
+       tile->data = tex;
+       gis_opengl_redraw(self->opengl);
+       g_object_unref(pixbuf);
+}
+
+static void _free_tile(GisTile *tile, gpointer _self)
+{
+       GisPluginBmng *self = _self;
+       g_debug("GisPluginBmng: _free_tile: %p=%d", tile->data, *(guint*)tile->data);
+       guint *data = tile->data;
+       glDeleteTextures(1, data);
+       g_free(data);
+}
+
+static gpointer _update_tiles(gpointer _self)
 {
+       g_debug("GisPluginBmng: _update_tiles");
        GisPluginBmng *self = _self;
-       if (gtk_toggle_button_get_active(self->button)) {
-               self->rotation += 1.0;
-               gis_opengl_redraw(self->opengl);
-       }
-       return TRUE;
+       gdouble lat, lon, elev;
+       gis_view_get_location(self->view, &lat, &lon, &elev);
+       gis_tile_update(self->tiles,
+                       MAX_RESOLUTION, TILE_WIDTH, TILE_WIDTH,
+                       lat, lon, elev,
+                       _load_tile, self);
+       gis_tile_gc(self->tiles, time(NULL)-10,
+                       _free_tile, self);
+       return NULL;
 }
 
+/*************
+ * Callbacks *
+ *************/
+static void _on_location_changed(GisView *view, gdouble lat, gdouble lon, gdouble elev,
+               GisPluginBmng *self)
+{
+       _update_tiles(self);
+}
 
 /***********
  * Methods *
@@ -43,42 +101,24 @@ GisPluginBmng *gis_plugin_bmng_new(GisWorld *world, GisView *view, GisOpenGL *op
 {
        g_debug("GisPluginBmng: new");
        GisPluginBmng *self = g_object_new(GIS_TYPE_PLUGIN_BMNG, NULL);
+       self->view   = view;
        self->opengl = opengl;
 
-       return self;
-}
+       /* Load initial tiles */
+       _load_tile(self->tiles, self);
+       _update_tiles(self);
 
-static GtkWidget *gis_plugin_bmng_get_config(GisPlugin *_self)
-{
-       GisPluginBmng *self = GIS_PLUGIN_BMNG(_self);
-       return GTK_WIDGET(self->button);
+       /* Connect signals */
+       g_signal_connect(self->view, "location-changed", G_CALLBACK(_on_location_changed), self);
+
+       return self;
 }
 
 static void gis_plugin_bmng_expose(GisPlugin *_self)
 {
        GisPluginBmng *self = GIS_PLUGIN_BMNG(_self);
        g_debug("GisPluginBmng: expose");
-
-       glMatrixMode(GL_PROJECTION);
-       glLoadIdentity();
-       glOrtho(1,-1, -1,1, -10,10);
-
-       glMatrixMode(GL_MODELVIEW);
-       glLoadIdentity();
-
-       float light_ambient[]  = {0.1f, 0.1f, 0.0f, 1.0f};
-       float light_diffuse[]  = {0.9f, 0.9f, 0.9f, 1.0f};
-       float light_position[] = {-30.0f, 50.0f, 40.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_COLOR_MATERIAL);
-
-       glTranslatef(-0.5, -0.5, -2);
-       glRotatef(self->rotation, 1, 1, 0);
-       glColor4f(0.9, 0.9, 0.7, 1.0);
-       glDisable(GL_CULL_FACE);
-       gdk_gl_draw_teapot(TRUE, 0.25);
+       gis_opengl_render_tiles(self->opengl, self->tiles);
 }
 
 
@@ -94,24 +134,22 @@ static void gis_plugin_bmng_plugin_init(GisPluginInterface *iface)
 {
        g_debug("GisPluginBmng: plugin_init");
        /* Add methods to the interface */
-       iface->expose     = gis_plugin_bmng_expose;
-       iface->get_config = gis_plugin_bmng_get_config;
+       iface->expose = gis_plugin_bmng_expose;
 }
 /* Class/Object init */
 static void gis_plugin_bmng_init(GisPluginBmng *self)
 {
        g_debug("GisPluginBmng: init");
        /* Set defaults */
-       self->button    = GTK_TOGGLE_BUTTON(gtk_toggle_button_new_with_label("Rotate"));
-       self->rotate_id = g_timeout_add(1000/60, rotate, self);
-       self->rotation  = 30.0;
-       self->opengl    = NULL;
+       self->tiles = gis_tile_new(NULL, NORTH, SOUTH, EAST, WEST);
+       self->wms   = gis_wms_new(
+               "http://www.nasa.network.com/wms", "bmng200406", "image/jpeg",
+               "bmng", ".jpg", TILE_WIDTH, TILE_HEIGHT);
 }
 static void gis_plugin_bmng_dispose(GObject *gobject)
 {
        g_debug("GisPluginBmng: dispose");
        GisPluginBmng *self = GIS_PLUGIN_BMNG(gobject);
-       g_source_remove(self->rotate_id);
        /* Drop references */
        G_OBJECT_CLASS(gis_plugin_bmng_parent_class)->dispose(gobject);
 }
@@ -120,6 +158,8 @@ static void gis_plugin_bmng_finalize(GObject *gobject)
        g_debug("GisPluginBmng: finalize");
        GisPluginBmng *self = GIS_PLUGIN_BMNG(gobject);
        /* Free data */
+       gis_tile_free(self->tiles, _free_tile, self);
+       gis_wms_free(self->wms);
        G_OBJECT_CLASS(gis_plugin_bmng_parent_class)->finalize(gobject);
 
 }
index 73701d8..21c562a 100644 (file)
@@ -34,10 +34,10 @@ struct _GisPluginBmng {
        GObject parent_instance;
 
        /* instance members */
-       GtkToggleButton *button;
-       guint            rotate_id;
-       float            rotation;
-       GisOpenGL       *opengl;
+       GisView   *view;
+       GisOpenGL *opengl;
+       GisTile   *tiles;
+       GisWms    *wms;
 };
 
 struct _GisPluginBmngClass {
index 54ae2c6..9f2474f 100644 (file)
@@ -156,6 +156,13 @@ RoamTriangle *roam_triangle_new(RoamPoint *l, RoamPoint *m, RoamPoint *r)
        self->split->height_func = m->height_func;
        self->split->height_data = m->height_data;
        roam_point_update_height(self->split);
+       //if ((float)self->split->lat > 44 && (float)self->split->lat < 46)
+       //      g_debug("RoamTriangle: new - (l,m,r,split).lats = %7.2f %7.2f %7.2f %7.2f",
+       //                      l->lat, m->lat, r->lat, self->split->lat);
+       //if ((float)l->lat == (float)r->lat &&
+       //    (float)self->split->lat != (float)l->lat)
+       //      g_warning("RoamTriangle: new - split.lat=%f != (l=r).lat=%f",
+       //              self->split->lat, l->lat);
 
        /* Update normal */
        double pa[3];
@@ -180,6 +187,36 @@ RoamTriangle *roam_triangle_new(RoamPoint *l, RoamPoint *m, RoamPoint *r)
        self->norm[1] /= total;
        self->norm[2] /= total;
 
+       /* Store bounding box, for get_intersect */
+       RoamPoint *m1 = roam_point_new((l->x+m->x)/2, (l->y+m->y)/2, (l->z+m->z)/2);
+       RoamPoint *m2 = roam_point_new((m->x+r->x)/2, (m->y+r->y)/2, (m->z+r->z)/2);
+       RoamPoint *m3 = roam_point_new((r->x+l->x)/2, (r->y+l->y)/2, (r->z+l->z)/2);
+       RoamPoint *p[] = {l,m,r,m1,m2,m3};
+       self->edge.n =  -90; self->edge.s =  90;
+       self->edge.e = -180; self->edge.w = 180;
+       gboolean maxed = FALSE;
+       for (int i = 0; i < G_N_ELEMENTS(p); i++) {
+               self->edge.n = MAX(self->edge.n, p[i]->lat);
+               self->edge.s = MIN(self->edge.s, p[i]->lat);
+               if (p[i]->lat == 90 || p[i]->lat == -90)
+                       continue;
+               if (p[i]->lon == 180) {
+                       maxed = TRUE;
+                       continue;
+               }
+               self->edge.e = MAX(self->edge.e, p[i]->lon);
+               self->edge.w = MIN(self->edge.w, p[i]->lon);
+       }
+       if (maxed) {
+               if (self->edge.e < 0)
+                       self->edge.w = -180;
+               else
+                       self->edge.e =  180;
+       }
+       g_free(m1);
+       g_free(m2);
+       g_free(m3);
+
        //g_message("roam_triangle_new: %p", self);
        return self;
 }
@@ -289,10 +326,10 @@ void roam_triangle_split(RoamTriangle *self, RoamSphere *sphere)
 
        /* Add new triangles */
        RoamPoint *mid = self->split;
-       RoamTriangle *sl = roam_triangle_new(self->p.m, mid, self->p.l); // Self Left
-       RoamTriangle *sr = roam_triangle_new(self->p.r, mid, self->p.m); // Self Right
-       RoamTriangle *bl = roam_triangle_new(base->p.m, mid, base->p.l); // Base Left
-       RoamTriangle *br = roam_triangle_new(base->p.r, mid, base->p.m); // Base Right
+       RoamTriangle *sl = self->kids[0] = roam_triangle_new(self->p.m, mid, self->p.l); // Self Left
+       RoamTriangle *sr = self->kids[1] = roam_triangle_new(self->p.r, mid, self->p.m); // Self Right
+       RoamTriangle *bl = base->kids[0] = roam_triangle_new(base->p.m, mid, base->p.l); // Base Left
+       RoamTriangle *br = base->kids[1] = roam_triangle_new(base->p.r, mid, base->p.m); // Base Right
 
        /*                tri,l,  base,      r,  sphere */
        roam_triangle_add(sl, sr, self->t.l, br, sphere);
@@ -410,6 +447,10 @@ void roam_diamond_merge(RoamDiamond *self, RoamSphere *sphere)
        roam_triangle_remove(kids[2], sphere);
        roam_triangle_remove(kids[3], sphere);
 
+       /* Clear kids */
+       tri->kids[0]  = tri->kids[1]  = NULL;
+       base->kids[0] = base->kids[1] = NULL;
+
        /* Add/Remove triangles */
        if (tri->t.l->t.l == tri->t.r->t.r &&
            tri->t.l->t.l != tri && tri->parent) {
@@ -485,6 +526,12 @@ RoamSphere *roam_sphere_new()
                        self->roots[_triangles[i][1][1]],
                        self->roots[_triangles[i][1][2]],
                        self);
+       for (int i = 0; i < 8; i++)
+               g_debug("RoamSphere: new - %p edge=%f,%f,%f,%f", self->roots[i],
+                               self->roots[i]->edge.n,
+                               self->roots[i]->edge.s,
+                               self->roots[i]->edge.e,
+                               self->roots[i]->edge.w);
 
        return self;
 }
@@ -571,6 +618,63 @@ void roam_sphere_draw_normals(RoamSphere *self)
        g_debug("RoamSphere: draw_normal");
        g_pqueue_foreach(self->triangles, (GFunc)roam_triangle_draw_normal, NULL);
 }
+GList *_roam_sphere_get_leaves(RoamTriangle *tri, GList *list)
+{
+       if (tri->kids[0] && tri->kids[1]) {
+               list = _roam_sphere_get_leaves(tri->kids[0], list);
+               list = _roam_sphere_get_leaves(tri->kids[1], list);
+               return list;
+       } else {
+               return g_list_append(list, tri);
+       }
+}
+GList *_roam_sphere_get_intersect_rec(RoamTriangle *tri, GList *list,
+               gdouble n, gdouble s, gdouble e, gdouble w)
+{
+       gdouble tn = tri->edge.n;
+       gdouble ts = tri->edge.s;
+       gdouble te = tri->edge.e;
+       gdouble tw = tri->edge.w;
+
+       gboolean debug = FALSE  &&
+               n==90 && s==45 && e==-90 && w==-180 &&
+               ts > 44 && ts < 46;
+
+       if (debug)
+               g_message("t=%p: %f < %f || %f > %f || %f < %f || %f > %f",
+                           tri, tn,   s,   ts,   n,   te,   w,   tw,   e);
+       if (tn < s || ts > n || te < w || tw > e) {
+               /* No intersect */
+               if (debug) g_message("no intersect");
+               return list;
+       } else if (tn < n && ts > s && te < e && tw > w) {
+               /* Triangle is completely contained */
+               if (debug) g_message("contained");
+               return _roam_sphere_get_leaves(tri, list);
+       } else if (tri->kids[0] && tri->kids[1]) {
+               /* Paritial intersect with children */
+               if (debug) g_message("edge w/ child");
+               list = _roam_sphere_get_intersect_rec(tri->kids[0], list, n, s, e, w);
+               list = _roam_sphere_get_intersect_rec(tri->kids[1], list, n, s, e, w);
+               return list;
+       } else {
+               /* This triangle is an edge case */
+               if (debug) g_message("edge");
+               return g_list_append(list, tri);
+       }
+}
+GList *roam_sphere_get_intersect(RoamSphere *self,
+               gdouble n, gdouble s, gdouble e, gdouble w)
+{
+       /* I think this is correct..
+        * i_cost = cost for triangle-rectagnle intersect test
+        * time = n_tiles * 2*tris_per_tile * i_cost
+        * time = 30      * 2*333           * i_cost = 20000 * i_cost */
+       GList *list = NULL;
+       for (int i = 0; i < G_N_ELEMENTS(self->roots); i++)
+               list = _roam_sphere_get_intersect_rec(self->roots[i], list, n, s, e, w);
+       return list;
+}
 void roam_sphere_free_tri(RoamTriangle *tri)
 {
        if (--tri->p.l->tris == 0) g_free(tri->p.l);
index da4f508..94b205f 100644 (file)
@@ -71,6 +71,10 @@ struct _RoamTriangle {
        double norm[3];
        double error;
        GPQueueHandle handle;
+
+       /* For get_intersect */
+       struct { gdouble n,s,e,w; } edge;
+       RoamTriangle *kids[2];
 };
 RoamTriangle *roam_triangle_new(RoamPoint *l, RoamPoint *m, RoamPoint *r);
 void roam_triangle_add(RoamTriangle *triangle,
@@ -119,6 +123,8 @@ void roam_sphere_merge_one(RoamSphere *sphere);
 gint roam_sphere_split_merge(RoamSphere *sphere);
 void roam_sphere_draw(RoamSphere *sphere);
 void roam_sphere_draw_normals(RoamSphere *sphere);
+GList *roam_sphere_get_intersect(RoamSphere *sphere,
+               gdouble n, gdouble s, gdouble e, gdouble w);
 void roam_sphere_free(RoamSphere *sphere);
 
 #endif