Opened 20 years ago
Last modified 14 years ago
#647 closed defect
Reprojection for Lat/Long Wrong — at Initial Version
Reported by: | warmerdam | Owned by: | sdlime |
---|---|---|---|
Priority: | normal | Milestone: | FUTURE |
Component: | WCS Server | Version: | unspecified |
Severity: | normal | Keywords: | |
Cc: |
Description
I have a lat/long layer with the following definition: LAYER NAME GTOPO TYPE raster STATUS ON DATA /usr2/data/jpeg2000/gtopo.jp2 DUMP TRUE # OFFSITE 255 255 255 PROJECTION "init=epsg:4326" END METADATA "wcs_title" "GTOPO" "wcs_label" "GTOPO - label" "wcs_extent" "-100 -10 -60 40" "wcs_resolution" "0.0083333 0.0083333" END END When I get the capabilities from the server, it returns the wrong latlong region: <ContentMetadata> <CoverageOfferingBrief> <name>GTOPO</name> <label>GTOPO - label</label> <lonLatEnvelope srsName="WGS84(DD)"> <gml:pos>-178.682 -572.958</gml:pos> <gml:pos>178.479 2291.83</gml:pos> </lonLatEnvelope> </CoverageOfferingBrief> </ContentMetadata> The problem seems to be that the following code at the end of msWCSGetCoverageMetadata() calls msProjectRect() with a latlong projection as input and NULL as output but down in msProjectPoint() it doesn't know how to reproject to/from a latlong coordinate system. msInitProjection(&proj); // or bad things happen sprintf(projstring, "init=epsg:%.20s", cm->srs+5); if (msLoadProjectionString(&proj, projstring) != 0) return MS_FAILURE; msProjectRect(&proj, NULL, &(cm->llextent)); The mode in msProjectPoint() where the input or output coordinate system can be NULL to mean lat/long does not seem to work properly if the other coordinate system is also lat/long. In particular the following logic fails to apply the degree to radian conversion on the input but it is applied on the output. else { p.u = point->x; p.v = point->y; if(in==NULL || in->proj==NULL) { /* input coordinates are lat/lon */ p.u *= DEG_TO_RAD; /* convert to radians */ p.v *= DEG_TO_RAD; p = pj_fwd(p, out->proj); } else { if(out==NULL || out->proj==NULL) { /* output coordinates are lat/lon */ p = pj_inv(p, in->proj); p.u *= RAD_TO_DEG; /* convert to decimal degrees */ p.v *= RAD_TO_DEG; } else { /* need to go from one projection to another */ p = pj_inv(p, in->proj); p = pj_fwd(p, out->proj); } } point->x = p.u; point->y = p.v; }
Note:
See TracTickets
for help on using tickets.