## Writing Modis data to an hdf5 file

I want to consolidate channel 1, channel 31, the lats and the lons in a smaller file so that I don't have to store all 147 Mbytes of the original L1B file. In this notebook I'll write out a new h5 file following the h5py quick start tutorial: http://docs.h5py.org/en/latest/quick.html. To do this, I'll start from the satelliteI_h5 notebook and add a file writing section.



In [0]:

from __future__ import print_function
import os,site
import glob
import h5py
#
# add the lib folder to the path assuming it is on the same
# level as the notebooks folder
#
libdir=os.path.abspath('../lib')
from modismeta_h5 import parseMeta
from h5dump import dumph5




In [0]:

ls



the glob function finds a file using a wildcard to save typing (google: python glob wildcard)



In [0]:

l1b_filename=glob.glob('../data/MYD02*.h5')[0]
print("found l1b file {}".format(l1b_filename))
geom_filename=glob.glob('../data/MYD03*.h5')[0]
print("found geom file {}".format(geom_filename))




In [0]:

l1b_file=h5py.File(l1b_filename)
geom_file=h5py.File(geom_filename)





In [0]:

print(l1b_file['MODIS_SWATH_Type_L1B']['Data Fields']['Band_1KM_Emissive'].shape)
print(l1b_file['MODIS_SWATH_Type_L1B']['Data Fields']['Band_1KM_Emissive'][...])
print(l1b_file['MODIS_SWATH_Type_L1B']['Data Fields']['EV_1KM_Emissive'].shape)



note that channel 31 occurs at index value 10



In [0]:

index31=10



the data is stored as unsigned, 2 byte integers which can hold values from 0 to $2^{16}$ - 1 = 65,535



In [0]:

chan31=l1b_file['MODIS_SWATH_Type_L1B']['Data Fields']['EV_1KM_Emissive'][index31,:,:]
print(chan31.shape,chan31.dtype)




In [0]:

chan31[:3,:3]



we need to apply a scale and offset to convert to radiance (the netcdf module did this for us automatically

$Data = (RawData - offset) \times scale$

this information is included in the attributes of each variable.

(see page 36 of the Modis users guide )



In [0]:




In [0]:

chan31=(chan31 - offset)*scale




In [0]:

%matplotlib inline



histogram the calibrated radiances and show that they lie between 0-10 $W\,m^{-2}\,\mu m^{-1}\,sr^{-1}$



In [0]:

import matplotlib.pyplot as plt
out=plt.hist(chan31.flat)



Now do the same for Channel 1

From the Modis users guide (p. 25, table 3.3.2), we know that the 1 km version of channel 1 is called EV_250_Aggr1km_RefSB and is at index 0, so just get that channel, scale and offset for index 0



In [0]:

reflective=l1b_file['MODIS_SWATH_Type_L1B']['Data Fields']['EV_250_Aggr1km_RefSB'][0,:,:]




In [0]:

scale=l1b_file['MODIS_SWATH_Type_L1B']['Data Fields']['EV_250_Aggr1km_RefSB'].attrs['reflectance_scales']
offset=l1b_file['MODIS_SWATH_Type_L1B']['Data Fields']['EV_250_Aggr1km_RefSB'].attrs['reflectance_offsets']
chan1=(reflective - offset[0])*scale[0]




In [0]:

out=plt.hist(chan1.flat)





In [0]:

the_lon=geom_file['MODIS_Swath_Type_GEO']['Geolocation Fields']['Longitude'][...]
the_lat=geom_file['MODIS_Swath_Type_GEO']['Geolocation Fields']['Latitude'][...]



# Now write these four fields out to a new hdf file



In [0]:

out_name="swath_output.h5"
f=h5py.File(out_name, "w")
dset = f.create_dataset("lattitude", the_lat.shape, dtype=the_lat.dtype)
dset[...]=the_lat[...]
dset = f.create_dataset("longitude", the_lon.shape, dtype=the_lon.dtype)
dset[...]=the_lon[...]
dset = f.create_dataset("channel1", chan1.shape, dtype=chan1.dtype)
dset[...]=chan1[...]
dset = f.create_dataset("channel31", chan31.shape, dtype=chan31.dtype)
dset[...]=chan31[...]





In [0]:



Transfer all the metadata to the new hdf5 file as global attributes



In [0]:




In [0]:

f.close()



Now reopen the output file and dump it to make sure it has all our stuff



In [0]:

f=h5py.File(out_name,'r')
dumph5(f)
f.close()