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.

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

Add field-line traces through satellites #4

Open
wants to merge 1 commit into
base: main
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
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
25 changes: 25 additions & 0 deletions enlilviz/enlil.py
Original file line number Diff line number Diff line change
Expand Up @@ -106,6 +106,31 @@ def get_satellite_data(self, satellite, var, coord=None):
varname += '_' + coord
return self.ds.loc[{'satellite': satellite}][varname]

@_validate_satellite
def get_satellite_fieldline(self, satellite, time):
"""Returns the position of the fieldline passing through the satellite.

Parameters
----------
satellite : str
Satellite of interest.
time : datetime-like
Time of interest. Looks for nearest time.

Returns
-------
r : float
Radial positions of the fieldline (AU)
lat : float
Latitudinal positions of the fieldline (deg).
lon : float
Longitudinal positions of the fieldline (deg).
"""
# Fieldlines are index based on 't' not 'earth_t'
ds = self.ds.sel({'t': time}, method='nearest')
ds = ds.loc[{'satellite': satellite}]
return (ds['fld_pos_r'], ds['fld_pos_lat'], ds['fld_pos_lon'])

def get_slice(self, var, slice_plane, time=None):
"""Get a 2D slice of data.

Expand Down
6 changes: 3 additions & 3 deletions enlilviz/io.py
Original file line number Diff line number Diff line change
Expand Up @@ -250,13 +250,13 @@ def _transform_variable(da):

# Position validation is coordinate dependent so do those
# conversions first
if da.name == 'pos_r':
if 'pos_r' in da.name:
da /= _unit_conversion['AU']
da.attrs['units'] = 'AU'
elif da.name == 'pos_lat':
elif 'pos_lat' in da.name:
da.values = np.rad2deg(np.pi/2 - da)
da.attrs['units'] = 'degrees_north'
elif da.name == 'pos_lon':
elif 'pos_lon' in da.name:
da.values = np.rad2deg(da - np.pi)
da.attrs['units'] = 'degrees_east'

Expand Down
46 changes: 46 additions & 0 deletions enlilviz/plotting/plots.py
Original file line number Diff line number Diff line change
Expand Up @@ -222,6 +222,18 @@ def _init_plot(self):
color=SAT_COLORS[sat],
markeredgecolor='k', markersize=10, zorder=2)
self.plot_data[sat] = marker
# Fieldlines
pos = run.get_satellite_fieldline(sat, self.time)
mask = np.logical_or(pos[0].values <= RMIN,
pos[0].values >= RMAX)

lon_field = np.ma.masked_array(np.deg2rad(pos[2]), mask=mask)
r_field = np.ma.masked_array(pos[0], mask=mask)

line, = ax.plot(lon_field, r_field,
color=SAT_COLORS[sat],
linestyle='-.', zorder=1)
self.plot_data[sat + '_fieldline'] = line

# Circle on Earth (1AU)
# ax.axvline(0., c='k', zorder=1)
Expand Down Expand Up @@ -258,6 +270,17 @@ def update(self):
pos = run.get_satellite_position(sat, self.time)
self.plot_data[sat].set_data(np.deg2rad(pos[2]), pos[0])

# Fieldlines
pos = run.get_satellite_fieldline(sat, self.time)
mask = np.logical_or(pos[0].values <= RMIN,
pos[0].values >= RMAX)

lon_field = np.ma.masked_array(np.deg2rad(pos[2]), mask=mask)
r_field = np.ma.masked_array(pos[0], mask=mask)

self.plot_data[sat + '_fieldline'].set_data(
lon_field, r_field)


class LongitudeSlice(_BasePlot):
"""A longitude polar slice plot, which is in latitude-radial coordinates.
Expand Down Expand Up @@ -310,6 +333,19 @@ def _init_plot(self):
markeredgecolor='k', markersize=10, zorder=2)
self.plot_data[sat] = marker

# Fieldlines
pos = run.get_satellite_fieldline(sat, self.time)
mask = np.logical_or(pos[0].values <= RMIN,
pos[0].values >= RMAX)

lat_field = np.ma.masked_array(np.deg2rad(pos[1]), mask=mask)
r_field = np.ma.masked_array(pos[0], mask=mask)

line, = ax.plot(lat_field, r_field,
color=SAT_COLORS[sat],
linestyle='-.', zorder=1)
self.plot_data[sat + '_fieldline'] = line

# x == theta, y == r
circle_points = np.linspace(0, 2*np.pi)
ax.plot(circle_points, np.ones(len(circle_points)), c='k',
Expand Down Expand Up @@ -342,6 +378,16 @@ def update(self):
for sat in ['Earth']:
pos = run.get_satellite_position(sat, self.time)
self.plot_data[sat].set_data(np.deg2rad(pos[1]), pos[0])
# Fieldlines
pos = run.get_satellite_fieldline(sat, self.time)
mask = np.logical_or(pos[0].values <= RMIN,
pos[0].values >= RMAX)

lat_field = np.ma.masked_array(np.deg2rad(pos[1]), mask=mask)
r_field = np.ma.masked_array(pos[0], mask=mask)

self.plot_data[sat + '_fieldline'].set_data(
lat_field, r_field)


class RadialSlice(_BasePlot):
Expand Down