Angle interpolation for MODIS data missing

[AranJM]

Describe the bug
When producing an image from level 2 MODIS data, the processing fails for any resolution below 1 km because the software cannot find the sun and satellite angle. These angles are provided in a different file (MOD03) for 1 km resolution, but the software does not interpolate them when trying to produce a finer resolution image.

To Reproduce

files = ['MOD02HKM.A2018246.1320.061.NRT.hdf', 'MOD03.A2018246.1320.061.NRT.hdf']
scn = Scene(filenames = files, reader = 'hdfeos_l1b')
scn.load(['true_color'], resolution = 500)

Expected behavior
A 500m resolution image is loaded.

Actual results
KeyError: 'satellite_longitude'

Environment Info:

• OS: Ubuntu 14.04.3 LTS
• SatPy Version: 0.10.0a0.dev0
• PyResample Version: 1.10.1

[djhoese]

Thanks for the issue. Could you please give the full traceback for the KeyError you got. I know @mraspaud started tracking down the exact issue, but any additional information will be helpful.

[AranJM]

---

KeyError Traceback (most recent call last)
in ()
----> 1 scn.load(['true_color'])

/usr/local/lib/python2.7/dist-packages/satpy/scene.pyc in load(self, wishlist, calibration, resolution, polarization, level, generate, unload, **kwargs)
848 self.read(**kwargs)
849 if generate:
--> 850 keepables = self.generate_composites()
851 else:
852 # don't lose datasets we loaded to try to generate composites

/usr/local/lib/python2.7/dist-packages/satpy/scene.pyc in generate_composites(self, nodes)
763 nodes = set(self.dep_tree.trunk(nodes=required_nodes)) -
764 set(self.datasets.keys())
--> 765 return self._read_composites(nodes)
766
767 def _remove_failed_datasets(self, keepables):

/usr/local/lib/python2.7/dist-packages/satpy/scene.pyc in _read_composites(self, compositor_nodes)
737 keepables = set()
738 for item in compositor_nodes:
--> 739 self._generate_composite(item, keepables)
740 return keepables
741

/usr/local/lib/python2.7/dist-packages/satpy/scene.pyc in _generate_composite(self, comp_node, keepables)
710 composite = compositor(prereq_datasets,
711 optional_datasets=optional_datasets,
--> 712 **self.attrs)
713
714 cid = DatasetID.from_dict(composite.attrs)

/usr/local/lib/python2.7/dist-packages/satpy/composites/init.pyc in call(self, projectables, optional_datasets, **info)
412 if not optional_datasets or len(optional_datasets) != 4:
413 vis, red = self.check_areas(projectables)
--> 414 sata, satz, suna, sunz = self.get_angles(vis)
415 red.data = da.rechunk(red.data, vis.data.chunks)
416 else:

/usr/local/lib/python2.7/dist-packages/satpy/composites/init.pyc in get_angles(self, vis)
396 sunz = sun_zenith_angle(vis.attrs['start_time'], lons, lats)
397 sata, satel = get_observer_look(
--> 398 vis.attrs['satellite_longitude'],
399 vis.attrs['satellite_latitude'],
400 vis.attrs['satellite_altitude'],

KeyError: 'satellite_longitude'

[mraspaud]

Thanks @AranJM

So the error comes from satpy trying to compute the angle with pyorbital as they are not available as dataset. Computing with pyorbital, the way it is implemented in satpy at the moment only works for geostationary satellites, and expect a satellite_longitude etc... to be provided as a dataset attribute, hence the error.

The deaper cause is that when loading angles for resolutions other than the ones provided in the file (1000m in this case), it doesn't work because interpolation of angles isn't implemented in the reader.

[AranJM]

Hi, would it be possible to get any update on this or any estimated time when it might be resolved, please?

[mraspaud]

I might have time to look at it next week. Of course any help would be highly appreciated!

[honnorat]

I have the very same problem here with Himawari-8 data from EumetCast (hrit_jma reader).

[djhoese]

@honnorat What problem are you exactly running in to? As far as I know the hrit_jma readers don't provide the angle datasets so they should be generated dynamically by pyorbital in any compositors that need them. Since hrit_jma is usually handling Himawari-8 AHI data (right?) you should have the satellite_longitude information as @mraspaud mentioned above.

[honnorat]

I came across this issue because I googled "KeyError: 'satellite_longitude'" and it appears that my problem produces exactly the same traceback and is deeply linked with this one.

After digging into the code, it seems linked to satpy.composites.PSPRayleighReflectance modifier since "green" and "true_color" composites both produce the error, but not the other ones.

[djhoese]

Interesting. This must be a bug with the hrit_jma reader not adding satellite_longitude to the data it returns since @adybbroe has an example making AHI true colors from the ahi_hsd reader here: https://github.com/pytroll/pytroll-examples/blob/master/satpy/ahi_true_color_pyspectral.ipynb

@honnorat Would you mind filing a separate bug for this with a minimal example and the traceback you receive? Thanks.

[honnorat]

Of course. However not before next week, I fear.

[mraspaud]

@AranJM As you can see I just submitted a PR that should fix this issue. It's still work in progress, but if you could try is out and see if it is a step in the right direction, it would be great.

[AranJM]

Hi @mraspaud, thanks! I hope to have time to look at it tomorrow.