Reprojecting from MODIS Sinusoidal to Web Mercator: how to calculate output raster dimensions?

Question

I have a MODIS image from AWS in sinusoidal projection and I want to calculate output raster size after reprojecting it into Web Mercator (EPSG:3857).

$ gdalinfo MCD43A4.A2017185.h11v03.006.2017194033101_B01.TIF Driver: GTiff/GeoTIFF
Files: MCD43A4.A2017185.h11v03.006.2017194033101_B01.TIF
Size is 2400, 2400
Coordinate System is:
PROJCS["unnamed",
    GEOGCS["Unknown datum based upon the custom spheroid",
        DATUM["Not_specified_based_on_custom_spheroid",
            SPHEROID["Custom spheroid",6371007.181,0]],
        PRIMEM["Greenwich",0],
        UNIT["degree",0.0174532925199433]],
    PROJECTION["Sinusoidal"],
    PARAMETER["longitude_of_center",0],
    PARAMETER["false_easting",0],
    PARAMETER["false_northing",0],
    UNIT["metre",1,
        AUTHORITY["EPSG","9001"]]]
Origin = (-7783653.637667000293732,6671703.117999999783933)
Pixel Size = (463.312716527916905,-463.312716527916507)
Metadata:
  add_offset=0
  add_offset_err=0
  ALBEDOFILEID=06121997
  AREA_OR_POINT=Area
  ASSOCIATEDINSTRUMENTSHORTNAME=MODIS
  ASSOCIATEDPLATFORMSHORTNAME=Terra
  ASSOCIATEDSENSORSHORTNAME=MODIS
  AUTOMATICQUALITYFLAG=Passed
  AUTOMATICQUALITYFLAGEXPLANATION=Passed was set as a default value. More algorithm will be developed
  AVERAGENUMBEROBS=7
  BRDFCODEID=AMBRALS_V4.0R1
  BRDFDATABASEVERSION=v1.0500m
  calibrated_nt=5
  CHARACTERISTICBINANGULARSIZE=15.0
  CHARACTERISTICBINSIZE=463.312716527778
  COVERAGECALCULATIONMETHOD=volume
  DATACOLUMNS=2400
  DATAROWS=2400
  DAYNIGHTFLAG=Day
  Description=Mandatory QA:
  0 = processed, good quality (full BRDF inversions)
  1 = processed, see other QA (magnitude BRDF inversions)

  DESCRREVISION=6.1
  EASTBOUNDINGCOORDINATE=-93.3304652293029
  EXCLUSIONGRINGFLAG=N
  GEOANYABNORMAL=False
  GEOESTMAXRMSERROR=75.0
  GLOBALGRIDCOLUMNS=86400
  GLOBALGRIDROWS=43200
  GRINGPOINTLATITUDE=49.5789113103597, 59.7052699100447, 60.1231268939801, 49.7821369581876
  GRINGPOINTLONGITUDE=-109.970765585882, -143.13163662668, -120.324853219622, -92.9933532661929
  GRINGPOINTSEQUENCENO=1, 2, 3, 4
  HDFEOSVersion=HDFEOS_V2.19
  HORIZONTALTILENUMBER=11
  identifier_product_doi=10.5067/MODIS/MCD43A4.006
  identifier_product_doi_authority=http://dx.doi.org
  INPUTPOINTER=MYD09GA.A2017177.h11v03.006.2017179030854.hdf, MYD09GA.A2017178.h11v03.006.2017180024736.hdf, MYD09GA.A2017179.h11v03.006.2017181181115.hdf, MYD09GA.A2017180.h11v03.006.2017182033059.hdf, MYD09GA.A2017181.h11v03.006.2017183034951.hdf, MYD09GA.A2017182.h11v03.006.2017184032100.hdf, MYD09GA.A2017183.h11v03.006.2017185025111.hdf, MYD09GA.A2017184.h11v03.006.2017191185535.hdf, MYD09GA.A2017185.h11v03.006.2017193215521.hdf, MYD09GA.A2017186.h11v03.006.2017188155338.hdf, MYD09GA.A2017187.h11v03.006.2017189035003.hdf, MYD09GA.A2017188.h11v03.006.2017190032826.hdf, MYD09GA.A2017189.h11v03.006.2017191032733.hdf, MYD09GA.A2017190.h11v03.006.2017192142155.hdf, MYD09GA.A2017191.h11v03.006.2017193031227.hdf, MYD09GA.A2017192.h11v03.006.2017194032826.hdf, MOD09GA.A2017177.h11v03.006.2017179030736.hdf, MOD09GA.A2017178.h11v03.006.2017180025015.hdf, MOD09GA.A2017179.h11v03.006.2017181030953.hdf, MOD09GA.A2017180.h11v03.006.2017187175807.hdf, MOD09GA.A2017181.h11v03.006.2017183034952.hdf, MOD09GA.A2017182.h11v03.006.2017184030419.hdf, MOD09GA.A2017183.h11v03.006.2017185025111.hdf, MOD09GA.A2017184.h11v03.006.2017186031131.hdf, MOD09GA.A2017185.h11v03.006.2017187025242.hdf, MOD09GA.A2017186.h11v03.006.2017193221446.hdf, MOD09GA.A2017187.h11v03.006.2017189031337.hdf, MOD09GA.A2017188.h11v03.006.2017190030447.hdf, MOD09GA.A2017189.h11v03.006.2017191030919.hdf, MOD09GA.A2017190.h11v03.006.2017192030346.hdf, MOD09GA.A2017191.h11v03.006.2017193030756.hdf, MOD09GA.A2017192.h11v03.006.2017194030750.hdf, MCD43DB.A2017184.h11v03.006.2017193032239.hdf
  LOCALGRANULEID=MCD43A4.A2017185.h11v03.006.2017194033101.hdf
  LOCALVERSIONID=6.1.34
  LONGNAME=MODIS/Terra+Aqua BRDF/Albedo Nadir BRDF-Adjusted Ref Daily L3 Global - 500m
  long_name=Nadir_Reflectance_Band1
  MAXIMUMOBSERVATIONS=0
  NADIRDATARESOLUTION=500m
  NORTHBOUNDINGCOORDINATE=59.9999999946118
  NUMBEROFGRANULES=1
  PARAMETERNAME=NOT SET
  PERCENTLANDINTILE=78
  PERCENTNEWBRDFS=61
  PERCENTPROCESSEDINTILE=82
  PERCENTSHAPEFIXEDBRDFS=38
  PERCENTSUBSTITUTEBRDFS=0
  PGEVERSION=6.0.39
  PROCESSINGCENTER=MODAPS
  PROCESSINGENVIRONMENT=Linux minion7039 3.10.0-514.10.2.el7.x86_64 #1 SMP Fri Mar 3 00:04:05 UTC 2017 x86_64 x86_64 x86_64 GNU/Linux
  PRODUCTIONDATETIME=2017-07-13T03:33:02.000Z
  QAPERCENTGOODQUALITY=50
  QAPERCENTNOTPRODUCEDCLOUD=7
  QAPERCENTNOTPRODUCEDOTHER=10
  QAPERCENTOTHERQUALITY=31
  RANGEBEGINNINGDATE=2017-06-26
  RANGEBEGINNINGTIME=00:00:00.000000
  RANGEENDINGDATE=2017-07-11
  RANGEENDINGTIME=23:59:59.999999
  REPROCESSINGACTUAL=processed once
  REPROCESSINGPLANNED=further update is anticipated
  scale_factor=0.0001
  scale_factor_err=0
  SCIENCEQUALITYFLAG=Not Investigated
  SCIENCEQUALITYFLAGEXPLANATION=See http://landweb.nascom/nasa.gov/cgi-bin/QA_WWW/qaFlagPage.cgi?sat=aqua the product Science Quality status.
  SETUPFILEID=06121997
  SHORTNAME=MCD43A4
  SOUTHBOUNDINGCOORDINATE=49.9999999955098
  SPSOPARAMETERS=2015
  TileID=51011003
  units=reflectance, no units
  valid_range=0, 32767
  VERSIONID=6
  VERTICALTILENUMBER=3
  WESTBOUNDINGCOORDINATE=-139.999999964624
  _FillValue=32767
Image Structure Metadata:
  COMPRESSION=DEFLATE
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left  (-7783653.638, 6671703.118) (140d 0' 0.00"W, 60d 0' 0.00"N)
Lower Left  (-7783653.638, 5559752.598) (108d54' 2.40"W, 50d 0' 0.00"N)
Upper Right (-6671703.118, 6671703.118) (120d 0' 0.00"W, 60d 0' 0.00"N)
Lower Right (-6671703.118, 5559752.598) ( 93d20'36.35"W, 50d 0' 0.00"N)
Center      (-7227678.378, 6115727.858) (113d19'26.55"W, 55d 0' 0.00"N)
Band 1 Block=512x512 Type=Int16, ColorInterp=Gray
  Description = Nadir_Reflectance_Band1
  NoData Value=32767
  Offset: 0,   Scale:0.0001

Reprojecting:

$ gdalwarp -t_srs "+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" MCD43A4.A2017185.h11v03.006.2017194033101_B01.TIF out.tif
Creating output file that is 3177P x 1195L.
Processing input file MCD43A4.A2017185.h11v03.006.2017194033101_B01.TIF.
Using internal nodata values (e.g. 32767) for image MCD43A4.A2017185.h11v03.006.2017194033101_B01.TIF.
Copying nodata values from source MCD43A4.A2017185.h11v03.006.2017194033101_B01.TIF to destination out.tif.
0...10...20...30...40...50...60...70...80...90...100 - done.

In the reprojected image info I see that the geographical extent of the image has been changed. Why did that happen?

Pixel size is not 1/cos(lat) * original pixel size. Does it mean that MODIS pixel size is not in 'real meters'? How 1634.89 pixel size had been calculated?

What is the way to calculate output raster size 3177P x 1195L? I've tried reprojecting the image extent polygon into Web Mercator, but got stuck what to do next

$ gdalinfo out.tif 
Driver: GTiff/GeoTIFF
Files: out.tif
Size is 3177, 1195
Coordinate System is:
PROJCS["unnamed",
    GEOGCS["unnamed ellipse",
        DATUM["unknown",
            SPHEROID["unnamed",6378137,0]],
        PRIMEM["Greenwich",0],
        UNIT["degree",0.0174532925199433]],
    PROJECTION["Mercator_1SP"],
    PARAMETER["central_meridian",0],
    PARAMETER["scale_factor",1],
    PARAMETER["false_easting",0],
    PARAMETER["false_northing",0],
    UNIT["metre",1,
        AUTHORITY["EPSG","9001"]]]
Origin = (-15584728.707120897248387,8399737.888618739321828)
Pixel Size = (1634.890398670264403,-1634.890398670264403)
Metadata:
  add_offset=0
  add_offset_err=0
  ALBEDOFILEID=06121997
  AREA_OR_POINT=Area
  ASSOCIATEDINSTRUMENTSHORTNAME=MODIS
  ASSOCIATEDPLATFORMSHORTNAME=Terra
  ASSOCIATEDSENSORSHORTNAME=MODIS
  AUTOMATICQUALITYFLAG=Passed
  AUTOMATICQUALITYFLAGEXPLANATION=Passed was set as a default value. More algorithm will be developed
  AVERAGENUMBEROBS=7
  BRDFCODEID=AMBRALS_V4.0R1
  BRDFDATABASEVERSION=v1.0500m
  calibrated_nt=5
  CHARACTERISTICBINANGULARSIZE=15.0
  CHARACTERISTICBINSIZE=463.312716527778
  COVERAGECALCULATIONMETHOD=volume
  DATACOLUMNS=2400
  DATAROWS=2400
  DAYNIGHTFLAG=Day
  Description=Mandatory QA:
  0 = processed, good quality (full BRDF inversions)
  1 = processed, see other QA (magnitude BRDF inversions)

  DESCRREVISION=6.1
  EASTBOUNDINGCOORDINATE=-93.3304652293029
  EXCLUSIONGRINGFLAG=N
  GEOANYABNORMAL=False
  GEOESTMAXRMSERROR=75.0
  GLOBALGRIDCOLUMNS=86400
  GLOBALGRIDROWS=43200
  GRINGPOINTLATITUDE=49.5789113103597, 59.7052699100447, 60.1231268939801, 49.7821369581876
  GRINGPOINTLONGITUDE=-109.970765585882, -143.13163662668, -120.324853219622, -92.9933532661929
  GRINGPOINTSEQUENCENO=1, 2, 3, 4
  HDFEOSVersion=HDFEOS_V2.19
  HORIZONTALTILENUMBER=11
  identifier_product_doi=10.5067/MODIS/MCD43A4.006
  identifier_product_doi_authority=http://dx.doi.org
  INPUTPOINTER=MYD09GA.A2017177.h11v03.006.2017179030854.hdf, MYD09GA.A2017178.h11v03.006.2017180024736.hdf, MYD09GA.A2017179.h11v03.006.2017181181115.hdf, MYD09GA.A2017180.h11v03.006.2017182033059.hdf, MYD09GA.A2017181.h11v03.006.2017183034951.hdf, MYD09GA.A2017182.h11v03.006.2017184032100.hdf, MYD09GA.A2017183.h11v03.006.2017185025111.hdf, MYD09GA.A2017184.h11v03.006.2017191185535.hdf, MYD09GA.A2017185.h11v03.006.2017193215521.hdf, MYD09GA.A2017186.h11v03.006.2017188155338.hdf, MYD09GA.A2017187.h11v03.006.2017189035003.hdf, MYD09GA.A2017188.h11v03.006.2017190032826.hdf, MYD09GA.A2017189.h11v03.006.2017191032733.hdf, MYD09GA.A2017190.h11v03.006.2017192142155.hdf, MYD09GA.A2017191.h11v03.006.2017193031227.hdf, MYD09GA.A2017192.h11v03.006.2017194032826.hdf, MOD09GA.A2017177.h11v03.006.2017179030736.hdf, MOD09GA.A2017178.h11v03.006.2017180025015.hdf, MOD09GA.A2017179.h11v03.006.2017181030953.hdf, MOD09GA.A2017180.h11v03.006.2017187175807.hdf, MOD09GA.A2017181.h11v03.006.2017183034952.hdf, MOD09GA.A2017182.h11v03.006.2017184030419.hdf, MOD09GA.A2017183.h11v03.006.2017185025111.hdf, MOD09GA.A2017184.h11v03.006.2017186031131.hdf, MOD09GA.A2017185.h11v03.006.2017187025242.hdf, MOD09GA.A2017186.h11v03.006.2017193221446.hdf, MOD09GA.A2017187.h11v03.006.2017189031337.hdf, MOD09GA.A2017188.h11v03.006.2017190030447.hdf, MOD09GA.A2017189.h11v03.006.2017191030919.hdf, MOD09GA.A2017190.h11v03.006.2017192030346.hdf, MOD09GA.A2017191.h11v03.006.2017193030756.hdf, MOD09GA.A2017192.h11v03.006.2017194030750.hdf, MCD43DB.A2017184.h11v03.006.2017193032239.hdf
  LOCALGRANULEID=MCD43A4.A2017185.h11v03.006.2017194033101.hdf
  LOCALVERSIONID=6.1.34
  LONGNAME=MODIS/Terra+Aqua BRDF/Albedo Nadir BRDF-Adjusted Ref Daily L3 Global - 500m
  long_name=Nadir_Reflectance_Band1
  MAXIMUMOBSERVATIONS=0
  NADIRDATARESOLUTION=500m
  NORTHBOUNDINGCOORDINATE=59.9999999946118
  NUMBEROFGRANULES=1
  PARAMETERNAME=NOT SET
  PERCENTLANDINTILE=78
  PERCENTNEWBRDFS=61
  PERCENTPROCESSEDINTILE=82
  PERCENTSHAPEFIXEDBRDFS=38
  PERCENTSUBSTITUTEBRDFS=0
  PGEVERSION=6.0.39
  PROCESSINGCENTER=MODAPS
  PROCESSINGENVIRONMENT=Linux minion7039 3.10.0-514.10.2.el7.x86_64 #1 SMP Fri Mar 3 00:04:05 UTC 2017 x86_64 x86_64 x86_64 GNU/Linux
  PRODUCTIONDATETIME=2017-07-13T03:33:02.000Z
  QAPERCENTGOODQUALITY=50
  QAPERCENTNOTPRODUCEDCLOUD=7
  QAPERCENTNOTPRODUCEDOTHER=10
  QAPERCENTOTHERQUALITY=31
  RANGEBEGINNINGDATE=2017-06-26
  RANGEBEGINNINGTIME=00:00:00.000000
  RANGEENDINGDATE=2017-07-11
  RANGEENDINGTIME=23:59:59.999999
  REPROCESSINGACTUAL=processed once
  REPROCESSINGPLANNED=further update is anticipated
  scale_factor=0.0001
  scale_factor_err=0
  SCIENCEQUALITYFLAG=Not Investigated
  SCIENCEQUALITYFLAGEXPLANATION=See http://landweb.nascom/nasa.gov/cgi-bin/QA_WWW/qaFlagPage.cgi?sat=aqua the product Science Quality status.
  SETUPFILEID=06121997
  SHORTNAME=MCD43A4
  SOUTHBOUNDINGCOORDINATE=49.9999999955098
  SPSOPARAMETERS=2015
  TileID=51011003
  units=reflectance, no units
  valid_range=0, 32767
  VERSIONID=6
  VERTICALTILENUMBER=3
  WESTBOUNDINGCOORDINATE=-139.999999964624
  _FillValue=32767
Image Structure Metadata:
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left  (-15584728.707, 8399737.889) (140d 0' 0.00"W, 60d 0' 0.00"N)
Lower Left  (-15584728.707, 6446043.862) (140d 0' 0.00"W, 49d59'55.18"N)
Upper Right (-10390681.911, 8399737.889) ( 93d20'27.90"W, 60d 0' 0.00"N)
Lower Right (-10390681.911, 6446043.862) ( 93d20'27.90"W, 49d59'55.18"N)
Center      (-12987705.309, 7422890.875) (116d40'13.95"W, 55d18'47.52"N)
Band 1 Block=3177x1 Type=Int16, ColorInterp=Gray
  Description = Nadir_Reflectance_Band1
  NoData Value=32767

What mdsumner said, plus you're converting from an equal area projection to a conformal (maintains shape, not area) projection that has dramatic north-south distortion. — mkennedy, Jul 18, 2017 at 21:38

PolyGeo · Accepted Answer · 2017-12-27 08:09:24Z

Reprojecting like this really means the target ends up very different from the source - consider what the Sinusoidal projection looks like: http://proj4.org/projections/sinu.html

The reason it's in this projection is because the L3-bins - daily, equal-area, statistical bins - from which this mapping is generated are not a "grid" in their native form, they are a kind of ragged array in this Sinusoidal projection. I.e. there's as many bins per latitude as required, not the cosine latitude hack.

Difference between MODIS sinusoidal and integerized sinusoidal grid

So, it's very "brute force" the transformation you are applying, it might be fine but it's not exactly light weight. If you are surprised by the shape of the outcome it might be a warning that the modification of the data is no suitable either.

Stack Exchange Network

Reprojecting from MODIS Sinusoidal to Web Mercator: how to calculate output raster dimensions?

1 Answer 1

Your Answer

Not the answer you're looking for? Browse other questions tagged
coordinate-system
gdal
gdalwarp
modis
web-mercator
or ask your own question.

Linked

Hot Network Questions

Reprojecting from MODIS Sinusoidal to Web Mercator: how to calculate output raster dimensions?

1 Answer 1

Your Answer

Sign up or log in

Post as a guest

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

Linked

Related

Hot Network Questions

Not the answer you're looking for? Browse other questions tagged
coordinate-system
gdal
gdalwarp
modis
web-mercator
or ask your own question.