satpy.readers.hrit_jma module
HRIT format reader for JMA data.
Introduction
The JMA HRIT format is described in the JMA HRIT - Mission Specific Implementation. There are three readers for this format in Satpy:
jami_hrit
: For data from the JAMI instrument on MTSAT-1Rmtsat2-imager_hrit
: For data from the Imager instrument on MTSAT-2ahi_hrit
: For data from the AHI instrument on Himawari-8/9
Although the data format is identical, the instruments have different characteristics, which is why there is a dedicated reader for each of them. Sample data is available here:
Example:
Here is an example how to read Himwari-8 HRIT data with Satpy:
from satpy import Scene
import glob
filenames = glob.glob('data/IMG_DK01B14_2018011109*')
scn = Scene(filenames=filenames, reader='ahi_hrit')
scn.load(['B14'])
print(scn['B14'])
Output:
<xarray.DataArray (y: 5500, x: 5500)>
dask.array<concatenate, shape=(5500, 5500), dtype=float64, chunksize=(550, 4096), ...
Coordinates:
acq_time (y) datetime64[ns] 2018-01-11T09:00:20.995200 ... 2018-01-11T09:09:40.348800
crs object +proj=geos +lon_0=140.7 +h=35785831 +x_0=0 +y_0=0 +a=6378169 ...
* y (y) float64 5.5e+06 5.498e+06 5.496e+06 ... -5.496e+06 -5.498e+06
* x (x) float64 -5.498e+06 -5.496e+06 -5.494e+06 ... 5.498e+06 5.5e+06
Attributes:
orbital_parameters: {'projection_longitude': 140.7, 'projection_latitud...
standard_name: toa_brightness_temperature
level: None
wavelength: (11.0, 11.2, 11.4)
units: K
calibration: brightness_temperature
file_type: ['hrit_b14_seg', 'hrit_b14_fd']
modifiers: ()
polarization: None
sensor: ahi
name: B14
platform_name: Himawari-8
resolution: 4000
start_time: 2018-01-11 09:00:20.995200
end_time: 2018-01-11 09:09:40.348800
area: Area ID: FLDK, Description: Full Disk, Projection I...
ancillary_variables: []
JMA HRIT data contain the scanline acquisition time for only a subset of scanlines. Timestamps of
the remaining scanlines are computed using linear interpolation. This is what you’ll find in the
acq_time
coordinate of the dataset.
Compression
Gzip-compressed MTSAT files can be decompressed on the fly using
FSFile
:
import fsspec
from satpy import Scene
from satpy.readers import FSFile
filename = "/data/HRIT_MTSAT1_20090101_0630_DK01IR1.gz"
open_file = fsspec.open(filename, compression="gzip")
fs_file = FSFile(open_file)
scn = Scene([fs_file], reader="jami_hrit")
scn.load(["IR1"])
- class satpy.readers.hrit_jma.HRITJMAFileHandler(filename, filename_info, filetype_info, use_acquisition_time_as_start_time=False)[source]
Bases:
HRITFileHandler
JMA HRIT format reader.
By default, the reader uses the start time parsed from the filename. To use exact time, computed from the metadata, the user can define a keyword argument:
scene = Scene(filenames=filenames, reader='ahi_hrit', reader_kwargs={'use_acquisition_time_as_start_time': True})
As this time is different for every channel, time-dependent calculations like SZA correction can be pretty slow when multiple channels are used.
The exact scanline times are always available as coordinates of an individual channels:
scene.load(["B03"]) print(scene["B03].coords["acq_time"].data)
would print something similar to:
array(['2021-12-08T06:00:20.131200000', '2021-12-08T06:00:20.191948000', '2021-12-08T06:00:20.252695000', ..., '2021-12-08T06:09:39.449390000', '2021-12-08T06:09:39.510295000', '2021-12-08T06:09:39.571200000'], dtype='datetime64[ns]')
The first value represents the exact start time, and the last one the exact end time of the data acquisition.
Initialize the reader.
- _check_sensor_platform_consistency(sensor)[source]
Make sure sensor and platform are consistent.
- Parameters:
sensor (str) – Sensor name from YAML dataset definition
- Raises:
ValueError if they don't match –
- _get_acq_time()[source]
Get the acquisition times from the file.
Acquisition times for a subset of scanlines are stored in the header as follows:
b’LINE:=1rTIME:=54365.022558rLINE:=21rTIME:=54365.022664r…’
Missing timestamps in between are computed using linear interpolation.
- _get_line_offset()[source]
Get line offset for the current segment.
Read line offset from the file and adapt it to the current segment or half disk scan so that
y(l) ~ l - loff
because this is what get_geostationary_area_extent() expects.
- _get_platform()[source]
Get the platform name.
The platform is not specified explicitly in JMA HRIT files. For segmented data it is not even specified in the filename. But it can be derived indirectly from the projection name:
GEOS(140.00): MTSAT-1R GEOS(140.25): MTSAT-1R # TODO: Check if there is more… GEOS(140.70): Himawari-8 GEOS(145.00): MTSAT-2
See [MTSAT], section 3.1. Unfortunately Himawari-8 and 9 are not distinguishable using that method at the moment. From [HIMAWARI]:
“HRIT/LRIT files have the same file naming convention in the same format in Himawari-8 and Himawari-9, so there is no particular difference.”
TODO: Find another way to distinguish Himawari-8 and 9.
References: [MTSAT] http://www.data.jma.go.jp/mscweb/notice/Himawari7_e.html [HIMAWARI] http://www.data.jma.go.jp/mscweb/en/himawari89/space_segment/sample_hrit.html
- property end_time
Get end time of the scan.
- property start_time
Get start time of the scan.