Atmospheric Radiation Measurement user facility (ARM)¶

TRacking Aerosol Convection interations ExpeRiment (TRACER)¶

July 16-19, 2022 Dust Event¶

Notebook for data exploration of ARM aerosol and lidar data.
Corresponding Author: Adam Theisen (atheisen@anl.gov)

Overview¶

The first notebook in this series was an introduction to some of the features in ACT. In this notebook, we will be exploring different types and data and bringing it all together to get a view of this dust event.

Micropulse Lidar (MPL) Data
Aerodynamic Particle Sizer (APS) Data
Putting it All Together
Questions for the User to Explore

Prerequisites¶

This notebook will rely heavily on Python and the Atmospheric data Community Toolkit (ACT). Don’t worry if you don’t have experience with either, this notebook will walk you though what you need to know.

You will also need an account and token to download data using the ARM Live webservice. Navigate to the webserive information page and log in to get your token. Your account username will be your ARM username.

Concepts	Importance	Notes
ACT	Helpful

Time to learn: 60 Minutes
System requirements:
- Python 3.11 or latest
- ACT v1.5.0 or latest
- numpy
- xarray
- matplotlib

Imports¶

Let’s get started with some data! But first, we need to import some libraries.

import act
import numpy as np
import xarray as xr
import matplotlib.pyplot as plt

Micropulse Lidar (MPL) Data¶

ARM has many value-added products (VAPs) that apply corrections, additional quality control, retrievals, and more. These are a great way to get the data if you don’t want to do the added processing yourself. In this case, there is a VAP for the MPL with the datastream hou30smplcmask1zwangM1.c1 that you could utilize. For this case though, let’s download the instrument data and correct it ourselves using ACT.

# Set your username and token here!
username = 'YourUserName'
token = 'YourToken'

# Set the datastream and start/enddates
datastream = 'houmplpolfsM1.b1'
startdate = '2022-07-16'
enddate = '2022-07-16'

# Use ACT to easily download the data.  Watch for the data citation!  Show some support
# for ARM's instrument experts and cite their data if you use it in a publication
result = act.discovery.download_arm_data(username, token, datastream, startdate, enddate)

# Let's read in the data using ACT and check out the data
ds_mpl = act.io.read_arm_netcdf(result)

# Now we can correct the data
ds_mpl = act.corrections.correct_mpl(ds_mpl)

# Let's plot up the data to see what it looks like. But first,
# if you look at the variables, you would see the the variable
# we're going to plot has range_bins as it's 2nd dimension.
# We want it to be height so we have to swap some coordinates around
ds_mpl.coords['height'] = ds_mpl.height
ds_mpl = ds_mpl.swap_dims({'range_bins': 'height'})

# It's a lot of data, so it will take some time!
display = act.plotting.TimeSeriesDisplay(ds_mpl, figsize=(10, 6))
display.plot('signal_return_co_pol', cb_friendly=True, vmin=-20, vmax=20)

# This instrument collects data over 30 km so let's constrain
# this to just the lower 5 km
display.set_yrng([0, 5])
plt.show()

# Let's look at the cross polarization to get some more information
# It's a lot of data, so it will take some time!
display = act.plotting.TimeSeriesDisplay(ds_mpl, figsize=(10, 6))
display.plot('signal_return_cross_pol', cb_friendly=True, vmin=-20, vmax=20)

# This instrument collects data over 30 km so let's constrain
# this to just the lower 5 km
display.set_yrng([0, 5])
plt.show()

Aerosol Particle Size Data¶

Dust particles are generally under 10 µm or 10,000 nm and ARM has a variety of instruments that measure across different ranges. The chart below shows the ranges for each instrument and where we can expect to see the dust particles.

# Download the data as before
datastream = 'houaosapsM1.b1'
result = act.discovery.download_arm_data(username, token, datastream, startdate, enddate)

# and read it in
ds_aps = act.io.arm.read_arm_netcdf(result)

# Let's plot up the size distribution and number concentration
# These data are taken every second so it will take a little bit
display = act.plotting.TimeSeriesDisplay(ds_aps, subplot_shape=(2,), figsize=(10,8))
display.plot('total_N_conc', subplot_index=(0,))
display.day_night_background(subplot_index=(0,))
display.plot('dN_dlogDp', subplot_index=(1,), vmax=40)

plt.show()

Putting it All Together¶

Let’s start pulling these data together into the same plots so we can see what’s going on.

# Let's start fresh and get all the data and maybe expand the dates
startdate = '2022-07-10'
enddate = '2022-07-20'

# MPL - Let's try the VAP this time
datastream = 'hou30smplcmask1zwangM1.c1'
result = act.discovery.download_arm_data(username, token, datastream, startdate, enddate)
ds_mpl = act.io.arm.read_arm_netcdf(result)

# APS
datastream = 'houaosapsM1.b1'
result = act.discovery.download_arm_data(username, token, datastream, startdate, enddate)
ds_aps = act.io.arm.read_arm_netcdf(result)

# SMPS
datastream = 'houaossmpsM1.b1'
result = act.discovery.download_arm_data(username, token, datastream, startdate, enddate)
ds_smps = act.io.arm.read_arm_netcdf(result)

#ACSM
datastream = 'houaosacsmM1.b2'
result = act.discovery.download_arm_data(username, token, datastream, startdate, enddate)
ds_acsm = act.io.arm.read_arm_netcdf(result)

#PSAP
datastream = 'houaospsap3w1mM1.b1'
result = act.discovery.download_arm_data(username, token, datastream, startdate, enddate)
ds_psap = act.io.arm.read_arm_netcdf(result)

#SP2
datastream = 'houaossp2bc60sM1.b1'
result = act.discovery.download_arm_data(username, token, datastream, startdate, enddate)
ds_sp2 = act.io.arm.read_arm_netcdf(result)

# AOSMET
datastream = 'houmetM1.b1'
result = act.discovery.download_arm_data(username, token, datastream, startdate, enddate)
ds_met = act.io.arm.read_arm_netcdf(result)

# Before we proceed to plotting, let's reduce the MPL data down a little bit
# This will remove all data where heights are greater than 5
ds_mpl = ds_mpl.where(ds_mpl.height <= 5, drop=True)

# This will resample to 1 minute
ds_mpl = ds_mpl.resample(time='1min').nearest()

# Let's not forget about QCing the data!
# We can remove all the bad data from each aerosol dataset
ds_aps.clean.cleanup()
ds_aps = act.qc.arm.add_dqr_to_qc(ds_aps)
ds_aps.qcfilter.datafilter(rm_assessments=['Bad'], del_qc_var=False)

ds_smps.clean.cleanup()
ds_smps = act.qc.arm.add_dqr_to_qc(ds_smps)
ds_smps.qcfilter.datafilter(rm_assessments=['Bad'], del_qc_var=False)

ds_acsm.clean.cleanup()
ds_acsm = act.qc.arm.add_dqr_to_qc(ds_acsm)
ds_acsm.qcfilter.datafilter(rm_assessments=['Bad'], del_qc_var=False)

ds_sp2.clean.cleanup()
ds_sp2 = act.qc.arm.add_dqr_to_qc(ds_sp2)
ds_sp2.qcfilter.datafilter(rm_assessments=['Bad'], del_qc_var=False)

ds_psap.clean.cleanup()
ds_psap = act.qc.arm.add_dqr_to_qc(ds_psap)
ds_psap.qcfilter.datafilter(rm_assessments=['Bad'], del_qc_var=False)

# We can pass a dictionary to the display objects with multiple datasets
# So let's plot all this up!
display = act.plotting.TimeSeriesDisplay({'aps': ds_aps, 'mpl': ds_mpl, 'smps': ds_smps, 'acsm': ds_acsm, 'psap': ds_psap, 'sp2': ds_sp2, 'met': ds_met},
                                         subplot_shape=(7,), figsize=(10,18))

# MPL Plot
# Variable names of interest linear_depol_ratio, linear_depol_snr, backscatter_snr
display.plot('linear_depol_ratio', dsname='mpl', subplot_index=(0,), cb_friendly=True)
display.set_yrng([0, 3], subplot_index=(0,))

# APS Plot
display.plot('total_N_conc', dsname='aps', subplot_index=(1,))
display.day_night_background(dsname='aps', subplot_index=(1,))

# SMPS Plot
display.plot('total_N_conc', dsname='smps', subplot_index=(2,))
display.day_night_background(dsname='smps', subplot_index=(2,))

# ACSM plot
display.plot('sulfate', dsname='acsm', subplot_index=(3,), label='sulfate')
display.plot('nitrate', dsname='acsm', subplot_index=(3,), label='nitrate')
display.plot('ammonium', dsname='acsm', subplot_index=(3,), label='ammonium')
display.plot('chloride', dsname='acsm', subplot_index=(3,), label='chloride')
display.plot('total_organics', dsname='acsm', subplot_index=(3,), label='total_organics')

display.day_night_background(dsname='acsm', subplot_index=(3,))

# SP2 Plot
display.plot('sp2_rbc_conc', dsname='sp2', subplot_index=(4,))
display.day_night_background(dsname='sp2', subplot_index=(4,))

# PSAP Plot
display.plot('Ba_B_Weiss', dsname='psap', subplot_index=(5,))
display.day_night_background(dsname='psap', subplot_index=(5,))

# MET Wind Plot
display.plot('wdir_vec_mean', dsname='met', subplot_index=(6,), linestyle='', marker='.')
#display.plot('wspd_vec_mean', dsname='met', subplot_index=(6,), marker=None, secondary_y=True, color='orange', label='Wind Speed')
display.day_night_background(dsname='met', subplot_index=(6,))

plt.subplots_adjust(hspace=0.3)
plt.legend()
plt.savefig('./images/output.png')
plt.show()

Check this out¶

Floutsi et al 2023 has a great paper that shows how different aerosols can have different ratios. In Figure 3, you can see that Sharan dust has a slightly different signal than Central Asian dust.

Try more data!¶

Go back and adjust the start/end date around this period and see how things change! You can also try adding in some of the other instruments like the OPC (houaosopcM1.b1: calc_dust_weight or total_N_conc) or the NEPHDRY (houaosnephdry1mM1.b1: Bs_B_Dry_Neph3W). You can find a full list of data products on the ARM TRACER website

Advanced Visualizations and Processing¶

Let’s try and dive into the data a little more and see if there are any patterns in the data based on the direction of the wind.

# We already should have the data loaded up so let's explore with some data roses
# First we need to combine data and to do that, we need to get it on the same time grid
ds_combined = xr.merge([ds_met.resample(time='30min').nearest(), ds_acsm.resample(time='30min').nearest()], compat='override')

# Plot out the data rose using the WindRose display object
display = act.plotting.WindRoseDisplay(ds_combined)
display.plot_data('wdir_vec_mean', 'wspd_vec_mean', 'sulfate', num_dirs=15, plot_type='line', line_plot_calc='mean')
plt.show()

# First we need to combine data and to do that, we need to get it on the same time grid
ds_combined = xr.merge([ds_met.resample(time='1min').nearest(), ds_sp2.resample(time='1min').nearest()], compat='override')

# Plot out the data rose using the WindRose display object
display = act.plotting.WindRoseDisplay(ds_combined)
# Let's try a different type of data rose that will show the mean Black Carbon Concentration
# depending on wind direction and speed
display.plot_data('wdir_vec_mean', 'wspd_vec_mean', 'sp2_rbc_conc', num_dirs=15, plot_type='contour', contour_type='mean')
plt.show()

Checkout the area¶

The AMF was deployed at La Porte Municipal Airport. Check out the google map and see if this mapes sense!

Now back to the dust event!¶

Let’s see if we can look at the lowest layers in the MPL data and see how they compare to the APS data. For this, we can use simple matplotlib plots to investigate as well.

# Let's take a slice of the data to figure out where the signal
# may be showing up in the profile
data = ds_mpl['linear_depol_ratio'].values[:, 0:10]

fig, ax = plt.subplots()
ax.plot(ds_mpl['time'], data)
plt.show()

# Let's do it again and add in the APS data
# Other variables to try: backscatter_snr, linear_depol_snr, linear_depol_ratio
data = ds_mpl['linear_depol_ratio'].values[:, 4]

fig, ax = plt.subplots(figsize=(12,6))
ax.plot(ds_mpl['time'], data, label='MPL LDR')

ax2 = ax.twinx()
ax2.plot(ds_aps['time'], ds_aps['total_N_conc'], color='purple', label='APS')
plt.legend()
plt.show()

# Let's try a scatter plot using matplotlib
# First we need the data on the same time scale
ds_mpl2 = ds_mpl.resample(time='30s').nearest()
ds_aps2 = ds_aps.resample(time='30s').nearest()

ds = xr.merge([ds_mpl2, ds_aps2], compat='override')

# Other variables to try: backscatter_snr, linear_depol_snr, linear_depol_ratio
data = ds['linear_depol_snr'].values[:, 4]

fig, ax = plt.subplots(figsize=(10, 8))
ax.scatter(data, ds['total_N_conc']) # , c=ds['linear_depol_snr'].values[:, 4])

plt.show()

Excercises¶

Here are a couple additional excercises that users can work through on their own time.

Time series of ratios of fine particles vs. coarse particles using SMPS and APS data¶

Using the SMPS (ds_smps) and APS (ds_aps) data, can you calculate a ratio of fine to coarse particles? Remember, that this data would need to be on the same timestep.

# First, let's resample to 1-minute averages.  Note, that this may take some time with this much data.
# The other option is to grab the nearest value
ds_smps_1m = ds_smps.resample(time='1min').nearest()
ds_aps_1m = ds_aps.resample(time='1min').nearest()

# But wait, the concentration variables have the same name!  We need to rename
ds_aps_1m = ds_aps_1m.rename_vars({'total_N_conc': 'total_N_conc_aps'})
ds_merge = xr.merge([ds_smps_1m, ds_aps_1m], compat='override')

# Next, let's calculate the ratio
ratio = ds_merge['total_N_conc'].values / ds_merge['total_N_conc_aps'].values

# Then we can add it back to the Dataset as a variable
atts = {'units': '1', 'long_name': 'Ration of SMPS to APS total_N_conc'}
da = xr.DataArray(ratio, coords=ds_merge['total_N_conc'].coords, dims=ds_merge['total_N_conc'].dims, attrs=atts)
ds_merge['smps_aps_ratio'] = da

# And then plot it wil ACT!
display = act.plotting.TimeSeriesDisplay(ds_merge, figsize=(12, 10), subplot_shape=(3,))
display.plot('smps_aps_ratio', subplot_index=(0,))
display.plot('total_N_conc', subplot_index=(1,))
display.plot('total_N_conc_aps', subplot_index=(2,))
plt.show()

Compare Chemical Composition Between Events¶

Make a chart of aerosol composition on July 12, 13, and compare it to that from 17-18. This will require subsetting the ds_acsm data and plotting it using Matplotlib

# First subset the data for both sets
ds_acsm_12 = ds_acsm.sel(time=slice('2022-07-12', '2022-07-13'))
ds_acsm_17 = ds_acsm.sel(time=slice('2022-07-17', '2022-07-18'))

# Run through and get data
labels = ['total_organics', 'sulfate', 'nitrate', 'ammonium', 'chloride']
data12 = []
data17 = []
for label in labels:
    data12.append(float(ds_acsm_12[label].mean().values.tolist()))
    data17.append(float(ds_acsm_17[label].mean().values.tolist()))
    
# Pie charts use values that add up to 100% so we have to adjust the data
#data12 = data12 / np.sum(data12)
#data17 = data17 / np.sum(data17)

fig, ax = plt.subplots(1,1, figsize=(12,6))
ax.bar(range(len(labels)), data12, label='July 12-13')
plt.xticks(range(len(labels)), labels)

ax.bar(range(len(labels)), data17, label='July 17-18')
plt.xticks(range(len(labels)), labels)
plt.legend()
plt.show()

Additional Questions for the User to Explore¶

Can you improve these matplotlib plots to make them more readable?
Are there other variables that you would be interested in for the data roses?
We looked at relationships between the MPL and APS, are there other sets of instruments that might have relationships? (SP2 and PSAP?)

More Advanced¶

Can we use a machine learning toolkit like scikit-learn to cluster the MPL/APS data?
Could we get relationships between those clusters or a classification of what it is (smoke, dust, pollution)?

Next Steps¶

This notebook showed users how they can explore multiple different datasets using ACT. The next steps are for the users to continue exploring and diving into the data and the questions posed at the end!

Data Used in this Notebook¶

Ermold, B., & Flynn, C. Particle Soot Absorption Photometer (AOSPSAP3W1M). Atmospheric Radiation Measurement (ARM) User Facility. Sedlacek et al. (2011)

Jackson, R., & Sedlacek, A. Single Particle Soot Photometer (AOSSP2BC60S). Atmospheric Radiation Measurement (ARM) User Facility. Jackson et al. (2021)

Kuang, C., & Singh, A. Aerodynamic Particle Sizer (AOSAPS). Atmospheric Radiation Measurement (ARM) User Facility. Singh et al. (2021)

Kyrouac, J., & Shi, Y. Surface Meteorological Instrumentation (MET). Atmospheric Radiation Measurement (ARM) User Facility. Kyrouac et al. (2021)

Muradyan, P., Cromwell, E., Koontz, A., & Coulter, R. Micropulse Lidar (MPLPOLFS). Atmospheric Radiation Measurement (ARM) User Facility. Muradyan et al. (2010)

Sivaraman, C., Flynn, D., Riihimaki, L., & Comstock, J. Cloud mask from Micropulse Lidar (30SMPLCMASK1ZWANG). Atmospheric Radiation Measurement (ARM) User Facility. Sivaraman et al. (2019)

Zawadowicz, M., & Howie, J. Aerosol Chemical Speciation Monitor (AOSACSM). Atmospheric Radiation Measurement (ARM) User Facility. Zawadowicz et al. (2021)

References¶

Sedlacek, A., Smith, S., Ermold, B., Flynn, C., Trojanowski, R., Andrews, E., Hayes, C., Salwen, C., & Flynn, C. (2011). aospsap3w1m.b1. Atmospheric Radiation Measurement User Facility. 10.5439/1225037
Jackson, R., Sedlacek, A., Salwen, C., & Hayes, C. (2021). Single Particle Soot Photometer (SP2) Black Carbon mass and number concentrations. Atmospheric Radiation Measurement (ARM) Archive, Oak Ridge National Laboratory (ORNL), Oak Ridge, TN (US); ARM Data Center, Oak Ridge National Laboratory (ORNL), Oak Ridge, TN (United States). 10.5439/1807910
Singh, A., Oliveira, D. C. D., Koontz, A., Sedlacek, A., & Kuang, C. (2021). aosaps.b1. Atmospheric Radiation Measurement (ARM) Archive, Oak Ridge National Laboratory (ORNL), Oak Ridge, TN (US); ARM Data Center, Oak Ridge National Laboratory (ORNL), Oak Ridge, TN (United States). 10.5439/1407135
Kyrouac, J., Shi, Y., & Tuftedal, M. (2021). met.b1. Atmospheric Radiation Measurement (ARM) Archive, Oak Ridge National Laboratory (ORNL), Oak Ridge, TN (US); ARM Data Center, Oak Ridge National Laboratory (ORNL), Oak Ridge, TN (United States). 10.5439/1786358
Muradyan, P., Cromwell, E., Koontz, A., Coulter, R., Flynn, C., Ermold, B., & OBrien, J. (2010). mplpolfs.b1. Atmospheric Radiation Measurement (ARM) Archive, Oak Ridge National Laboratory (ORNL), Oak Ridge, TN (US); ARM Data Center, Oak Ridge National Laboratory (ORNL), Oak Ridge, TN (United States). 10.5439/1320657
Sivaraman, C., Flynn, D., Riihimaki, L., Comstock, J., & Zhang, D. (2019). Cloud Mask from Micropulse Lidar. Atmospheric Radiation Measurement (ARM) Archive, Oak Ridge National Laboratory (ORNL), Oak Ridge, TN (US); ARM Data Center, Oak Ridge National Laboratory (ORNL), Oak Ridge, TN (United States). 10.5439/1508389
Zawadowicz, M., Howie, J., Hayes, C., Allain, M., Salwen, C., & Behrens, B. (2021). Aerosol Chemical Speciation Monitor, mentor processed, b2. Atmospheric Radiation Measurement (ARM) Archive, Oak Ridge National Laboratory (ORNL), Oak Ridge, TN (US); ARM Data Center, Oak Ridge National Laboratory (ORNL), Oak Ridge, TN (United States). 10.5439/1762267