Doesn't display feature if another (possibly invalid) feature is present

[semprebon]

In the attached file, I have two shape files. The "fail" dataset contains Angola, Ziare, and Antartica, and only Ziare is drawn on the map. The "works" dataset is the same, but Antarctica has been removed. Now both Ziare and Angola appear. Both should appear when the first shape file is drawn.

I'm running version 0.6.0 on OS X.

[semprebon]

Test files for generating problem

[springmeyer]

hey semprebon,

Yes, thanks for the report as this has troubled a number of users.

This issue should only occur when attempting to reproject data when the layers dataset contains coordinates out of the bounds the projection can handle, which happens in this case because mercator does not accept coordinates > 90 which occur in antarctica.

So, ultimately the problem is the data combined with the nature of the projection.

This issue is also described here: http://trac.osgeo.org/gdal/ticket/2588 (with that same dataset)

See also: http://www.mail-archive.com/mapnik-users@lists.berlios.de/msg00373.html

I think it may be worthwhile to add a layer option to do the equivalent of what gdal does - essentially skipping geometries that fail to reproject, or otherwise notifying the user that something is wrong in their input data.

-dane

[semprebon]

Patch to handle unprojectable extents in layers

[semprebon]

Yeah, the problem seems to occur when mapnik is building the bounding box in feature_style_processor.hpp, in that it is basically ignoring the return value from proj_trans.forward or proj_trans.backward. When the projection fails on the layer's envelope, it ends up with a completely bogus bounding box.

By checking the return values and using the corresponding map extent point if the projection fails, I at least got my original problem resolved, although it's still not correctly handling the test case I gave you.

[semprebon]

Here is another patch to coord_transform2 that filters out unprojectable points.

[semprebon]

Filters out unprojectable points

[springmeyer]

Replying to semprebon:

Yeah, the problem seems to occur when mapnik is building the bounding box in feature_style_processor.hpp, in that it is basically ignoring the return value from proj_trans.forward or proj_trans.backward. When the projection fails on the layer's envelope, it ends up with a completely bogus bounding box.

Really? I feel like the bounding box of the map has always been generally right. So if a layer's got messed up data/reprojection and shapes fail the bounding box of the map still ends up including at least what was able to be parsed of that shape. I've not look terribly closely however, so perhaps you mean that the messed up bbox at the layer level is why shapes are being aborted altogether. I think that ideal behavior is to skip rendering of features that have invalid coordinates but to expand the map extent by their set of valid coordinates. Does this match with your goals?

By checking the return values and using the corresponding map extent point if the projection fails, I at least got my original problem resolved, although it's still not correctly handling the test case I gave you.

What was your original problem?

[springmeyer]

Okay, I've tried your patches together. They nicely solve your testcase in the way I am imagining such that both angola and zaire show up while antarctica does not while the full map extent (zoom_all()) seems to include antarctica.

However, there was a typo that I fixed before recompiling:

         {
-             command == SEG_MOVETO;
+             command = SEG_MOVETO;
         }

However this approach is not working for another testcase here:

http://mapnik-utils.googlecode.com/svn/example_code/world_population/

The result of this command shows that the combined layer extents (nik2img by default calls m.zoom_all()) are not working:

nik2img.py -m population.xml -o messed_up_bbox.png -p epsg:900913

And forcing a larger bounding box reveals that antarctica as a feature is not being skipped and therefore is rendered oddly:

nik2img.py -m population.xml -o borked_shape.png -p epsg:900913 -e -179,-70,179.0000019073486,70.57026863098145

[springmeyer]

Marking this as decision_decision_needed as we need to take a closer look at this.

Pushing to 0.6.2 for now.

[springmeyer]

Just noticed that the first few shapefiles I tried from the newly released http://www.naturalearthdata.com/ are going to suffer enormously from this. e.g they contain coords outside the range of what mercator supports (and proj_transform can handle) and we need a way to drop just those coordinates without loosing the features themselves, otherwise rendering from files like this is going to take a lot of pre-processing, and headache.

While I was not satisfied with these existing patches, it is clear we need a solution that is is able to keep coords that are within the bounds of the map while avoiding those that are invalid, but which does not visibly corrupt shapes or slow mapnik down noticeably.