2

when I use:

gdalwarp -dstalpha -multi -t_srs epsg:900913 b172.tif b282.tif b770.tif b771.tif b772.tif b832.tif w_25k.tif

I can get one merged tif then go through gdal2tiles for tiling

But if I use gdalwarp it takes really long time to process and very slow.

So I use gdalbuildvrt to create vrt first, then use gdaltranslate to get merged tif.

gdalbuildvrt -overwrite -a_srs epsg:900913 b172.tif b282.tif b770.tif b771.tif b772.tif b832.tif test.vrt

gdal_translate -of GTiff -a_srs epsg:900913 test.vrt  25k.tif

So I get test.tif then go through gdal2tiles the result tif is on top of Africa not the the destination of my tif (Australia).

I use:

gdalinfo 25k.tif

Driver: GTiff/GeoTIFF
Files: 25k.tif
Size is 27063, 13066
Coordinate System is:
PROJCS["Google Maps Global Mercator",
    GEOGCS["WGS 84",
        DATUM["WGS_1984",
            SPHEROID["WGS 84",6378137,298.257223563,
                AUTHORITY["EPSG","7030"]],
            AUTHORITY["EPSG","6326"]],
        PRIMEM["Greenwich",0,
            AUTHORITY["EPSG","8901"]],
        UNIT["degree",0.01745329251994328,
            AUTHORITY["EPSG","9122"]],
        AUTHORITY["EPSG","4326"]],
    PROJECTION["Mercator_2SP"],
    PARAMETER["standard_parallel_1",0],
    PARAMETER["latitude_of_origin",0],
    PARAMETER["central_meridian",0],
    PARAMETER["false_easting",0],
    PARAMETER["false_northing",0],
    UNIT["Meter",1],
    EXTENSION["PROJ4","+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +wktext  +no_defs"],
    AUTHORITY["EPSG","900913"]]
Origin = (2511424.026000000100000,2339504.385999999900000)
Pixel Size = (0.846642320907504,-0.846706498927300)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left  ( 2511424.026, 2339504.386) ( 22d33'37.82"E, 20d33'36.44"N)
Lower Left  ( 2511424.026, 2328441.319) ( 22d33'37.82"E, 20d28' 1.35"N)
Upper Right ( 2534336.707, 2339504.386) ( 22d45'58.80"E, 20d33'36.44"N)
Lower Right ( 2534336.707, 2328441.319) ( 22d45'58.80"E, 20d28' 1.35"N)
Center      ( 2522880.367, 2333972.852) ( 22d39'48.31"E, 20d30'48.92"N)
Band 1 Block=27063x1 Type=Byte, ColorInterp=Red
  Mask Flags: PER_DATASET ALPHA
Band 2 Block=27063x1 Type=Byte, ColorInterp=Green
  Mask Flags: PER_DATASET ALPHA
Band 3 Block=27063x1 Type=Byte, ColorInterp=Blue
  Mask Flags: PER_DATASET ALPHA
Band 4 Block=27063x1 Type=Byte, ColorInterp=Alpha

gdalinfo w_25k.tif

Driver: GTiff/GeoTIFF
Files: w_25k.tif
Size is 27080, 13104
Coordinate System is:
PROJCS["Google Maps Global Mercator",
    GEOGCS["WGS 84",
        DATUM["WGS_1984",
            SPHEROID["WGS 84",6378137,298.257223563,
                AUTHORITY["EPSG","7030"]],
            AUTHORITY["EPSG","6326"]],
        PRIMEM["Greenwich",0,
            AUTHORITY["EPSG","8901"]],
        UNIT["degree",0.01745329251994328,
            AUTHORITY["EPSG","9122"]],
        AUTHORITY["EPSG","4326"]],
    PROJECTION["Mercator_2SP"],
    PARAMETER["standard_parallel_1",0],
    PARAMETER["latitude_of_origin",0],
    PARAMETER["central_meridian",0],
    PARAMETER["false_easting",0],
    PARAMETER["false_northing",0],
    UNIT["Meter",1],
    EXTENSION["PROJ4","+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +wktext  +no_defs"],
    AUTHORITY["EPSG","900913"]]
Origin = (16155902.033114368000000,-4642590.182159157500000)
Pixel Size = (1.080346179541381,-1.080346179541381)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left  (16155902.033,-4642590.182) (145d 7'51.37"E, 38d26'44.75"S)
Lower Left  (16155902.033,-4656747.038) (145d 7'51.37"E, 38d32'43.07"S)
Upper Right (16185157.808,-4642590.182) (145d23'37.49"E, 38d26'44.75"S)
Lower Right (16185157.808,-4656747.038) (145d23'37.49"E, 38d32'43.07"S)
Center      (16170529.920,-4649668.610) (145d15'44.43"E, 38d29'43.97"S)
Band 1 Block=27080x1 Type=Byte, ColorInterp=Red
  Mask Flags: PER_DATASET ALPHA
Band 2 Block=27080x1 Type=Byte, ColorInterp=Green
  Mask Flags: PER_DATASET ALPHA
Band 3 Block=27080x1 Type=Byte, ColorInterp=Blue
  Mask Flags: PER_DATASET ALPHA
Band 4 Block=27080x1 Type=Byte, ColorInterp=Alpha

I found the Corner Coordinates are changes can some one help me to keep the origin Corner Coordinates, please?

enter image description here

The origin Data has *.ers file

DatasetHeader Begin
Version     = "7.1"
Name        = "m734wh.ers"
SourceDataset   = "m734wh.tif"
LastUpdated = Thu May 15 16:07:08 GMT 2014
DataFile    = "m734wh.tif"
DataSetType = Translated
DataType    = Raster
ByteOrder   = LSBFirst
CoordinateSpace Begin
    Datum       = "GDA94"
    Projection  = "LOCAL"
    CoordinateType  = EN
    Rotation    = 0:0:0.0
CoordinateSpace End
RasterInfo Begin
    CellType    = Unsigned8BitInteger
    CellInfo Begin
        Xdimension  = 0.8465509849362
        Ydimension  = 0.846666666666667
    CellInfo End
    NrOfLines   = 1200
    NrOfCellsPerLine    = 1726
    RegistrationCoord Begin
        Eastings    = 2519823.203
        Northings   = 2336158.68
    RegistrationCoord End
    NrOfBands   = 3
    BandId Begin
        Value       = "Red Layer"
    BandId End
    BandId Begin
        Value       = "Green Layer"
    BandId End
    BandId Begin
        Value       = "Blue Layer"
    BandId End
    RegionInfo Begin
        Type        = Polygon
        RegionName  = "All"
        SubRegion   = {
            0   0
            0   1200
            1726    1200
            1726    0
            }
    RegionInfo End
RasterInfo End

DatasetHeader End

3
  • At least using -a_srs in the latter command is not needed because GDAL can read the projection from the vrt file. It should not make problems either but make still a new trial.
    – user30184
    Commented Jul 26, 2016 at 7:13
  • 1
    Are your input tiffs really in EPSG:3857 (900913)?
    – user30184
    Commented Jul 26, 2016 at 7:21
  • @user30184 I have update the question with origin Datum
    – Daqian
    Commented Jul 26, 2016 at 7:31

1 Answer 1

2

Your source data is in GDA94, supposed to be Vicgrid94 covering Philip Island. gdalwarp transforms correctly to a target CRS of epsg:3857, but gdal_translate will not do that.

-a_srs has a different purpose: just assign a CRS, no reprojection.

So you have to use:

gdalwarp -of GTiff -s_srs epsg:3111 -a_srs epsg:3857 test.vrt  25k.tif

BTW I used this batch to merge and tile a set of maps in a local projection:

for %%N in (D:\Karten\gdal\gdal2tiles\NL25\*.tif) DO gdal_translate -of vrt -expand rgba %%N D:\Karten\gdal\gdal2tiles\NL25\%%~nN.vrt
gdalbuildvrt -allow_projection_difference index25.vrt NL25\*.vrt
gdal2tiles --s_srs EPSG:28992 --zoom 15-16 index25.vrt tiles
pause

No big tif file needed, and gdal2tiles reprojects to EPSG:3857 itself. The -expand rgba is needed for paletted source files, so I wrapped those in a vrt for each source file before merging. That might fix your issue in GDAL merge warp has black block can't remove as well, since it adds an alpha channel too.

2
  • Thanks sooo much @AndreJ, you're life saver. I have been go through so much to find the answer, I wish I could met early. But a small question, if you don't merge the big tif how do you tile it, if you do gdal2tile on one tif a time it will overwrite the previous images. Can you please give a hint on that???
    – Daqian
    Commented Jul 27, 2016 at 1:24
  • I don't tile single tifs, but the merged index25.vrt file. Contains everything like a big tif, but does not need much disk space.
    – AndreJ
    Commented Jul 27, 2016 at 5:50

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.