Reprojecting GEOTIFF to epsg:3857 with scrambled CRS data (possibly GCJ-02)


toku2288@...
 

Hello, I am trying to reproject a population distribution raster file of China into 3857 space. I tried but the images are not aligning. Looking at the CRS it has a lot of different standards, I am not sure what to make of it.

The meta of the source file looks something like this:
{'driver': 'GTiff', 
'dtype': 'int32', 'nodata': 2147483647.0, 'width': 5049, 'height': 5582, 'count': 1, 'crs': CRS.from_wkt('PROJCS["Albers_Conic_Equal_Area",GEOGCS["GCS_Unknown_datum_based_upon_the_Krassowsky_1940_ellipsoid",DATUM["Not_specified_based_on_Krassowsky_1940_ellipsoid",SPHEROID["Krasovsky_1940",6378245,298.3,AUTHORITY["EPSG","7024"]],AUTHORITY["EPSG","6024"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433]],PROJECTION["Albers_Conic_Equal_Area"],PARAMETER["standard_parallel_1",25],PARAMETER["standard_parallel_2",47],PARAMETER["latitude_of_center",0],PARAMETER["longitude_of_center",105],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]]]'), 'transform': Affine(1000.0, 0.0, -2638055.0614613006, 0.0, -1000.0, 5922165.484718928)}
PopTIF.crs.data returns:
{'proj': 'aea',
 'lat_1': 25,
 'lat_2': 47,
 'lat_0': 0,
 'lon_0': 105,
 'x_0': 0,
 'y_0': 0,
 'ellps': 'krass',
 'units': 'm',
 'no_defs': True}
after running the reprojection algorithm:
{'driver': 'GTiff',
 'dtype': 'int32',
 'nodata': 2147483647.0,
 'width': 6708,
 'height': 6102,
 'count': 1,
 'crs': CRS.from_dict(init='epsg:3857'),
 'transform': Affine(1184.7466003556895, 0.0, 7551104.353792087,
        0.0, -1184.7466003556895, 7299154.069454485)}




I am then using PopTIF.dataset_mask()  to mask the area of a Heightfield GEOTIFF (I requested the raster area of this file to be the same as the lat/lon of the PopTIF data:  'lat_1': 25,'lat_2': 47, 'lat_0': 0, 'lon_0': 105
 {'driver': 'GTiff', 
'dtype': 'float32', 'nodata': -32767.0, 'width': 2816, 'height': 2288, 'count': 1, 'crs': CRS.from_dict(init='epsg:4326'), 'transform': Affine(0.02197265625, 0.0, 73.4765625, 0.0, -0.02197265625, 53.7890625)}

 EPSG:4326 to EPSG:3857

{'driver': 'GTiff',
 'dtype': 'float32',
 'nodata': -32767.0,
 'width': 2594,
 'height': 2537,
 'count': 1,
 'crs': CRS.from_dict(init='epsg:3857'),
 'transform': Affine(2655.775407243914, 0.0, 8179373.5227401415,
        0.0, -2655.775407243914, 7130308.041557059)}


So, After moving both GEOTIFFs into 3857 space. I then needed to upscale the Heightfield so the array matched the the size of the Population data:

import skimage.transform as st

pathToRaster  = r'./DEMTIF/cnl10_3857.tif'
newsize = reprojectedmask.shape


with rio.open(pathToRaster) as src:
    arr = src.read()
    arr[arr ==-9999] = np.nan
    arr  = arr[0,:,:]
    new_shape = newsize
    newarr= st.resize(arr, new_shape, mode='constant')

Resulting array can be masked, and looks something like this, Notice how there is no height data in areas such as Taiwan, There has been a mismatch with the projections, However I cant tell where/how. This may be due to the fact I specified the lat/lon in AEA(plus other CRS)? Would appreciate any guidance on how to find the correct area to sample and how to transform it 

Join main@rasterio.groups.io to automatically receive all group messages.