SWOT

This notebook demonstrates access to SWOT level 2 data. Broad information about the dataset can be found in the PODAAC website (see here)

Requirements to run this notebook

1. Have an Earth Data Login account

2. Have a Bearer Token.

Objectives

To demonstrate a workflow for remote access and plotting of Complex (Level 2 with Groups) SWOT Data via OPeNDAP

Author: Miguel Jimenez-Urias, ‘24

import matplotlib.pyplot as plt
import numpy as np
import requests
from pydap.client import open_url
import json
import cartopy.crs as ccrs

Access EARTHDATA

The access link can be found at PODACC. This may require to be logged on to EarthDataLogin. There is data for 2023 and 2024.

data_url1 = 'https://opendap.earthdata.nasa.gov/collections/C2799438313-POCLOUD/granules/SWOT_GPR_2PfP507_010_20230501_003247_20230501_012352'

Add to session’s headers Token Authorization

edl_token = ""

auth_hdr="Bearer " + edl_token

# pass Token Authorization to a new Session.

my_session = requests.Session()
my_session.headers={"Authorization": auth_hdr}

Create a dataset access via pydap

dataset1 = open_url(data_url1, session=my_session, protocol="dap4")

dataset1.tree()

.SWOT_GPR_2PfP507_010_20230501_003247_20230501_012352.nc
└──data_01
   ├──time
   ├──longitude
   ├──latitude
   ├──rad_side_2_rain_flag
   ├──rad_wet_tropo_cor
   ├──depth_or_elevation
   ├──rad_qual
   ├──rad_water_vapor
   ├──rad_side_2_sea_ice_flag
   ├──pole_tide
   ├──model_dry_tropo_cor_zero_altitude
   ├──solid_earth_tide
   ├──dac
   ├──rad_side_1_surface_type_flag
   ├──surface_classification_flag
   ├──time_tai
   ├──altitude
   ├──rad_side_1_rain_flag
   ├──rad_side_1_sea_ice_flag
   ├──ice_flag
   ├──rad_wet_tropo_cor_interp_qual
   ├──inv_bar_cor
   ├──ocean_tide_fes
   ├──mean_dynamic_topography
   ├──rad_side_2_surface_type_flag
   ├──alt_qual
   ├──internal_tide_hret
   ├──meteo_map_availability_flag
   ├──wind_speed_alt
   ├──wind_speed_alt_mle3
   ├──rain_flag
   ├──mean_sea_surface_cnescls
   ├──geo_qual
   ├──rad_cloud_liquid_water
   └──ku
      ├──iono_cor_alt_filtered_mle3
      ├──sea_state_bias_mle3
      ├──wvf_main_class
      ├──iono_cor_alt_filtered
      ├──range_ocean
      ├──ssha_mle3
      ├──sig0_ocean
      ├──sig0_ocean_mle3
      ├──sea_state_bias
      ├──swh_ocean
      ├──range_ocean_mle3
      ├──ssha
      └──swh_ocean_mle3

Note

PyDAP accesses the remote dataset’s metadata, and no data has been downloaded yet!

This is a dataset pointing to a remote data location

Data remains remote, no data has been downloaded.

dataset1['data_01/mean_dynamic_topography'].shape

(2880,)

dataset1['data_01/time'].shape

(2880,)

print('total array memory: ', dataset1.nbytes/1e9)

total array memory:  0.00031392

Inspect the values

• longitude

• latitude

• time

dataset1['data_01/time'].attributes

{'long_name': 'time in UTC',
 'standard_name': 'time',
 'calendar': 'gregorian',
 'units': 'seconds since 2000-01-01 00:00:00.0',
 'comment': 'Time of measurement in seconds in the UTC time scale since 1 Jan 2000 00:00:00 UTC. [tai_utc_difference] is the difference between TAI and UTC reference time (seconds) for the first measurement of the data set. If a leap second occurs within the data set, the attribute [leap_second] is set to the UTC time at which the leap second occurs',
 'tai_utc_difference': 37.0,
 'leap_second': '0000-00-00 00:00:00',
 'path': '/data_01'}

%%time
dyn_topo = dataset1['data_01/mean_dynamic_topography'][:] # downloads as BaseType - a thin wrapper for numpy arrays

CPU times: user 40 ms, sys: 3.86 ms, total: 43.9 ms
Wall time: 729 ms

Maps

Refers to the coverage of the Satellite track. This is, how the trajectory “maps” with time (i.e. the dimension)

dyn_topo.Maps

('/data_01/longitude', '/data_01/latitude')

longitude1 = dataset1[dyn_topo.Maps[0]][:]
latitude1 = dataset1[dyn_topo.Maps[1]][:]

longitude1.attributes

{'_FillValue': 2147483647,
 'long_name': 'longitude',
 'standard_name': 'longitude',
 'units': 'degrees_east',
 'scale_factor': 1e-06,
 'comment': 'East longitude relative to Greenwich meridian. See SWOT Nadir Altimeter User Handbook. Associated quality flag is orb_state_diode_flag for the OGDR products, orb_state_rest_flag for the IGDR and GDR products',
 'path': '/data_01',
 'checksum': array([4151025818], dtype=uint32),
 'Maps': ()}

latitude1.attributes

{'_FillValue': 2147483647,
 'long_name': 'latitude',
 'standard_name': 'latitude',
 'units': 'degrees_north',
 'scale_factor': 1e-06,
 'comment': 'Positive latitude is North latitude, negative latitude is South latitude. See SWOT Nadir Altimeter User Handbook. Associated quality flag is orb_state_diode_flag for the OGDR products, orb_state_rest_flag for the IGDR and GDR products',
 'path': '/data_01',
 'checksum': array([1388797074], dtype=uint32),
 'Maps': ()}

Decoding data values

xarray decodes time and spatial values internally by default, everytime one accesses the data values, whereas currently there is no such method within pydap to do so. But it is often useful to understand how this works internally, and what type of parameters are used for decoding. Because OPeNDAP is based on the NetCDF data model, it if a CF-compliant software. Below are some of the most used metadata attributes associated for decoding data:

CF - Conventions

In OPeNDAP’s metadata rich datasets, each contains standard attributes used to describe missing data, units in which the data is presented, and any stretching/scaling of the values.

• standard name

• units

• _FillValue

• scale_factor

• off_set

def decode(variable) -> np.ndarray:
    """Decodes the variable BaseType according with atributes:
        _FillValue
        scale_factor
    """
    scale_factor = 1
    _Fill_value = None

    if 'scale_factor' in variable.attributes:
        scale_factor = variable.scale_factor
    if '_FillValue' in variable.attributes:
        data = np.where(variable.data == variable._FillValue, np.nan, variable.data)    
    else:
        data = variable.data
    return scale_factor * data

Lets make some plots!

OPeNDAP does NOT include a plotting service, but OPeNDAP-served data integrates easily with plotting packages like

• Matplotlib

• Cartopy

plt.figure(figsize=(15, 5))
ax = plt.axes(projection=ccrs.PlateCarree())
ax.set_global()
ax.coastlines()
ax.stock_img() # comment this line if you do not want any background color
plt.scatter(x=decode(longitude1), y=decode(latitude1), c=decode(dyn_topo), marker='.',  cmap='jet')
plt.colorbar().set_label(dyn_topo.name + ' ['+dyn_topo.units+']')

Fig. 1 Global map showing the track (trajectory) of the sampled (satellite) data. Values refer to mean dynamic topography.

plt.figure(figsize=(10, 5))
plt.plot(decode(dyn_topo), 'k', lw=3)
plt.ylabel(dyn_topo.name + ' ['+dyn_topo.units+']', fontsize=15)
plt.xlabel('Along Track Samples', fontsize=15)
plt.show()

Fig 2. Along track values. The first values along the track represent the North Atlantic, whereas the steep dropoff in dynamic topography represent the Southern Ocean.