Skip to content

Commit 5d2b1fd

Browse files
authored
Merge pull request #2905 from pdebuyl/mcd12q1_draft
Mcd12q1 draft
2 parents 4aed70b + 9e3b342 commit 5d2b1fd

File tree

6 files changed

+329
-4
lines changed

6 files changed

+329
-4
lines changed

satpy/etc/readers/mcd12q1.yaml

Lines changed: 49 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,49 @@
1+
reader:
2+
name: mcd12q1
3+
short_name: MCD12Q1
4+
long_name: MODIS Level 3 (mcd12Q1) data in HDF-EOS format
5+
description: MODIS HDF-EOS MCD12Q1 L3 Reader
6+
status: Beta
7+
supports_fsspec: false
8+
reader: !!python/name:satpy.readers.yaml_reader.FileYAMLReader
9+
sensors: [modis]
10+
11+
file_types:
12+
modis_mcd12q1_hdf_eos:
13+
file_patterns: ['MCD12Q1.A{start_time:%Y%j}.{tile_id}.{collection:03d}.{production_time:%Y%j%H%M%S}.hdf']
14+
file_reader: !!python/name:satpy.readers.mcd12q1.MCD12Q1HDFFileHandler
15+
16+
datasets:
17+
LC_Type1:
18+
name: LC_Type1
19+
resolution: 500
20+
file_type: modis_mcd12q1_hdf_eos
21+
22+
LC_Type2:
23+
name: LC_Type2
24+
resolution: 500
25+
file_type: modis_mcd12q1_hdf_eos
26+
LC_Type3:
27+
name: LC_Type3
28+
resolution: 500
29+
file_type: modis_mcd12q1_hdf_eos
30+
LC_Type4:
31+
name: LC_Type4
32+
resolution: 500
33+
file_type: modis_mcd12q1_hdf_eos
34+
LC_Type5:
35+
name: LC_Type5
36+
resolution: 500
37+
file_type: modis_mcd12q1_hdf_eos
38+
LC_Prop1:
39+
name: LC_Prop1
40+
resolution: 500
41+
file_type: modis_mcd12q1_hdf_eos
42+
LC_Prop2:
43+
name: LC_Prop2
44+
resolution: 500
45+
file_type: modis_mcd12q1_hdf_eos
46+
LC_Prop3:
47+
name: LC_Prop3
48+
resolution: 500
49+
file_type: modis_mcd12q1_hdf_eos

satpy/readers/hdfeos_base.py

Lines changed: 8 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -89,6 +89,12 @@ def _find_and_run_interpolation(interpolation_functions, src_resolution, dst_res
8989
logger.debug("Interpolating from {} to {}".format(src_resolution, dst_resolution))
9090
return interpolation_function(*args)
9191

92+
def _modis_date(date):
93+
"""Transform a date and time string into a datetime object."""
94+
if len(date) == 19:
95+
return dt.datetime.strptime(date, "%Y-%m-%d %H:%M:%S")
96+
else:
97+
return dt.datetime.strptime(date, "%Y-%m-%d %H:%M:%S.%f")
9298

9399
class HDFEOSBaseFileReader(BaseFileHandler):
94100
"""Base file handler for HDF EOS data for both L1b and L2 products."""
@@ -183,7 +189,7 @@ def start_time(self):
183189
try:
184190
date = (self.metadata["INVENTORYMETADATA"]["RANGEDATETIME"]["RANGEBEGINNINGDATE"]["VALUE"] + " " +
185191
self.metadata["INVENTORYMETADATA"]["RANGEDATETIME"]["RANGEBEGINNINGTIME"]["VALUE"])
186-
return dt.datetime.strptime(date, "%Y-%m-%d %H:%M:%S.%f")
192+
return _modis_date(date)
187193
except KeyError:
188194
return self._start_time_from_filename()
189195

@@ -196,7 +202,7 @@ def end_time(self):
196202
try:
197203
date = (self.metadata["INVENTORYMETADATA"]["RANGEDATETIME"]["RANGEENDINGDATE"]["VALUE"] + " " +
198204
self.metadata["INVENTORYMETADATA"]["RANGEDATETIME"]["RANGEENDINGTIME"]["VALUE"])
199-
return dt.datetime.strptime(date, "%Y-%m-%d %H:%M:%S.%f")
205+
return _modis_date(date)
200206
except KeyError:
201207
return self.start_time
202208

satpy/readers/mcd12q1.py

Lines changed: 144 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,144 @@
1+
#!/usr/bin/env python
2+
# -*- coding: utf-8 -*-
3+
# Copyright (c) 2024 Satpy developers
4+
#
5+
# This file is part of satpy.
6+
#
7+
# satpy is free software: you can redistribute it and/or modify it under the
8+
# terms of the GNU General Public License as published by the Free Software
9+
# Foundation, either version 3 of the License, or (at your option) any later
10+
# version.
11+
#
12+
# satpy is distributed in the hope that it will be useful, but WITHOUT ANY
13+
# WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
14+
# A PARTICULAR PURPOSE. See the GNU General Public License for more details.
15+
#
16+
# You should have received a copy of the GNU General Public License along with
17+
# satpy. If not, see <http://www.gnu.org/licenses/>.
18+
"""MCD12Q1 hdf-eos format reader.
19+
20+
Introduction
21+
------------
22+
23+
The ``mcd12q1`` reader reads MCD12Q1 products in HDF-EOS format.
24+
25+
The 500m product is provided on a sinusoidal grid.
26+
27+
Reference documents and links:
28+
- MODIS land products grid: https://modis-land.gsfc.nasa.gov/MODLAND_grid.html
29+
- User guide: https://lpdaac.usgs.gov/documents/101/MCD12_User_Guide_V6.pdf
30+
- MCD12Q1 v061: MODIS/Terra+Aqua Land Cover Type Yearly L3 Global 500 m SIN Grid
31+
32+
The reader has been tested with:
33+
- MCD12Q1: Land cover data.
34+
35+
To get a list of the available datasets for a given file refer to the "Load data" section in :doc:`../reading`.
36+
37+
"""
38+
import logging
39+
from typing import Iterable
40+
41+
from pyresample import geometry
42+
43+
from satpy.readers.hdfeos_base import HDFEOSBaseFileReader
44+
45+
logger = logging.getLogger(__name__)
46+
47+
48+
class MCD12Q1HDFFileHandler(HDFEOSBaseFileReader):
49+
"""File handler for MCD12Q1 HDF-EOS 500m granules."""
50+
def available_datasets(self, configured_datasets=None):
51+
"""Automatically determine datasets provided by this file."""
52+
# Initialise set of variable names to carry through code
53+
handled_var_names = set()
54+
55+
ds_dict = self.sd.datasets()
56+
57+
for is_avail, ds_info in (configured_datasets or []):
58+
file_key = ds_info.get("file_key", ds_info["name"])
59+
# we must add all variables here even if another file handler has
60+
# claimed the variable. It could be another instance of this file
61+
# type, and we don't want to add that variable dynamically if the
62+
# other file handler defined it by the YAML definition.
63+
handled_var_names.add(file_key)
64+
if is_avail is not None:
65+
# some other file handler said it has this dataset
66+
# we don't know any more information than the previous
67+
# file handler so let's yield early
68+
yield is_avail, ds_info
69+
continue
70+
if self.file_type_matches(ds_info["file_type"]) is None:
71+
# this is not the file type for this dataset
72+
yield None, ds_info
73+
continue
74+
yield file_key in ds_dict.keys(), ds_info
75+
76+
yield from self._dynamic_variables_from_file(handled_var_names)
77+
78+
def _dynamic_variables_from_file(self, handled_var_names: set) -> Iterable[tuple[bool, dict]]:
79+
res = self._get_res()
80+
for var_name in self.sd.datasets().keys():
81+
if var_name in handled_var_names:
82+
# skip variables that YAML had configured
83+
continue
84+
common = {"file_type": "mcd12q1_500m_hdf",
85+
"resolution": res,
86+
"name": var_name}
87+
yield True, common
88+
89+
90+
def _get_res(self):
91+
"""Compute the resolution from the file metadata."""
92+
gridname = self.metadata["GridStructure"]["GRID_1"]["GridName"]
93+
if "MCD12Q1" not in gridname:
94+
raise ValueError("Only MCD12Q1 grids are supported")
95+
96+
resolution_string = self.metadata["ARCHIVEDMETADATA"]["NADIRDATARESOLUTION"]["VALUE"]
97+
if resolution_string[-1] == "m":
98+
return int(resolution_string.removesuffix("m"))
99+
else:
100+
raise ValueError("Cannot parse resolution of MCD12Q1 grid")
101+
102+
103+
def get_dataset(self, dataset_id, dataset_info):
104+
"""Get DataArray for specified dataset."""
105+
dataset_name = dataset_id["name"]
106+
# xxx
107+
dataset = self.load_dataset(dataset_name, dataset_info.pop("category", False))
108+
109+
self._add_satpy_metadata(dataset_id, dataset)
110+
111+
return dataset
112+
113+
def _get_area_extent(self):
114+
"""Get the grid properties."""
115+
# Now compute the data extent
116+
upperleft = self.metadata["GridStructure"]["GRID_1"]["UpperLeftPointMtrs"]
117+
lowerright = self.metadata["GridStructure"]["GRID_1"]["LowerRightMtrs"]
118+
119+
return upperleft[0], lowerright[1], lowerright[0], upperleft[1]
120+
121+
def get_area_def(self, dsid):
122+
"""Get the area definition.
123+
124+
This is fixed, but not defined in the file. So we must
125+
generate it ourselves with some assumptions.
126+
127+
The proj_param string comes from https://lpdaac.usgs.gov/documents/101/MCD12_User_Guide_V6.pdf
128+
"""
129+
proj_param = "proj=sinu +a=6371007.181 +b=6371007.181 +units=m"
130+
131+
# Get the size of the dataset
132+
nrows = self.metadata["GridStructure"]["GRID_1"]["YDim"]
133+
ncols = self.metadata["GridStructure"]["GRID_1"]["XDim"]
134+
135+
# Construct the area definition
136+
area = geometry.AreaDefinition("sinusoidal_modis",
137+
"Tiled sinusoidal L3 MODIS area",
138+
"sinusoidal",
139+
proj_param,
140+
ncols,
141+
nrows,
142+
self._get_area_extent())
143+
144+
return area

satpy/tests/reader_tests/modis_tests/_modis_fixtures.py

Lines changed: 66 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -215,6 +215,30 @@ def _get_l1b_geo_variable_info(filename: str,
215215
variables_info.update(_get_angles_variable_info(geo_resolution))
216216
return variables_info
217217

218+
def _get_l3_land_cover_info() -> dict:
219+
220+
lc_data = np.zeros((2400, 2400), dtype=np.uint8)
221+
222+
variables_info = \
223+
{
224+
"LC_Type1": {"data": lc_data,
225+
"type": SDC.UINT8,
226+
"fill_value": 255,
227+
"attrs": {
228+
"dim_labels": ["YDim:MCD12Q1", "XDim:MCD12Q1"],
229+
},
230+
},
231+
"LC_Type2": {"data": lc_data,
232+
"type": SDC.UINT8,
233+
"fill_value": 255,
234+
"attrs": {
235+
"dim_labels": ["YDim:MCD12Q1", "XDim:MCD12Q1"],
236+
},
237+
},
238+
}
239+
240+
return variables_info
241+
218242

219243
def generate_nasa_l1b_filename(prefix):
220244
"""Generate a filename that follows NASA MODIS L1b convention."""
@@ -331,19 +355,30 @@ def _create_struct_metadata_cmg(ftype: str) -> str:
331355
gridline = 'GridName="MOD09CMG"\n'
332356
upleft = "UpperLeftPointMtrs=(-180000000.000000,90000000.000000)\n"
333357
upright = "LowerRightMtrs=(180000000.000000,-90000000.000000)\n"
358+
XDim=7200
359+
YDim=3600
360+
# Case of a MCD12Q1 file
361+
elif ftype == "MCD12Q1":
362+
gridline = 'GridName="MCD12Q1"\n'
363+
upleft = "UpperLeftPointMtrs=(-8895604.157333,-1111950.519667)\n"
364+
upright = "LowerRightMtrs=(-7783653.637667,-2223901.039333)\n"
365+
XDim=2400
366+
YDim=2400
334367
# Case of a MCD43 file
335368
else:
336369
gridline = 'GridName="MCD_CMG_BRDF_0.05Deg"\n'
337370
upleft = "UpperLeftPointMtrs=(-180.000000,90.000000)\n"
338371
upright = "LowerRightMtrs=(180.000000,-90.000000)\n"
372+
XDim=7200
373+
YDim=3600
339374

340375
struct_metadata_header = ("GROUP=SwathStructure\n"
341376
"END_GROUP=SwathStructure\n"
342377
"GROUP=GridStructure\n"
343378
"GROUP=GRID_1\n"
344379
f"{gridline}\n"
345-
"XDim=7200\n"
346-
"YDim=3600\n"
380+
f"XDim={XDim}\n"
381+
f"YDim={YDim}\n"
347382
f"{upleft}\n"
348383
f"{upright}\n"
349384
"END_GROUP=GRID_1\n"
@@ -598,6 +633,10 @@ def generate_nasa_l2_filename(prefix: str) -> str:
598633
now = dt.datetime.now()
599634
return f"{prefix}_L2.A{now:%Y%j.%H%M}.061.{now:%Y%j%H%M%S}.hdf"
600635

636+
def generate_nasa_l3_tile_filename(prefix: str) -> str:
637+
"""Generate a file name that follows MODIS sinusoidal grid tile pattern."""
638+
now = dt.datetime.now()
639+
return f"{prefix}.A{now:%Y}001.h34v07.061.{now:%Y%j%H%M%S}.hdf"
601640

602641
@pytest.fixture(scope="session")
603642
def modis_l2_nasa_mod35_file(tmpdir_factory) -> list[str]:
@@ -613,6 +652,31 @@ def modis_l2_nasa_mod35_file(tmpdir_factory) -> list[str]:
613652
_create_header_metadata())
614653
return [full_path]
615654

655+
@pytest.fixture(scope="session")
656+
def modis_l3_nasa_mcd12q1_file(tmpdir_factory) -> list[str]:
657+
"""Create a single MOD35 L2 HDF4 file with headers."""
658+
filename = generate_nasa_l3_tile_filename("MCD12Q1")
659+
full_path = str(tmpdir_factory.mktemp("modis_l2").join(filename))
660+
variable_infos = _get_l3_land_cover_info()
661+
archive_header = \
662+
"""GROUP = ARCHIVEDMETADATA
663+
GROUPTYPE = MASTERGROUP
664+
665+
OBJECT = NADIRDATARESOLUTION
666+
NUM_VAL = 1
667+
VALUE = "500m"
668+
END_OBJECT = NADIRDATARESOLUTION
669+
670+
END_GROUP = ARCHIVEDMETADATA
671+
672+
END
673+
"""
674+
create_hdfeos_test_file(full_path,
675+
variable_infos,
676+
_create_struct_metadata_cmg("MCD12Q1"),
677+
_create_core_metadata("MCD12Q1"),
678+
archive_header)
679+
return [full_path]
616680

617681
def generate_nasa_l3_filename(prefix: str) -> str:
618682
"""Generate a file name that follows MODIS 09 L3 convention in a temporary directory."""

satpy/tests/reader_tests/modis_tests/conftest.py

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -32,6 +32,7 @@
3232
modis_l2_nasa_mod06_file,
3333
modis_l2_nasa_mod35_file,
3434
modis_l2_nasa_mod35_mod03_files,
35+
modis_l3_nasa_mcd12q1_file,
3536
modis_l3_nasa_mod09_file,
3637
modis_l3_nasa_mod43_file,
3738
)

0 commit comments

Comments
 (0)