Problems with mapserver and oblique projections

[Steve Roberts]

Hello,

I've been experiencing this same problem using a standard Polar 
Stereographic projection. The problem has not gone away with the latest 
release(4.10.0). Below is a url to an example mapserver site where I am 
trying to plot ship tracks and coring sites that cross the North Pole:

http://mapserver.eol.ucar.edu/cgi-bin/bug/mapserv-4.10.0?mode=browse&map=/export/web/mapserver/data/bug/test.map&imgxy=300+300&imgext=-2000000+-2000000+2000000+2000000&zoomsize=2

I'm using mapserver 4.10.0 in this example. As you can see not much shows 
up with the initial map extent which is centered on the North Pole. But 
when you recenter the map away from the pole more of the ship tracks and 
core sites will show up.

As Frank mentioned the core of the problem is based on the false 
assumption that the lat/long equivalent of the edge of the map bounding 
box can be used to represent a rectangle bounding box to filter lat/long 
data that is outside of your map view. In the above example url with the 
map centered on the North Pole the edges of the bounding box never gets 
above 71.5 degrees so we end up with a circular filter "box" of 
180W,180E,64N,71.5N where everything outside of this "box" is not plotted. 
This is obviously not what we want. There is an excellent article on the 
web which discusses these issues at:

  http://chukchi.colorado.edu/PAPER/

There was mention of MapServer "sampling" through the destination 
rectangle if there is a failure to reproject any of the edge points. In 
this example and many others this never occurs so this does not solve the 
problem. And doing a course sampling through the interior of the rectangle 
is not an optimal solution either. There is a high probability of missing 
the pole. In my example url above I have a coring site right at the North 
Pole and when I modify the code to force "sampling" most of the data will 
now show up but the North Pole core site still fails to plot. A proposed 
simple solution would be to add some code to detect if the poles are 
inside the map bounding box and adjust the extents accordingly. Adding the 
following code to the subroutine msProjectRect in mapproject.c just before 
the line "if( failure > 0 )" would do the trick:


   /* Check if North or South Pole is in search rectangle */
   if(!pj_is_latlong(in->proj) && pj_is_latlong(out->proj) ) {
       /* Check if North Pole is in search rectangle */
       prj_point.x = 0;
       prj_point.y = 90;
       if( msProjectPoint( out, in, &prj_point ) == MS_FAILURE ) {
           msDebug( "msProjectRect(): Failed to project North Pole.\n" );
           pole_failure++;
       } else {;
           if(msPointInRect(&prj_point,rect)) {
              prj_rect.maxy = 90;
           }
       }

       /* Check if South Pole is in search rectangle */
       prj_point.x = 0;
       prj_point.y = -90;
       if( msProjectPoint( out, in, &prj_point ) == MS_FAILURE ) {
           msDebug( "msProjectRect(): Failed to project South Pole.\n" );
           pole_failure++;
       } else {;
           if(msPointInRect(&prj_point,rect)) {
              prj_rect.miny = -90;
           }
       }
       if(pole_failure > 1) {
           msDebug( "msProjectRect(): Failed to project both North and South Pole.\n" );
           failure++;
       }
   }


And also add:

   int     pole_failure=0;

at the beginning of the subroutine. The above code works by projecting the 
north and south poles from lat/longs to the map projection and then using 
the routine "msPointInRect" to check if these points are inside the map 
bounding box. With the above addition the above example will now generate 
the bounds 180W,180E,64N,90N and all the data should now show up.

However, I've noticed that msProjectRect still has issues when the 
lat/long bounding box crosses the dateline. Regardless of how small the 
bounding box is, if it crosses the dateline it will produce a longitudinal 
extent of -180 to 180. Using the terminology of the above referenced 
article this will potentially generate a large number of "False 
Positives". Normally the user would not be aware of this but if you are 
dealing with a large dataset I would imagine this could potentially slow 
things down significantly. I'm not sure how to deal with this but just 
thought I would point this out.

-Steve
-- 
Steve Roberts
NCAR/EOL
P.O. Box 3000 
Boulder, CO 80307-3000		eol:  http://www.eol.ucar.edu