<jupyter_start><jupyter_text>CoastWatch Python Exercises   Extracting data at ship stations tutorial> history | updated August 2021  > owner | NOAA CoastWatch In this exercise, you will import a ship track and extract satellite SST values at locations along the track.The exercise demonstrates the following skills:*    Getting data from a remote CSV file*    Downloading satellite data from ERDDAP in netCDF format*    Extracting data with Python  Look for python modules you might not have installed* We will be using the xarray, numpy, and pandas modules for this exercise. * Make sure that they are installed in your Python 3 environment. * A quick way to do this is with the script below. * If you find missing modules, please go to the appendix at the bottom of the page for installation instructions.<jupyter_code>import pkg_resources

# Create a set 'curly brackets' of the modules to look for
# You can put any modules that you want to in the set
required = {'xarray', 'numpy', 'pandas'}

# look for the installed packages
installed = {pkg.key for pkg in pkg_resources.working_set}

# Find which modules are missing
missing = required - installed
if len(missing)==0:
    print('All modules are installed')
else:
    print('These modules are missing', ', '.join(missing))<jupyter_output>All modules are installed<jupyter_text>Import the primary modules used in this tutorial* numpy is used for matrix operations* numpy.ma specifically is used for masked arrays* pandas is used for tabular data* xarray is used for opening the gridded dataset<jupyter_code>import numpy as np  
import numpy.ma as ma 
import pandas as pd  
import xarray as xr<jupyter_output><empty_output><jupyter_text>Generate means for a ship station track In the notebook repo you will see a .csv file file in the /data directory named my_monterey_track.csv. Load a station location file (.csv) with the time, latitude, and longitude coordinates for ship station locations.For this example the CSV file has three columns, with each row locating a station in time and space:  `sample_date	sample_time	sample_lat	sample_lon`  `6/3/20	        18:00	    37.3	    -123.2` Load the CSV file data using the "read_csv" function from the Pandas module.  To use the date and time, we will need to create a Pandas date object. Pandas has a nice option (parse_dates) to create a column date object from the date and time fields as the CSV file is loaded. We will also keep the original sample_date and sample_time by setting the keep_date_col to True. The Xarray module requires dates be formatted in a ISO 8601 date format, so we will also create an iso_date column<jupyter_code># Load the file into a Pandas data frame
# This creates a date object column named 'sample_date_sample_time' 
track_DF = pd.read_csv('../data/my_monterey_track.csv', parse_dates=[['sample_date', 'sample_time']], keep_date_col=True) 


# track_DF['iso_date'] = track_DF.sample_date_sample_time.dt.strftime('%Y-%m-%dT12:%M:%S')
track_DF['iso_date'] = track_DF.sample_date_sample_time.dt.strftime('%Y-%m-%d')
# Display the data
display(track_DF)<jupyter_output><empty_output><jupyter_text>Load the satellite SST data from ERDDAP * We use two functions to get the data, these are described in greater detail in Tutotial 1* point_to_dataset sets up the request url and returns a dataset array* get_data makes a subset request and returns the satellite data<jupyter_code>def point_to_dataset(dataset_id, base_url='https://coastwatch.pfeg.noaa.gov/erddap/griddap'):
    base_url = base_url.rstrip('/')
    full_url = '/'.join([base_url, dataset_id])
    return xr.open_dataset(full_url)

def get_data(my_da, my_var,
             my_lt_min, my_lt_max, 
             my_ln_min, my_ln_max, 
             my_tm_min, my_tm_max
            ):
    
    my_data = my_da[my_var].sel(
                                latitude=slice(my_lt_min, my_lt_max), 
                                longitude=slice(my_ln_min, my_ln_max), 
                                time=slice(my_tm_min, my_tm_max)
                               )
    return my_data
da = point_to_dataset('nesdisGeoPolarSSTN5SQNRT')
da<jupyter_output><empty_output><jupyter_text>Prepare the track data Check that longitudes are in the same formatEnsure the longitude format of the track data is compatible with the SST dataset. Some datasets in ERDDAP use a -180°E to +180°E format, others use 0-360°.* We can see from the metadata printed above that the satellite SST data product has longitudes in a -180 to +180 format.      longitude range: -179.975 to 179.975We can use the following function to check the longitude formats in the satellite and ship track datasets and adjust the track longitudes if needed. Here is how the function works:  * If both longitude lists have the same format, then the station location longitudes are returned unchanged.* If the SST longitudes are -180 to 180 and the station location longitudes are 0-360, then 360 is subtracted from any station location longitude > 180* If the SST longitudes are 0 to 360 and the station location longitudes are -180 to 180, then 360 is added to any station location longitude < 0    * `dataset_lon` is the longitude column from the xarray data array ( sst['longitude'] )  * `location_lon` is the longitude column from the Pandas data frame, ( track_DF['sample_lon'] )<jupyter_code>def is_it_180_to_180(dataset_lon, location_lon):
    datasat_is = np.any((dataset_lon < 0))
    location_is = np.any((location_lon < 0))
    if datasat_is and location_is:
        return location_lon    
    elif not datasat_is and location_is:
        return [ x if x>0 else x+360 for x in location_lon ]
    elif datasat_is and not location_is:
        return [ x if x<=180 else x-360 for x in location_lon ]<jupyter_output><empty_output><jupyter_text>__Call the `is_it_180_to_180` function__* Pass the dataset longitudes and the station location longitudes to the function * Add the results to the Pandas track data frame with the ship station information<jupyter_code>track_DF['sample_lon_corrected'] = is_it_180_to_180(
                                                    sst['longitude'], 
                                                    track_DF['sample_lon']
                                                   )

track_DF<jupyter_output><empty_output><jupyter_text>Select an area around each ship station to get matching satellite dataHere we have a function for extracting satellite (or any gridded data) from ERDDAP at a ship track point. There are options to return only the exact match as a single point or to extract all satellite data within a distance from the ship location.* Use a single latitude and longitude pair to extract SST data, then the data from a single point will be returned. * To increase to chances of getting data (if there are gaps due to clouds for example) and improve your statistics, you can request a box of SST data around the station location, e.g. a latitude range and a longitude range.  * The box can increase by multiples of the spatial resolution of the dataset.     > For this dataset the spatial resolution was 0.05 degrees for the latitude and longitude. That is approximately 5 km resolution. So, a 5 pixel X 5 pixel box around each station location would result in a 25 x 25 km box (625 km^2) of data.  The following `data_box` function generates a latitude and longitude range box for each station location. The function asks for:* Station latitude and longitude* Spatial resolution of satellite dataset. (`da.geospatial_lat_resolution`, `da.geospatial_lon_resolution`)  * Number of pixels you want added to each side of the latitude and longitude center point      * Choosing `2` for latitude, results in a range of five pixels: the station latitude for the center plus two pixels on each side. If you also choose `2` for the longitude, then you would get a 5x5 box around the station location. The default value for the function is `1`, which is a 3x3 box<jupyter_code>def data_box(
             st_lat, st_lon, 
             lat_res, lon_res, 
             lat_pixels = 1, lon_pixels = 1
            ):
    
    lt_min = st_lat - lat_res * lat_pixels
    lt_max = st_lat + lat_res * lat_pixels
    ln_min = st_lon - lon_res * lon_pixels
    ln_max = st_lon + lon_res * lon_pixels

    return lt_min, lt_max, ln_min, ln_max<jupyter_output><empty_output><jupyter_text>Call the data_box function  Add the satellite data request areas to the ship track data frame:```new_Pandas_dataframe_column = data_box(st_lat, st_lon, lat_res, lon_res, lat_pixels = 1, lon_pixels = 1)```Results in the following new columns in the pandas data frame:`track_DF['lat_min']`, `track_DF['lat_min']`, `track_DF['lon_min']`, `track_DF['lon_min']`Inputs:  * st_lat = the latitude column from the Pandas station location data frame, track_DF['sample_lat']* st_lon = the longitude column from the Pandas station location data frame, track_DF['sample_lon_corrected']    * make sure to use the corrected longitude column* lat_res = the spatial resolution for satellite latitude, da.geospatial_lat_resolution* lon_res = the spatial resolution for satellite longitude, da.geospatial_lon_resolution * lat_pixels = the number of pixels to add to each side of the latitude center point (default = 1) * lon_pixels = the number of pixels to add to each side of the longitude center point (default = 1)<jupyter_code>[
 track_DF['lat_min'], track_DF['lat_max'],
 track_DF['lon_min'], track_DF['lon_max']
] = data_box(
             track_DF['sample_lat'], track_DF['sample_lon_corrected'],
             da.geospatial_lat_resolution, da.geospatial_lon_resolution,
             2, 2
            )
                                                                                              
track_DF<jupyter_output><empty_output><jupyter_text>Extract SST data for each station location and generate some basic statistics* We now have the input parameters needed to obtain satellite data for the cruise track* Extract data for each station in the ship track   * Build a __for__ loop to use the `sample_date_sample_time`, `lat_min`, `lat_max`, `lon_min`, and `lon_max` columns from `track_DF`    * Use the `get_data()` function to download satellite data for each station.* Calculate the mean, standard deviation, and n stats, then add them to `track_DF`<jupyter_code>means = []
stdevs = []
ns = []

for index, row in track_DF.iterrows():
    
    my_data = get_data(da, my_var,
                       row['lat_min'], row['lat_max'], 
                       row['lon_min'], row['lon_max'], 
                       row['iso_date'], row['iso_date']
                      )
    
    means.append(np.nanmean(my_data, axis=(1,2))[0])
    stdevs.append(np.nanstd(my_data, axis=(1,2))[0])
    ns.append(np.count_nonzero(my_data, axis=(1,2))[0])

print('means', means)
print('stdevs', stdevs)
print('ns', ns)<jupyter_output>means [14.251244, 14.836244, 16.224995, 16.861244, 17.026245, 17.41437]
stdevs [0.3615224, 0.12873788, 0.09499985, 0.048589565, 0.09873407, 0.10271109]
ns [16, 16, 16, 16, 16, 16]<jupyter_text>>means = the average SST over the selected area at the ship location>>stdevs = An indicator of variability over the selected satellite data cells>>ns = The number of non-zero data points used to calculate the averageAdd the statistics to the track_DF dataframe and export as a csv* Uncomment the last line to save the dataframe as a CSV file called 'mbnms_results.csv'<jupyter_code>track_DF['mean'] = means
track_DF['stdev'] = stdevs
track_DF['n'] = ns
display(track_DF)

# uncomment to save the data
# track_DF.to_csv('crusetrack_results.csv', index = False)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.07129375428963623, 'fraction_numerical': 0.02299245024021963, 'mean_word_length': 3.756017951856385, 'pattern_counts': {'":': 0, '<': 43, '<?xml version=': 0, '>': 49, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 13, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '10', 'dir': 'jupyter-structured-clean-dedup', 'id': '212716', 'license': 'MIT', 'n_tokens_mistral': 3593, 'n_tokens_neox': 3291, 'n_words': 1383, 'path': '2_extract_at_ship_stations_edmw.ipynb', 'repo_name': 'mwengren/CoastWatch_EDMW_2021_tutorial', 'size': '11656'}