developmentseed / developmentseed/geolambda
createcopy causing coordinate system to change in geo tiff
- Dominant language
- Dockerfile
- Stars
- 306
- Forks
- 84
- PR merge metrics
- No merged PRs in 30d
Description
I'm attempting to correct a compression issue with a geotiff by creating a copy or using translate which works locally with gdal_translate. The operation below using geolambda oddly causes the projected coordinate system to change from WGS 84 / Pseudo-Mercator to WGS_1984_Web_Mercator_Auxiliary_Sphere. I'm using the latest public geolambda layer referenced from arn:aws:lambda:us-east-1:552188055668:layer:geolambda:4 and geolambda-python layer from arn:aws:lambda:us-east-1:552188055668:layer:geolambda-python:3
Why is the PROJCS changing along with other values? Perhaps there is a better way to accomplish the task but still would like to know why the issue below happens.
**Lambda Code**
```
import logging
from osgeo import gdal
import boto3
from botocore.exceptions import ClientError
def upload_file(file_name, bucket, object_name=None):
```
# If S3 object_name was not specified, use file_name
if object_name is None:
object_name = file_name
# Upload the file
s3_client = boto3.client('s3')
try:
response = s3_client.upload_file(file_name, bucket, object_name)
except ClientError as e:
logging.error(e)
return False
return True
```
def lambda_handler(event, context=None):
s3 = boto3.client('s3')
```
s3_filename = event["s3file"]
fname = event.get('filename', s3_filename)
fname = fname.replace('s3://', '/vsis3/')
bucket_name = 'bucket-name-here'
key = fname.replace('/vsis3/bucket-name-here/', '')
key = key.replace('.tif', '-gdal.tif')
dataset = gdal.Open(fname, gdal.GA_ReadOnly)
#gdal.Translate('/tmp/output.tif', dataset)
# Load the dataset into the virtual filesystem
temp_name = '/vsimem/current.tif'
tiff_driver = gdal.GetDriverByName('GTiff')
tiff_driver.CreateCopy(temp_name, dataset, False, [ 'COMPRESS=LZW' ])
#tiff_driver.Warp(temp_name,dataset,outputSRS='EPSG:4326')
# Read the raw data from the virtual filesystem
f = gdal.VSIFOpenL(temp_name, 'rb')
gdal.VSIFSeekL(f, 0, 2) # seek to end
size = gdal.VSIFTellL(f)
gdal.VSIFSeekL(f, 0, 0) # seek to beginning
data = gdal.VSIFReadL(1, size, f)
gdal.VSIFCloseL(f)
# Upload the raw data to s3
if s3.put_object(Key=key, Bucket=bucket_name, Body=data, ContentLength=size):
details = 'GDAL conversion success. Input: ' + s3_filename + ' and Output: ' + key
status = '200'
summary_msg = 'success'
else:
details = 'GDAL conversion failed. Input: ' + s3_filename + ' and Output: ' + key
status = '400'
summary_msg = 'error'
response = {
'status':status,
'message':{
'summary' : summary_msg,
'details' : details
}
}
return response
gdal.Unlink(temp_name)
dataset = None
f = None
**Input**
`{
"s3file": "s3://bucket-name/project/f228713b-cb88-41aa-8d0d-b728c3022d0e_2019:12:10T04:05:36:343_dsm.tif"
}`
**Result**
`{
"status": "200",
"message": {
"summary": "success",
"details": "GDAL conversion success. Input: s3://bucket-name/project/f228713b-cb88-41aa-8d0d-b728c3022d0e_2019:12:10T04:05:36:343_dsm.tif and Output: s3://bucket-name/project/f228713b-cb88-41aa-8d0d-b728c3022d0e_2019:12:10T04:05:36:343_dsm-gdal.tif"
}
}`
**Original gdalinfo**
`Driver: GTiff/GeoTIFF
Files: f228713b-cb88-41aa-8d0d-b728c3022d0e_2019:12:10T04:05:36:343_dsm.tif
Size is 15974, 14281
Coordinate System is:
PROJCS["WGS 84 / Pseudo-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.0174532925199433,
AUTHORITY["EPSG","9122"]],
AUTHORITY["EPSG","4326"]],
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"]],
AXIS["X",EAST],
AXIS["Y",NORTH],
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","3857"]]
Origin = (-8637199.014870001003146,4665140.998520000837743)
Pixel Size = (0.020000000000000,-0.020000000000000)
Metadata:
AREA_OR_POINT=Area
TIFFTAG_SOFTWARE=pix4dmapper
Image Structure Metadata:
COMPRESSION=LZW
INTERLEAVE=BAND
Corner Coordinates:
Upper Left (-8637199.015, 4665140.999) ( 77d35'21.40"W, 38d36'15.29"N)
Lower Left (-8637199.015, 4664855.379) ( 77d35'21.40"W, 38d36' 8.07"N)
Upper Right (-8636879.535, 4665140.999) ( 77d35'11.07"W, 38d36'15.29"N)
Lower Right (-8636879.535, 4664855.379) ( 77d35'11.07"W, 38d36' 8.07"N)
Center (-8637039.275, 4664998.189) ( 77d35'16.24"W, 38d36'11.68"N)
Band 1 Block=15974x1 Type=Float32, ColorInterp=Gray
NoData Value=-10000`
**Output gdalinfo**
`Driver: GTiff/GeoTIFF
Files: f228713b-cb88-41aa-8d0d-b728c3022d0e_2019:12:10T04:05:36:343_dsm-gdal.tif
Size is 15974, 14281
Coordinate System is:
PROJCS["WGS_1984_Web_Mercator_Auxiliary_Sphere",
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"]],
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"]],
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","3857"]]
Origin = (-8637199.014870001003146,4665140.998520000837743)
Pixel Size = (0.020000000000000,-0.020000000000000)
Metadata:
AREA_OR_POINT=Area
Image Structure Metadata:
INTERLEAVE=BAND
Corner Coordinates:
Upper Left (-8637199.015, 4665140.999) ( 77d35'21.40"W, 38d36'15.29"N)
Lower Left (-8637199.015, 4664855.379) ( 77d35'21.40"W, 38d36' 8.07"N)
Upper Right (-8636879.535, 4665140.999) ( 77d35'11.07"W, 38d36'15.29"N)
Lower Right (-8636879.535, 4664855.379) ( 77d35'11.07"W, 38d36' 8.07"N)
Center (-8637039.275, 4664998.189) ( 77d35'16.24"W, 38d36'11.68"N)
Band 1 Block=15974x1 Type=Float32, ColorInterp=Gray`
**gdal_translate local command that works**
gdal_translate -q f228713b-cb88-41aa-8d0d-b728c3022d0e_2019:12:10T04:05:36:343_dsm.tif f228713b-cb88-41aa-8d0d-b728c3022d0e_2019:12:10T04:05:36:343_dsm-gdal.tif -co COMPRESS=LZW
Contributor guide
No contributing guide indexed for this repository
Assessment
This issue has not been assessed yet.