Skip to content

Commit

Permalink
UVRaster: add support for query mode / WMS GetFeatureInfo
Browse files Browse the repository at this point in the history
  • Loading branch information
@rouault
rouault committed Jan 6, 2020
1 parent 91782df commit 70aeb86
Show file tree
Hide file tree
Showing 5 changed files with 165 additions and 94 deletions.
89 changes: 2 additions & 87 deletions mapdraw.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -974,94 +974,9 @@ int msDrawVectorLayer(mapObj *map, layerObj *layer, imageObj *image)
/* original area of interest is (map->extent, map->projection)... */
/* Useful when dealin with UVRASTER that extend beyond 180 deg */
msUVRASTERLayerUseMapExtentAndProjectionForNextWhichShapes( layer, map );
}

/* For UVRaster, it is important that the searchrect is not too large */
/* to avoid insufficient intermediate raster resolution, which could */
/* happen if we use the default code path, given potential reprojection */
/* issues when using a map extent that is not in the validity area of */
/* the layer projection. */
if( layer->connectiontype == MS_UVRASTER &&
!layer->projection.gt.need_geotransform &&
!(msProjIsGeographicCRS(&(map->projection)) &&
msProjIsGeographicCRS(&(layer->projection))) ) {
rectObj layer_ori_extent;

if( msLayerGetExtent(layer, &layer_ori_extent) == MS_SUCCESS ) {
projectionObj map_proj;

double map_extent_minx = map->extent.minx;
double map_extent_miny = map->extent.miny;
double map_extent_maxx = map->extent.maxx;
double map_extent_maxy = map->extent.maxy;
rectObj layer_extent = layer_ori_extent;

/* Create a variant of map->projection without geotransform for */
/* conveniency */
msInitProjection(&map_proj);
msCopyProjection(&map_proj, &map->projection);
map_proj.gt.need_geotransform = MS_FALSE;
if( map->projection.gt.need_geotransform ) {
map_extent_minx = map->projection.gt.geotransform[0]
+ map->projection.gt.geotransform[1] * map->extent.minx
+ map->projection.gt.geotransform[2] * map->extent.miny;
map_extent_miny = map->projection.gt.geotransform[3]
+ map->projection.gt.geotransform[4] * map->extent.minx
+ map->projection.gt.geotransform[5] * map->extent.miny;
map_extent_maxx = map->projection.gt.geotransform[0]
+ map->projection.gt.geotransform[1] * map->extent.maxx
+ map->projection.gt.geotransform[2] * map->extent.maxy;
map_extent_maxy = map->projection.gt.geotransform[3]
+ map->projection.gt.geotransform[4] * map->extent.maxx
+ map->projection.gt.geotransform[5] * map->extent.maxy;
}

/* Reproject layer extent to map projection */
msProjectRect(&layer->projection, &map_proj, &layer_extent);

if( layer_extent.minx <= map_extent_minx &&
layer_extent.miny <= map_extent_miny &&
layer_extent.maxx >= map_extent_maxx &&
layer_extent.maxy >= map_extent_maxy ) {
/* do nothing special if area to map is inside layer extent */
}
else {
if( layer_extent.minx >= map_extent_minx &&
layer_extent.maxx <= map_extent_maxx &&
layer_extent.miny >= map_extent_miny &&
layer_extent.maxy <= map_extent_maxy ) {
/* if the area to map is larger than the layer extent, then */
/* use full layer extent and add some margin to reflect the */
/* proportion of the useful area over the requested bbox */
double extra_x =
(map_extent_maxx - map_extent_minx) /
(layer_extent.maxx - layer_extent.minx) *
(layer_ori_extent.maxx - layer_ori_extent.minx);
double extra_y =
(map_extent_maxy - map_extent_miny) /
(layer_extent.maxy - layer_extent.miny) *
(layer_ori_extent.maxy - layer_ori_extent.miny);
searchrect.minx = layer_ori_extent.minx - extra_x / 2;
searchrect.maxx = layer_ori_extent.maxx + extra_x / 2;
searchrect.miny = layer_ori_extent.miny - extra_y / 2;
searchrect.maxy = layer_ori_extent.maxy + extra_y / 2;
}
else
{
/* otherwise clip the map extent with the reprojected layer */
/* extent */
searchrect.minx = MS_MAX( map_extent_minx, layer_extent.minx );
searchrect.maxx = MS_MIN( map_extent_maxx, layer_extent.maxx );
searchrect.miny = MS_MAX( map_extent_miny, layer_extent.miny );
searchrect.maxy = MS_MIN( map_extent_maxy, layer_extent.maxy );
/* and reproject into the layer projection */
msProjectRect(&map_proj, &layer->projection, &searchrect);
}
bDone = MS_TRUE;
}

msFreeProjection(&map_proj);
}
searchrect = msUVRASTERGetSearchRect( layer, map );
bDone = MS_TRUE;
}

if( !bDone )
Expand Down
1 change: 1 addition & 0 deletions mapserver.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2534,6 +2534,7 @@ void msPopulateTextSymbolForLabelAndString(textSymbolObj *ts, labelObj *l, char

int LayerDefaultGetShapeCount(layerObj *layer, rectObj rect, projectionObj *rectProjection);
void msUVRASTERLayerUseMapExtentAndProjectionForNextWhichShapes(layerObj* layer, mapObj* map);
rectObj msUVRASTERGetSearchRect( layerObj* layer, mapObj* map );

/* ==================================================================== */
/* Prototypes for functions in mapdraw.c */
Expand Down
136 changes: 129 additions & 7 deletions mapuvraster.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -31,6 +31,7 @@

#include <assert.h>
#include <math.h>
#include "mapows.h"
#include "mapresample.h"
#include "mapthread.h"

Expand Down Expand Up @@ -175,6 +176,28 @@ static void msUVRasterLayerInfoInitialize(layerObj *layer)
uvlinfo->width = 0;
uvlinfo->height = 0;

/* Set attribute type to Real, unless the user has explicitly set */
/* something else. */
{
const char* const items[] = {
MSUVRASTER_ANGLE,
MSUVRASTER_MINUS_ANGLE,
MSUVRASTER_LENGTH,
MSUVRASTER_LENGTH_2,
MSUVRASTER_U,
MSUVRASTER_V,
};
size_t i;
for( i = 0; i < sizeof(items)/sizeof(items[0]); ++i ) {
char szTmp[100];
snprintf(szTmp, sizeof(szTmp), "%s_type", items[i]);
if (msOWSLookupMetadata(&(layer->metadata), "OFG", szTmp) == NULL) {
snprintf(szTmp, sizeof(szTmp), "gml_%s_type", items[i]);
msInsertHashTable(&(layer->metadata), szTmp, "Real");
}
}
}

/* uvlinfo->query_result_hard_max = 1000000; */

/* if( CSLFetchNameValue( layer->processing, "RASTER_QUERY_MAX_RESULT" ) */
Expand Down Expand Up @@ -342,6 +365,104 @@ static char **msUVRASTERGetValues(layerObj *layer, float *u, float *v)
return values;
}

rectObj msUVRASTERGetSearchRect( layerObj* layer, mapObj* map )
{
rectObj searchrect = map->extent;
int bDone = MS_FALSE;

/* For UVRaster, it is important that the searchrect is not too large */
/* to avoid insufficient intermediate raster resolution, which could */
/* happen if we use the default code path, given potential reprojection */
/* issues when using a map extent that is not in the validity area of */
/* the layer projection. */
if( !layer->projection.gt.need_geotransform &&
!(msProjIsGeographicCRS(&(map->projection)) &&
msProjIsGeographicCRS(&(layer->projection))) ) {
rectObj layer_ori_extent;

if( msLayerGetExtent(layer, &layer_ori_extent) == MS_SUCCESS ) {
projectionObj map_proj;

double map_extent_minx = map->extent.minx;
double map_extent_miny = map->extent.miny;
double map_extent_maxx = map->extent.maxx;
double map_extent_maxy = map->extent.maxy;
rectObj layer_extent = layer_ori_extent;

/* Create a variant of map->projection without geotransform for */
/* conveniency */
msInitProjection(&map_proj);
msCopyProjection(&map_proj, &map->projection);
map_proj.gt.need_geotransform = MS_FALSE;
if( map->projection.gt.need_geotransform ) {
map_extent_minx = map->projection.gt.geotransform[0]
+ map->projection.gt.geotransform[1] * map->extent.minx
+ map->projection.gt.geotransform[2] * map->extent.miny;
map_extent_miny = map->projection.gt.geotransform[3]
+ map->projection.gt.geotransform[4] * map->extent.minx
+ map->projection.gt.geotransform[5] * map->extent.miny;
map_extent_maxx = map->projection.gt.geotransform[0]
+ map->projection.gt.geotransform[1] * map->extent.maxx
+ map->projection.gt.geotransform[2] * map->extent.maxy;
map_extent_maxy = map->projection.gt.geotransform[3]
+ map->projection.gt.geotransform[4] * map->extent.maxx
+ map->projection.gt.geotransform[5] * map->extent.maxy;
}

/* Reproject layer extent to map projection */
msProjectRect(&layer->projection, &map_proj, &layer_extent);

if( layer_extent.minx <= map_extent_minx &&
layer_extent.miny <= map_extent_miny &&
layer_extent.maxx >= map_extent_maxx &&
layer_extent.maxy >= map_extent_maxy ) {
/* do nothing special if area to map is inside layer extent */
}
else {
if( layer_extent.minx >= map_extent_minx &&
layer_extent.maxx <= map_extent_maxx &&
layer_extent.miny >= map_extent_miny &&
layer_extent.maxy <= map_extent_maxy ) {
/* if the area to map is larger than the layer extent, then */
/* use full layer extent and add some margin to reflect the */
/* proportion of the useful area over the requested bbox */
double extra_x =
(map_extent_maxx - map_extent_minx) /
(layer_extent.maxx - layer_extent.minx) *
(layer_ori_extent.maxx - layer_ori_extent.minx);
double extra_y =
(map_extent_maxy - map_extent_miny) /
(layer_extent.maxy - layer_extent.miny) *
(layer_ori_extent.maxy - layer_ori_extent.miny);
searchrect.minx = layer_ori_extent.minx - extra_x / 2;
searchrect.maxx = layer_ori_extent.maxx + extra_x / 2;
searchrect.miny = layer_ori_extent.miny - extra_y / 2;
searchrect.maxy = layer_ori_extent.maxy + extra_y / 2;
}
else
{
/* otherwise clip the map extent with the reprojected layer */
/* extent */
searchrect.minx = MS_MAX( map_extent_minx, layer_extent.minx );
searchrect.maxx = MS_MIN( map_extent_maxx, layer_extent.maxx );
searchrect.miny = MS_MAX( map_extent_miny, layer_extent.miny );
searchrect.maxy = MS_MIN( map_extent_maxy, layer_extent.maxy );
/* and reproject into the layer projection */
msProjectRect(&map_proj, &layer->projection, &searchrect);
}
bDone = MS_TRUE;
}

msFreeProjection(&map_proj);
}
}

if( !bDone )
msProjectRect(&map->projection, &layer->projection, &searchrect); /* project the searchrect to source coords */

return searchrect;
}

int msUVRASTERLayerWhichShapes(layerObj *layer, rectObj rect, int isQuery)
{
uvRasterLayerInfo *uvlinfo = (uvRasterLayerInfo *) layer->layerinfo;
Expand All @@ -366,13 +487,6 @@ int msUVRASTERLayerWhichShapes(layerObj *layer, rectObj rect, int isQuery)
if( uvlinfo == NULL )
return MS_FAILURE;

/* QUERY NOT SUPPORTED YET */
if (isQuery == MS_TRUE) {
msSetError( MS_MISCERR, "Query is not supported for UV layer.", "msUVRASTERLayerWhichShapes()" );
return MS_FAILURE;
}


if( CSLFetchNameValue( layer->processing, "BANDS" ) == NULL ) {
msSetError( MS_MISCERR, "BANDS processing option is required for UV layer. You have to specified 2 bands.",
"msUVRASTERLayerWhichShapes()" );
Expand Down Expand Up @@ -563,6 +677,12 @@ int msUVRASTERLayerWhichShapes(layerObj *layer, rectObj rect, int isQuery)
}
}

if( isQuery )
{
/* For query mode, use layer->map->extent reprojected rather than */
/* the provided rect. Generic query code will filter returned features. */
rect = msUVRASTERGetSearchRect(layer, layer->map);
}

map_cellsize = MS_MAX(MS_CELLSIZE(rect.minx, rect.maxx,layer->map->width),
MS_CELLSIZE(rect.miny,rect.maxy,layer->map->height));
Expand Down Expand Up @@ -758,6 +878,8 @@ int msUVRASTERLayerGetShape(layerObj *layer, shapeObj *shape, resultObj *record)

shape->numvalues = layer->numitems;
shape->values = msUVRASTERGetValues(layer, &uvlinfo->u[x][y], &uvlinfo->v[x][y]);
shape->index = shapeindex;
shape->resultindex = shapeindex;

return MS_SUCCESS;

Expand Down
23 changes: 23 additions & 0 deletions msautotest/wxs/expected/wms_uvraster_map_reprojection_getfeatureinfo.xml
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,23 @@
<?xml version="1.0" encoding="UTF-8"?>

<msGMLOutput
xmlns:gml="http://www.opengis.net/gml"
xmlns:xlink="http://www.w3.org/1999/xlink"
xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance">
<test_layer>
<gml:name>test</gml:name>
<test_feature>
<gml:boundedBy>
<gml:Box srsName="EPSG:4326">
<gml:coordinates>82.689740,80.025912 82.689740,80.025912</gml:coordinates>
</gml:Box>
</gml:boundedBy>
<uv_angle>0.000000</uv_angle>
<uv_minus_angle>180.000000</uv_minus_angle>
<uv_length>2.378277</uv_length>
<uv_length_2>1.189139</uv_length_2>
<u>2.378277</u>
<v>0.000000</v>
</test_feature>
</test_layer>
</msGMLOutput>
10 changes: 10 additions & 0 deletions msautotest/wxs/wms_uvraster_map_reprojection.map
View file
Original file line number Diff line number Diff line change
Expand Up @@ -8,6 +8,8 @@
# RUN_PARMS: wms_uvraster_map_reprojection_extent_smaller_than_layer.png [MAPSERV] QUERY_STRING="map=[MAPFILE]&SERVICE=WMS&VERSION=1.3.0&REQUEST=GetMap&BBOX=30,-170,60,-100&CRS=EPSG:4326&WIDTH=400&HEIGHT=200&LAYERS=test&STYLES=&FORMAT=image/png&TRANSPARENT=TRUE" > [RESULT_DEMIME]
#
# RUN_PARMS: wms_uvraster_map_reprojection_extent_intersecting_layer.png [MAPSERV] QUERY_STRING="map=[MAPFILE]&SERVICE=WMS&VERSION=1.3.0&REQUEST=GetMap&BBOX=30,-200,60,-100&CRS=EPSG:4326&WIDTH=400&HEIGHT=200&LAYERS=test&STYLES=&FORMAT=image/png&TRANSPARENT=TRUE" > [RESULT_DEMIME]
#
# RUN_PARMS: wms_uvraster_map_reprojection_getfeatureinfo.xml [MAPSERV] QUERY_STRING="map=[MAPFILE]&SERVICE=WMS&VERSION=1.3.0&REQUEST=GetFeatureInfo&BBOX=-30.28229908443537965,-68.01500508646998355,123.14801627670399853,120.98499491353001645&CRS=EPSG:4326&WIDTH=983&HEIGHT=798&LAYERS=test&STYLES=&FORMAT=image/png&QUERY_LAYERS=test&INFO_FORMAT=application/vnd.ogc.gml&I=798&J=232&FEATURE_COUNT=10" > [RESULT_DEMIME]

MAP
NAME test
Expand Down Expand Up @@ -39,6 +41,7 @@ MAP
WEB
METADATA
"ows_enable_request" "*"
"wms_onlineresource" "http://localhost:8080/mapserv.cgi?map=/home/even/mapserver/mapserver/msautotest/wxs/wms_uvraster_map_reprojection.map&"
END
END

Expand All @@ -48,6 +51,13 @@ MAP
CONNECTIONTYPE uvraster
PROCESSING "BANDS=1,2"
PROCESSING "UV_SPACING=20"
TEMPLATE "ttt"
TOLERANCE 20

METADATA
"wms_title" "test"
"gml_include_items" "all"
END

PROJECTION
"proj=stere"
Expand Down

0 comments on commit 70aeb86

Please sign in to comment.