Geographic Information Systems Asked on December 14, 2020
I am using DBDB-V v6.2 and GDAL, specifically the dbdbv6_level0c.h5 HDF5 file and the GDAL 2.3.3 C# binaries.
Using these lines of code I can easily find the full list of Subdatasets used in the raster:
var path = @"F:testdbdbv6_level0c.h5";
var raster = Gdal.Open(path, Access.GA_ReadOnly);
var metadata = raster.GetMetadata("SUBDATASETS");
This produces an array of strings, pairs of which contain subdataset names and descriptions, for example:
SUBDATASET_150_NAME=HDF5:"F:testdbdbv6_level0c.h5"://geo/0.05000/229_543/depth
SUBDATASET_150_DESC=[600x600] //geo/0.05000/229_543/depth (32-bit floating-point)
This then allows me to apply Gdal.Open(...)
to any of Subdataset NAMEs and open up that particular raster for processing.
It can also be correlated with the basic metadata, which has information such as:
geo_0.05000_229_543_depth_compile_date=20091230
geo_0.05000_229_543_depth_datatype=depth
geo_0.05000_229_543_depth_ellipsoid_offset=0
geo_0.05000_229_543_depth_units_of_measure=meters
geo_0.05000_229_543_depth_vertical_datum=MSL
geo_0.05000_229_543_depth_zoo_value=0
geo_0.05000_229_543_east=272
geo_0.05000_229_543_north=25
geo_0.05000_229_543_security_key=100
geo_0.05000_229_543_south=24.5
geo_0.05000_229_543_west=271.5
Note that 229_543 is the shared identifier between the Subdataset Name and the metadata, which provides the geographic coordinates of interest.
My question is, can I start with known geo coordinates (e.g., 24.75N 271.75E) and somehow determine that the correct Subdataset to open is the one identified with 229_543?
Please note, I am able to open the Subdataset just fine, but what I want is to be able to determine which Subdataset to open based on geo coordinates of interest, without having to iterate the entire metadata collection for every new coordinate of interest… ideally, some sort of algorithm that takes in geo coordinates and spits out DBDB-V Subdataset identifiers.
I sorted this out, and wrote a C# Fiddle to show the answer with examples, if it may help someone in the future. Sample output:
Lat 24.75 / Lon 271.75
For resolution 0.05... id = 229_543
For resolution 0.50... id = 22_54
For resolution 1.00... id = 11_27
For resolution 2.00... id = 5_13
The gist is that an identifier like 229_543 is row index 229, column index 543. The rows and columns are determined by the arc-minute resolution (e.g., 0.05) and the size of the subdatasets (e.g., 600x600).
So, using the example of 24.75N 271.75E as the desired geo point, we can find the identifiers for any of the available resolutions. Each resolution will use a different "scaling factor," determined by
(Subdataset Size / 60 arc-minutes per degree) * Resolution in arc-minutes
For example, resolution 2.0 arc-minutes has scaling factor 20.0.
Each row index is determined by taking the latitude + 90, dividing by the scaling factor, and taking the floor of the result. For 2.0 arc-minute data, the index would be Floor((24.75 + 90) / 20) = 5
.
Each column index is determined by taking the longitude, dividing by the scaling factor, and taking the floor of the result. For 2.0 arc-minute data, the index would be Floor(271.75 / 20) = 13
.
We can then just combine these index values with the desired resolution to create the Subdataset name that includes the desired geo coordinates: "geo/2.00000/5_13/depth".
Answered by OrangeWombat on December 14, 2020
Get help from others!
Recent Questions
Recent Answers
© 2024 TransWikia.com. All rights reserved. Sites we Love: PCI Database, UKBizDB, Menu Kuliner, Sharing RPP