1

I have an HDF5 file of radiation data for the Indian Ocean region. I plotted the data in Python using the following code

import h5py as h5

import matplotlib.pyplot as plt

f = h5.File("D:\Nikka\Geo Data\Downwelling Surface shortwave Flux\New folder\HDF5_LSASAF_MSG_DSSF_MSG-Disk_201710290130.h5", "r", driver='core')
data=f['DSSF'][:]

print(type(data))

print(data.shape)

plt.imshow(data)

plt.show()

f.close()

The plot has columns and lines on the X and Y axis. I'm new to this and I don't even know how to get the projection information out of an HDF file.

I want to extract the value of radiation for a given coordinate from the whole data first. When I'm done for one data, I want to get the radiation value from the data downloaded for the whole year for the same coordinate.

1

TL;DR

You have to find out what is the index in the DSSF array that corresponds to the latitude and longitude you want, and you do so by using the latitude and longitude arrays.


I don't have your file to use here in this answer, so I will use some other file I have at hand as an example.

Let's say I'm interested in soil moisture, so I will use a file from SMAP. I chose a Level-2 product and some date, such that I downloaded the file SMAP_L2_SM_P_E_08013_D_20160801T112156_R14010_001.h5.

HDF Format

First, note that HDF5 (*.h5 files) is not a data format per se, rather it is a container that defines how data is stored in memory. One can organize the data inside HDF files in different ways. This means you have to know how the data is structured for the particular data set you are using.

Opening HDF files

There are a couple ways to figure that out. One is to look on the website you downloaded the data for some documentation. Likely they have some description of how the data is organized and the meaning of each variable inside it. Another easy approach is to open the file in some program that reads HDF files and explore how it is organized. Here are two examples:

Here is how my file looks like when I open it in HDFView. The columns you see in the middle are three variables that I double clicked after opening the file. I'll get there in a minute.

SMAP example

Now, let me try to answer all your questions.

Projection

I think you may not be able to find out right away what the projection is. Again, that is something that you have to check on the documentation. But you can try to guess. When you open the file with viewer, most likely you will see two main groups (like directories inside the HDF file): metadata (Metadata in my example) and the data itself (Soil_Moisture_Retrieval_Data in my example). Explore all the information inside the metadata group, depending on how your data provider organized it, you may find that information there as well. In case you can't find it, try to open it in some GIS tool (e.g., QGIS or ArcGIS should be able to read your file) and try different options, you may be able to narrow down which projection is it if you can make some educated guesses.

Extract data for a given coordinate

In the figure above, I opened three columns, all under Soil_Moisture_Retrieval_Data. The first is the variable soil_moisture, which corresponds to the data I want to extract. The second is latitude. I didn't open longitude but it is available as well. The third column is time (see next section).

Now, I can't tell for sure for your dataset (again, this is something you have to confirm in the documentation), but most likely the following is true for your case as well. What happens is that a given row number in all the vectors correspond to the same data point. For example, the soil_moisture entry in row number 201591 (highlighted in the figure) has the latitude int the same row number in the latitude vector. You see where this goes: it is the same case for longitude, time, and any other variable in the data set.

So, to find the soil moisture value of a given location, I can do something like this:

import h5py
import numpy as np

hdffile = 'SMAP_L2_SM_P_E_08013_D_20160801T112156_R14010_001.h5'
# define a box
box_lat = [39.9767, 40.0559]
box_lon = [-79.6315, -79.5231]

f = h5py.File(hdffile, mode='r')
# this is the full path within the HDF file
name = '/Soil_Moisture_Retrieval_Data/soil_moisture'
# get the geolocation data
latitude = f['/Soil_Moisture_Retrieval_Data/latitude'][:]
longitude = f['/Soil_Moisture_Retrieval_Data/longitude'][:]
lat_index = np.logical_and(latitude > box_lat[0], latitude < box_lat[1])
lon_index = np.logical_and(longitude > box_lon[0], longitude < box_lon[1])
box_index = np.logical_and(lat_index, lon_index)
data = f[name][box_index]

data
>>> array([ 0.2939842], dtype=float32)

Note that I provided the location as a box and not just a point. This is because the coordinates in latitude and longitude refer to the center of a cell (in my example), and if you just use equality (==) you may not find any point with the coordinates you are looking for. In this example, I loaded the HDF file in QGIS and visually inspected what were the coordinates of the point that I was interested on.

Extract data for the whole year

In this case the process is pretty much similar. In my example, each file only has data for one day. If it is the same with your data set, then you have to iterate over your data files, and extract the data value as above for each day/file. If your HDF file has a year worth of data, then when you extract the index corresponding to the given coordinates, you should get multiple indices instead of one (same point, different times), and you can use the time variable (tb_time_utc in my example) to confirm the date.

More code examples

See more examples here, including code to load and plot. The example in my answer was inspired by code available in that link.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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