[mapserver-users] Use Ground Control Points from GeoTIFF

Rahkonen Jukka (MML) jukka.rahkonen at maanmittauslaitos.fi
Thu Sep 3 04:31:09 PDT 2020


Hi,

I would appreciate if you could have a try with
"gdalwarp -of vrt -t_srs "EPSG:3857" sat.tif sat_3857.vrt"

-Jukka Rahkonen-

Lähettäjä: mapserver-users <mapserver-users-bounces at lists.osgeo.org> Puolesta Jesper Kempe
Lähetetty: torstai 3. syyskuuta 2020 13.57
Vastaanottaja: Jeff McKenna <jmckenna at gatewaygeomatics.com>
Kopio: Mapserver <mapserver-users at lists.osgeo.org>
Aihe: Re: [mapserver-users] Use Ground Control Points from GeoTIFF

shp2img creates an empty blank image with a mapfile like this:

MAP
    NAME "map"
    WEB
        METADATA
            "wms_title"            "WMS mapserver"
            "wms_srs"              "EPSG:4326 EPSG:3021 EPSG:9999 EPSG:3006 EPSG:900913"
            "wms_extent"           "1225000 6130000 1910000 7675000"
            "wms_enable_request"   "*"
            "wfs_enable_request"   "*"
        END
    END

    PROJECTION
        "init=epsg:4326"
    END

    EXTENT 15 60 30 68
    SIZE 200 200

    LAYER
      NAME "sat"
      DATA "/data/temp/Jesper/sat.tif"
      TYPE RASTER
      STATUS DEFAULT

      METADATA
                "wms_title" "sat"
      END
    END


The output from shp2img looks like this:

$ shp2img -m ttt.map -o /data/temp/Jesper/ttt.png -all_debug 5
msLoadMap(): 0.003s
msDrawMap(): rendering using outputformat named png (AGG/PNG).
msDrawRasterLayerLow(sat): entering.
msDrawRasterLayerGDAL(): Entering transform.
msDrawRasterLayerGDAL(): src=15,10645,15,15, dst=0,0,200,200
msDrawRasterLayerGDAL(): source raster PL (14.962,10645.462) for dst PL (0,0).
msDrawRasterLayerGDAL(): red,green,blue,alpha bands = 1,2,3,4
msDrawMap(): Layer 0 (sat), 0.011s
msDrawMap(): Drawing Label Cache, 0.000s
msDrawMap() total time: 0.011s
msSaveImage(/data/temp/Jesper/ttt.png) total time: 0.003s
msFreeMap(): freeing map at 0x1d14f50.
freeLayer(): freeing layer at 0x1d1d080.
freeLayer(): freeing layer at 0x1d1e1b0.
shp2img total time: 0.017s


If I add:
      PROJECTION
        "init=epsg:4326"
      END
      EXTENT 17.0473154408269387 62.1606759790530177 27.5792082876237146 66.7744046276809797

shp2img outputs:

msLoadMap(): 0.003s
msDrawMap(): rendering using outputformat named png (AGG/PNG).
msDrawRasterLayerLow(sat): entering.
msDrawRasterLayerGDAL(): Entering transform.
msDrawRasterLayerGDAL(): src=0,0,10474,10717, dst=27,63,140,61
msDrawRasterLayerGDAL(): source raster PL (-49.554,-33.774) for dst PL (27,63).
msDrawRasterLayerGDAL(): red,green,blue,alpha bands = 1,2,3,4
msDrawMap(): Layer 0 (sat), 0.025s
msDrawMap(): Drawing Label Cache, 0.000s
msDrawMap() total time: 0.025s
msSaveImage(/data/temp/Jesper/ttt.png) total time: 0.004s
msFreeMap(): freeing map at 0x2453f50.
freeLayer(): freeing layer at 0x245c080.
freeLayer(): freeing layer at 0x245d430.
shp2img total time: 0.033s

Then the image is shown in ttt.png but is not in the correct place and not rotated correctly.

Any help or ideas are appreciated!

Thanks
Jesper

Den ons 2 sep. 2020 kl 18:23 skrev Jeff McKenna <jmckenna at gatewaygeomatics.com<mailto:jmckenna at gatewaygeomatics.com>>:
I'm curious what happens when you try to generate a map image from the
commandline: shp2img -m test.map -o ttt.png -all_debug 5

-jeff



--
Jeff McKenna
MapServer Consulting and Training Services
co-founder of FOSS4G
http://gatewaygeo.com/



On 2020-09-02 4:51 a.m., Jesper Kempe wrote:
> Hi
>
> This is the mapfile:
>
> MAP
>      NAME "map"
>      WEB
>          METADATA
>              "wms_title"            "WMS mapserver"
>              "wms_srs"              "EPSG:4326 EPSG:3021 EPSG:9999
> EPSG:3006 EPSG:900913"
>              "wms_extent"           "1225000 6130000 1910000 7675000"
>              "wms_enable_request"   "*"
>              "wfs_enable_request"   "*"
>          END
>      END
>
>      PROJECTION
>          "init=epsg:4326"
>      END
>
>      EXTENT -180 -90 180 90
>
>      LAYER
>        NAME "sat"
>        DATA "sat.tif"
>        TYPE RASTER
>        STATUS ON
>
>        METADATA
>                  "wms_title" "sat"
>        END
>
>        PROJECTION
>          "init=epsg:4326"
>        END
>        EXTENT 17.0473154408269387 62.1606759790530177
> 27.5792082876237146 66.7744046276809797
>
>      END
>
> END
>
> I have also tried both without PROJECTION in LAYER and with
> PROJECTION
>       AUTO
> END
>
> None of them works.
>
> EXTENT in the mapfile is taken from Qgis metadata of the layer.
>
>
> Thanks
> Jesper
>
> Den tis 1 sep. 2020 kl 23:41 skrev Steve Lime <sdlime at gmail.com<mailto:sdlime at gmail.com>
> <mailto:sdlime at gmail.com<mailto:sdlime at gmail.com>>>:
>
>     What does the mapfile look like?
>
>     On Tue, Sep 1, 2020 at 10:26 AM Jesper Kempe <jesper.kempe at gmail.com<mailto:jesper.kempe at gmail.com>
>     <mailto:jesper.kempe at gmail.com<mailto:jesper.kempe at gmail.com>>> wrote:
>
>         Hi
>
>         I have a GeoTiff file that gdalinfo presents like this:
>
>         gdalinfo sat.tif
>         Driver: GTiff/GeoTIFF
>         Files: sat.tif
>         Size is 10474, 10717
>         Coordinate System is `'
>         GCP Projection =
>         GEOGCS["WGS 84",
>              DATUM["WGS_1984",
>                  SPHEROID["WGS 84",6378137,298.257223563,
>                      AUTHORITY["EPSG","7030"]],
>                  AUTHORITY["EPSG","6326"]],
>              PRIMEM["Greenwich",0],
>              UNIT["degree",0.0174532925199433],
>              AUTHORITY["EPSG","4326"]]
>         GCP[  0]: Id=1, Info=
>                    (0,0) ->
>         (27.5783052509812,65.9184735836942,229.992792810313)
>         GCP[  1]: Id=2, Info=
>                    (524,0) ->
>         (27.1416007200269,65.9731480840222,172.994255514815)
>         GCP[  2]: Id=3, Info=
>                    (1048,0) ->
>         (26.7005564512604,66.0268770755367,152.994545626454)
>         GCP[  3]: Id=4, Info=
>                    (1572,0) ->
>         (26.257537545879,66.0793551447287,131.994995356537)
>         GCP[  4]: Id=5, Info=
>                    (2096,0) ->
>         (25.8102107755995,66.1308402121171,163.993212723173)
>         ....
>
>         with 482 GCPs.
>
>
>         The file opens ok in Qgis and is located in the right location.
>         But I can´t in any way get it to show correctly in Mapserver.
>
>         Is this possible to do? Or do I have to gdalwarp the file to a
>         known projection and set the PROJECTION of the layer?
>
>         Thanks
>         Jesper
>
>



_______________________________________________
mapserver-users mailing list
mapserver-users at lists.osgeo.org<mailto:mapserver-users at lists.osgeo.org>
https://lists.osgeo.org/mailman/listinfo/mapserver-users
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/mapserver-users/attachments/20200903/9ea97930/attachment.htm>


More information about the MapServer-users mailing list