Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Approximate angles without TLE #65

Merged
merged 17 commits into from
Jun 18, 2021
Merged

Approximate angles without TLE #65

Show file tree
Hide file tree
Changes from 14 commits
Commits
Show all changes
17 commits
Select commit Hold shift + click to select a range
0b32a88
Approximate angles without tle
Jun 4, 2020
d6b59fe
Created function for satellite angles without TLEs
Jun 25, 2020
8294e90
Do not assume 409 columns
Sep 9, 2020
a650a16
Test the get angles without tle
Sep 9, 2020
6346817
Fixing style errors.
stickler-ci Sep 9, 2020
22a1c79
style (codefactor)
Nov 2, 2020
55dec03
Merge branch 'get_angles_without_tle' of github.com:/ninahakansson/py…
Nov 2, 2020
b3b535d
Fix problem at nadir
Jan 26, 2021
efaa8e0
Only do fix for older versions of pyorbital
May 26, 2021
b273b0a
Added comment about pyorbital issue
May 26, 2021
46c5c36
Added test that trigger pyorbital <= 1.6.1 bug
May 27, 2021
00d5bca
Expanded and renamed test of nadir sat elevation
May 27, 2021
676a4c2
Added a warning about satellite height which is approximated to 850km
Jun 18, 2021
23c1e0f
Remove pass introduced by mistake
Jun 18, 2021
9f5c66a
Merge branch 'main' into get_angles_without_tle
mraspaud Jun 18, 2021
f08178f
Fixing the angle tests
Jun 18, 2021
c53b737
Fix one more reader_test broken at merge
Jun 18, 2021
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
34 changes: 28 additions & 6 deletions pygac/reader.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -34,6 +34,7 @@
import six
import types
import warnings
import pyorbital

from pygac.configuration import get_config
from pygac.utils import (centered_modulus,
Expand All @@ -43,6 +44,7 @@
from pyorbital import astronomy
from pygac.calibration import calibrate_solar, calibrate_thermal
from pygac import gac_io
from distutils.version import LooseVersion

LOG = logging.getLogger(__name__)

Expand Down Expand Up @@ -723,6 +725,24 @@ def get_tle_lines(self):
self.tle_lines = tle1, tle2
return tle1, tle2

def get_sat_angles_without_tle(self):
"""Get satellite angles using lat/lon from data to approximate satellite postition instead of TLE."""
from pyorbital.orbital import get_observer_look as get_observer_look_no_tle
LOG.warning('Approximating satellite height to 850km (TIROS-N OSCAR)!')
sat_alt = 850.0 # km TIROS-N OSCAR
mid_column = int(0.5*self.lons.shape[1])
sat_azi, sat_elev = get_observer_look_no_tle(
self.lons[:, mid_column][:, np.newaxis],
self.lats[:, mid_column][:, np.newaxis], # approximate satellite position
sat_alt, # approximate satellite altitude
self.times[:, np.newaxis],
self.lons, self.lats, 0)
# Sometimes (pyorbital <= 1.6.1) the get_observer_look_not_tle returns nodata instead of 90.
# Problem solved with https://github.com/pytroll/pyorbital/pull/77
if LooseVersion(pyorbital.__version__) <= LooseVersion('1.6.1'):
sat_elev[:, mid_column] = 90
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Could you issue a warning here saying eg that for best results, pyorbital > 1.6.1 is needed?

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Results are not better with pyorbital > 1.6.1. With the fix they are the same (= 90). Without the fix they sometimes get rounded to 89.999999 instead of 90 but I don't think that is better. I added a test to check sat_elev at nadir are 90 degrees.

return sat_azi, sat_elev

def get_angles(self):
"""Get azimuth and zenith angles.

Expand All @@ -743,15 +763,17 @@ def get_angles(self):
"""
self.get_times()
self.get_lonlat()
tle1, tle2 = self.get_tle_lines()
orb = Orbital(self.spacecrafts_orbital[self.spacecraft_id],
line1=tle1, line2=tle2)

sat_azi, sat_elev = orb.get_observer_look(self.times[:, np.newaxis],
self.lons, self.lats, 0)
try:
tle1, tle2 = self.get_tle_lines()
orb = Orbital(self.spacecrafts_orbital[self.spacecraft_id],
line1=tle1, line2=tle2)
sat_azi, sat_elev = orb.get_observer_look(self.times[:, np.newaxis],
self.lons, self.lats, 0)
except IndexError:
sat_azi, sat_elev = self.get_sat_angles_without_tle()

sat_zenith = 90 - sat_elev

sun_zenith = astronomy.sun_zenith_angle(self.times[:, np.newaxis],
self.lons, self.lats)

Expand Down
63 changes: 62 additions & 1 deletion pygac/tests/test_reader.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -298,6 +298,67 @@ def test_get_tle_lines(self, read_tle_file, *mocks):
self.assertEqual(tle1, tle_data[tle_idx])
self.assertEqual(tle2, tle_data[tle_idx + 1])

def test_get_sat_angles_without_tle_at_nadir(self):
"""Test that the get satellite angles without tle at nadir."""
rng = np.random.RandomState(125)
self.reader.lons = rng.rand(100, 409) * 90
self.reader.lats = rng.rand(100, 409) * 90
self.reader.times = np.array([datetime.datetime(1980, 1, 3, 11, 47, 15, 469000) for date in range(100)])
sat_azi, sat_elev = self.reader.get_sat_angles_without_tle()
self.assertEqual(np.sum(np.isnan(sat_elev)), 0)
np.testing.assert_allclose(sat_elev[:, 204], 90., atol=0.01)

def test_get_sat_angles_without_tle(self):
"""Test the get satellite angles without tle."""
# Test data correspond to columns 0:2, 201:208 and 407:409. Extracted like this:
# self.lons[0:5, [0, 1, 201, 202, 203, 204, 205, 206, 207, -2, -1]]
self.reader.lons = np.array([[69.41555135, 68.76815744, 28.04133742, 27.94671757, 27.85220562,
27.7578125, 27.66354783, 27.5694182, 27.47542957, 2.66416611,
2.39739436],
[69.41409536, 68.76979616, 28.00228658, 27.9076628, 27.8131467,
27.71875, 27.62448295, 27.53035209, 27.43636312, 2.61727408,
2.35049275],
[69.42987929, 68.78543423, 27.96407251, 27.86936406, 27.77457923,
27.6796875, 27.58467527, 27.48959853, 27.39453053, 2.5704025,
2.30362323],
[69.44430772, 68.80064104, 27.91910034, 27.82340242, 27.72794715,
27.6328125, 27.53805662, 27.44366144, 27.34959008, 2.53088093,
2.26729486],
[69.47408815, 68.8259859, 27.87666513, 27.78224611, 27.68795045,
27.59375, 27.49962326, 27.40557682, 27.31162435, 2.48359319,
2.21976689]])
self.reader.lats = np.array([[71.62830288, 71.67081539, 69.90976034, 69.89297223, 69.87616536,
69.859375, 69.84262997, 69.82593315, 69.80928089, 61.61334632,
61.466097],
[71.65222644, 71.69455292, 69.9331981, 69.91640991, 69.89960295,
69.8828125, 69.86606741, 69.84937054, 69.83271828, 61.62903602,
61.48182119],
[71.68324681, 71.7256586, 69.96444365, 69.94765637, 69.93085095,
69.9140625, 69.89731955, 69.8806245, 69.86397337, 61.65247459,
61.50526043],
[71.71489701, 71.75715746, 69.98789144, 69.97110185, 69.9542928,
69.9375, 69.92075278, 69.90405431, 69.88740114, 61.66489162,
61.51597801],
[71.7390608, 71.78108857, 70.01126173, 69.99449392, 69.97770882,
69.9609375, 69.9442054, 69.92751555, 69.91086544, 61.67769951,
61.52802445]])
self.reader.times = np.array([datetime.datetime(1980, 1, 3, 11, 47, 15, 469000),
datetime.datetime(1980, 1, 3, 11, 47, 15, 969000),
datetime.datetime(1980, 1, 3, 11, 47, 16, 469000),
datetime.datetime(1980, 1, 3, 11, 47, 16, 969000),
datetime.datetime(1980, 1, 3, 11, 47, 17, 469000)])
expected_sat_azi_0 = np.array([283.09872924, 283.12775589, 283.13951497, 283.14786413, 283.19638805])
expected_sat_azi_201 = np.array([272.85051989, 273.79847634, 272.04794616, 273.01363377, 274.00055397])
expected_sat_azi_408 = np.array([39.77021472, 39.71516966, 39.68104134, 39.60503726, 39.5403431])
expected_sat_elev_0 = np.array([20.94889204, 20.96041284, 20.96368521, 20.96826309, 20.95849212])
expected_sat_elev_204 = np.array([89.24533884, 89.22663677, 89.25079817, 89.24938043, 89.23004118])
Comment on lines +362 to +366
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Were these expected values computed with TLE and then for the tests they are approximated without TLE?

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes, that is how I made theml.

sat_azi, sat_elev = self.reader.get_sat_angles_without_tle()
np.testing.assert_allclose(sat_azi[:, 0], expected_sat_azi_0, atol=1.0)
np.testing.assert_allclose(sat_azi[:, 2], expected_sat_azi_201, atol=35.0) # Azi bad close to center!
np.testing.assert_allclose(sat_azi[:, -1], expected_sat_azi_408, atol=1.0)
np.testing.assert_allclose(sat_elev[:, 0], expected_sat_elev_0, atol=1.0)
np.testing.assert_allclose(sat_elev[:, 5], expected_sat_elev_204, atol=1.0)

@mock.patch('pygac.gac_reader.GACReader._get_corrupt_mask')
def test_get_angles(self, get_corrupt_mask):
"""Test get_angles function of the reader."""
Expand Down Expand Up @@ -346,7 +407,7 @@ class MockConfigParser(object):
def get(self, section, option, **kwargs):
if option == 'tledir':
return 'path/to/TLEs'
elif option == 'tlename':
if option == 'tlename':
return 'tle_file.txt'
get_config.return_value = MockConfigParser()
tle_file = self.reader.get_tle_file()
Expand Down