We are excited to release a new version of the OpenGeo Suite! In order to capture the many improvements and bug fixes happening in the open source community, we are moving toward a more rapid release cycle. For example, GeoServer now has JDBC datastore session startup/teardown SQL comments, as well as support for paletted PNG images with alpha transparency.

In GeoExplorer (which really is pretty amazing, if you haven’t seen it recently) there is now smoother tile display, including fade-in. Also, the map legend has now been integrated directly into the layer tree. Finally, we have changed the default base layer to be MapQest OSM, moving away from Google (though Google base layers are still available).

All of these new features are available in the Community Edition, Enterprise Edition (which includes a free 30 day trial of our support), and all Cloud Editions! Try any version you’d like and contact us to purchase the support you need to put your project into production!

Another feature that has received relatively scant attention in the upcoming release is the ability to index in higher dimensions than the usual two.

Initially I had hoped to made n-dimensional indexing the default, so that multi-dimensional indexes would just automagically be available. Unfortunately it turned out that the overhead involved in having a n-dimensional index key caused the index to be slower.

Rather than slow down the standard use case (2D searching), I re-instated fixed two-dimensional indexing as the default for geometry, and added a new N-D index type for those who want the thrill of extra dimensions.

To build an N-D index, use the following CREATE INDEX statement:

CREATE INDEX my_nd_index
ON my_table USING GIST (geom gist_geometry_ops_nd);

To query your N-D index, use the &&& operator (note! three ampersands!), like so:

SELECT * FROM my_table
WHERE geom &&& 'LINESTRING(0 0 0, 2 2 2)';

(“Hey, what’s with the LINESTRING?” Because the index-only query is bounding-volume based, a linestring from the lower-left to the upper-right of my desired bounding volume will have exactly the bounding volume I want.)

Some of the fancy new 3D functions like ST_3DDistance can make use of the 3D index for fast searches in high dimensional spaces. I imagine folks working with XYT data may also find the new indexes useful.

Welcome to another edition of OpenGeo Connections. Today we’re sitting down with David Winslow, technical lead on GeoNode. David works out of OpenGeo’s NYC headquarters, is a member of the GeoNode Project Steering Committee, and a core commiter to GeoServer. David has been with OpenGeo for over four years and is an integral part of our GeoNode team.

OpenGeo’s David Dubovsky: Hi David, thanks for sitting down to chat.
David Winslow: It’s my pleasure.

(DD): Let’s jump right in –  where are you from, where are you now, and how long have you been with OpenGeo?
(DW): I’m a Duke University comp-sci alum, I grew up in Virginia and now I live in Brooklyn. I’ve been with OpenGeo for about four and a half years now.

(DD): Four and a half years? Does that make you an elder statesman at OpenGeo?
(DW): We’ve grown a lot over that time but I wouldn’t say I’m an elder anything.

(DD): So, what brought you to join OpenGeo?
(DW): Before I came here I was studying computer science at Duke. I was attending a carrer fair when  I met Vanessa Hamer, now our Vice President of Operations. The rest is history, I found her pitch intriguing as I wanted to work at a place that was devoted to open source and open data.

(DD): So you studied computer science? What about web mapping and open source?
(DW): At the time mapping was new to me. When I joined OpenPlans, which OpenGeo is a division of, we were called The Open Planning Project (TOPP). I was hired to work on a OpenCore, a tool that provided basic project infrastructure for those rallying around a cause – it included things like a wiki, calendar, blog etc.

(DD): And how did web mapping get onto your resume?
(DW): I really consider myself more of a nerd than a geo-nerd. Chris Holmes had helped develop GeoServer while at TOPP. OpenCore was one of his first implementations of GeoServer. The idea was for users would to utilize maps to share and disperse information. This would enable volunteers to identify points of interest, concern, in need of attention etc.

(DD): And Chris pulled you in to help with that?
(DW): Yes, I helped Chris integrate GeoServer, and that was my first experience with open source geospatial.

(DD): That’s great. So what are you primarily working on these days?
(DW): For OpenGeo I work mostly on GeoNode, both supporting the community work – making sure releases come out, helping with new development etc, and for OpenGeo I perform custom deployments and bringing new features into the project. Recent deployments required customization that would be beneficial to the broader community; right now I’m working on brigning those into main GeoNode project

In my spare time I’ve been working on GeoServer styling with CSS. Right now you can only edit styles for GeoServer in an XML format called SLD. CSS is a lot less complex to edit and would would really be a benefit to the community.

(DD): Interesting, the word is that your GeoNode work has had you traveling a lot, where have you been?
(DW): You can say that, last year I traveled to Fiji, Portland, Zurich, Jakarta, Denver for FOSS4G, the list goes on.

(DD): Wow, what were you doing in the South Pacific?
(DW): I was in Fiji for a GeoNode deployment at the Secretariat of the Pacific Community (SOPAC). The SPC acts as an information broker for all the South Pacific nations. They’re using GeoNode to serve up data they’ve collected and the results of a large study on disater management

In Jakarta I worked with The Australia-Indonesia Facility for Disaster Reduction (AIFDR) on a project called Risk in a box. They’re helping inform people who don’t have the ability to make disaster models what would happen if a disaster did strike their area. For example, what would happen to their community in the event of a four meter storm surge. The online calculator is based on GeoNode and all of the calculations are published as GeoNode data sets and layers.

(DD): Very interesting, has the travel slowed down at all
(DW): Well OpenGeo will be in Washington D.C. for FOSS4G North America and then in Portland for an organization wide get-together - so no, but the flights are shorter.

(DD): So what’s your claim to fame in OpenGeo?
(DW): I’ve touched at least one line of code in every piece of the OpenGeo stack – PostGIS, GeoExt, GeoServer, REST.API, KML etc. I don’t think many others in our organization can say that.

(DD): Well, I’m sure you’ll get some competition for that now. Thanks so much for letting us pick your brain.
(DW): I’m happy to talk any time.

If you want to get in touch with David look out for dwins in #geoserver or #geonode. He’s also active on the the GeoServer and GeoNode mailing lists.

The new raster type gets most of the press in the PostGIS 2.0 features round-ups, but one of my personal favorites is the support for “typmod” on geometry objects.

A “typmod” is a “type modifier”. If you’ve used SQL databases you’re probably familiar with type modifiers on variable length character columns, like this one:

CREATE TABLE test (
  description VARCHAR(255)
)

The “255″ after the “VARCHAR” is a modifier, specifying more information about the column. For geometry columns in PostGIS 2.0, it’s now possible to define a column like this:

CREATE TABLE geotbl (
  geom GEOMETRY(PointZ, 4326)
)

Note the modifiers on the geometry, specifying the specific type “Point”, the dimensionality “Z” and the spatial reference identifier “4326″. And the table description reflects the typmod information:

           Table "public.geotbl"
 Column |         Type          | Modifiers
--------+-----------------------+-----------
 geom   | geometry(PointZ,4326) | 

So what, you say? So, now the “geometry_columns” table can automatically be kept up to date as a view on the system tables. Once a spatial table is created, it now automatically appears in “geometry_columns”.

 f_table_catalog | f_table_schema | f_table_name | f_geometry_column | coord_dimension | srid  |      type
-----------------+----------------+--------------+-------------------+-----------------+-------+-----------------
 cded            | public         | geotbl       | geom              |               3 |  4326 | POINT

Still not impressed? The classic problem, “how do I reproject the data in my table” used to require multiple steps: drop spatial reference constraint, update the column, re-add the constraint, refresh the geometry_columns table. Now, it’s a one-liner, converting the contraints, global metadata, and table data in one step.

ALTER TABLE geotbl
  ALTER COLUMN geom
  SET DATA TYPE geometry(Point,26910)
  USING ST_Transform(ST_Force_2D(geom),26910)

And as a bonus, we also converted the type from 3D to 2D. And all the metadata is up-to-date automagically.

           Table "public.geotbl"
 Column |         Type          | Modifiers
--------+-----------------------+-----------
 geom   | geometry(Point,26910) | 

The manual management of the “geometry_columns” table has bothered me ever since it was introduced. After only 10 years, I’m pleased to see it finally fixed! This update required the new extended typmod support brought in to PostgreSQL 8.4 and was prototyped on the “geography” type in PostGIS 1.5.

A PostGIS user asked me an interesting question this week, looking for the intersection of these two polygons:

POLYGON((170 50,170 72,-130 72,-130 50,170 50))
POLYGON((-170 68,-170 90,-141 90,-141 68,-170 68))

Now, from the coordinate sizes, you can see these shapes are described in longitude/latitude. So naturally you should consider using the PostGIS “geography” type, particularly since the shapes are quite big.

However, when you plug the shape into SQL and try out the calculation:

SELECT ST_AsText(
  ST_Intersection(
    ST_GeogFromText('POLYGON((170 50,170 72,-130 72,-130 50,170 50))'),
    ST_GeogFromText('POLYGON((-170 68,-170 90,-141 90,-141 68,-170 68))')
  ));

You get a bad result:

ERROR: transform: couldn't project point (-170 90 0): tolerance condition error (-20)

Hmm. Transform? This sounds like reprojection. We’re working in geography though, what’s going on? Check the definition of the ST_Intersects(geography, geography) function.

CREATE OR REPLACE FUNCTION ST_Intersection(geography, geography)
  RETURNS geography
  AS 'SELECT
    geography(
      ST_Transform(
        ST_Intersection(
          ST_Transform(
            geometry($1),
            _ST_BestSRID($1, $2)),
          ST_Transform(
            geometry($2),
            _ST_BestSRID($1, $2))
          ), 4326))'
  LANGUAGE 'SQL' IMMUTABLE STRICT;

Crazy stuff. Geographies are being cast to geometry, transformed into a planar projection (magically selected by _ST_BestSRID) intersected, then transformed back into geographics and cast back to geography.

The problem is that _ST_BestSRID is failing to pick an appropriate projection to do the calculation in. Now, there are some shapes that are just too big to pick any projection for, but in this case the shapes are a good deal smaller than a hemisphere. However, they are inconveniently located:

One of them goes all the way up to the pole, but the other extends moderately far south. The “best SRID” routine evaluates first whether the shapes can fit into a single UTM zone (they can’t), if they are sufficiently far north to use a polar stereographic (only one is), and finally falls back to a world mercator projection. But, one of the shapes goes up to the pole (which mercator places at infinity), so in the end even the mercator projection fails.

What to do? In order to accurately derive a new shape in planar space that reflects the result on the sphere, we need a projection that represents great circles (used to connect points on the sphere) as straight lines (used to connect points on the plane). We need a Gnomic projection.

Fortunately the projection library underneath PostGIS — Proj4 — supports the gnomic projection. Now, we need a gnomic projection that is parameterized to optimally fit our geography shapes. The approximate center of the two shapes falls at (-160, 70) so we’ll use that as the center of the projection.

Now we insert a record into the PostGIS “spatial_ref_sys” table, giving our new projection an identifier number of 55000:

INSERT INTO spatial_ref_sys (srid, proj4text)
VALUES (55000, '+proj=gnom +lat_0=70 +lon_0=-160');

Now we can run the intersection again, this time in our specially chosen gnomic:

WITH shapes AS
(
SELECT
  ST_GeogFromText('POLYGON((170 50,170 72,-130 72,-130 50,170 50))') AS p1,
  ST_GeogFromText('POLYGON((-170 68,-170 90,-141 90,-141 68,-170 68))') AS p2
)
SELECT ST_AsText(
  geography(ST_Transform(
    ST_Intersection(
      ST_Transform(geometry(p1), 55000),
      ST_Transform(geometry(p2), 55000)
  ),4326)))
FROM shapes;

And this time we get a result:

POLYGON((-169.908 74.035,-141.155 73.424,-141 68,-170 68,-169.908 74.035))

As you can see eyeballing the globe above the resultant shape should have four corners (check), the southern two corners should be the southern two corners of the second input shape (check), and the northern two corners should be derived from the intersection of the northern edge of the first shape and the eastern/western edges of the second (check).

Astute readers of the geometry text will note that the northern edge of the first shape runs from (170 72) to (-130 72): why doesn’t the north edge of the derived shape fall on the 72d parallel? Because the shortest path between those points doesn’t run along the parallel, it’s a great circle arc that stretches to the north along its path.

This example has made it clear to me that the “best SRID” code needs a little upgrade, so that the logic first checks if UTM is acceptable, then switches into finding a best gnomic, and only breaks down and tries mercator if the shape is drastically too large to fit in a hemisphere.

“Best SRID” finding  will always remain a hack, and direct math on the sphere is preferable, but for now we only have a few operations implemented on the sphere (area, length, distance) so we need some hacks to get our work done.