Just as in this question, I download a sentinel tif file:
from pystac_client import Client
catalog = Client.open("https://earth-search.aws.element84.com/v0")
mysearch = catalog.search(collections=['sentinel-s2-l2a-cogs'],
bbox=[-72.5,40.5,-72,41],
query = {"eo:cloud_cover":{"lt":1}},
datetime='2020-01-01/2021-03-31',
max_items=10)
print(f"{mysearch.matched()} items found")
resdict = mysearch.get_all_items_as_dict()
Select one:
resdict['features'][3]['assets']['B08']
This is a dict
:
{'eo:bands': [{'center_wavelength': 0.8351,
'common_name': 'nir',
'full_width_half_max': 0.145,
'name': 'B08'}],
'gsd': 10,
'href': 'https://sentinel-cogs.s3.us-west-2.amazonaws.com/sentinel-s2-l2a-cogs/18/T/YL/2021/3/S2A_18TYL_20210330_0_L2A/B08.tif',
'proj:shape': [10980, 10980],
'proj:transform': [10, 0, 699960, 0, -10, 4600020, 0, 0, 1],
'roles': ['data'],
'title': 'Band 8 (nir)',
'type': 'image/tiff; application=geotiff; profile=cloud-optimized'}
Let's see what Rasterio can tell us about the file:
src = rasterio.open('https://sentinel-cogs.s3.us-west-2.amazonaws.com/sentinel-s2-l2a-cogs/18/T/YL/2021/3/S2A_18TYL_20210330_0_L2A/B08.tif')
Then, if I do:
src.crs.wkt
I get:
PROJCS["WGS 84 / UTM zone 18N",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["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",-75],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Easting",EAST],AXIS["Northing",NORTH],AUTHORITY["EPSG","32618"]]
This seem to suggest that for Latitude - Longitude conversion, I need to use UTM zone 18N. Lets plot the image:
plt.figure(figsize=(8,8))
src = rasterio.open(resdict['features'][3]['assets']['B08']['href'])
plt.imshow(src.read(1), cmap='pink')
This is the tip of Long Island. On Google Maps, lets plot a point 41, -71.5
on this area:
The coordinates 41, -71.5
seem to be on my image as well. I then use this website to convert 41, -71.5
to UTM, I get this result:
Ie it seem to think that UTM Zone 18T is used.
Using Wikipedia's map, the 18T Zone seems to be correct.
Shouldn't an image showing the tip of Long Island have the same UTM Zone as a coordinate corresponding to the tip of Long Island? Why is this discrepancy?
I am concerned about this because it was suggested in the comment section of this question that to find pixel coordinates from Latitude-Longitude, first I need to convert to UTM. My UTM Easting and UTM Northing value might be wrong if I use the wrong Zone, at least that is my understanding.
(I am new to GIS.)