Bring your laptops to class.  Reagan McKinney leads us in the hysplit tutorials.

1. Start by viewing the first item linked in the Resources section of the Hysplit web page for this course.
2. Next, Read the two items listed in the Theory section of the Hysplit web page for this course.
   
3. Then read the BAMS journal paper about Hysplit, which is the third link in the Resources section of the Hysplit web page for this course.  Also, read the Supplement to the BAMS paper.

1. Present to the class on 26 Feb a section of the BAMS Supplement, as was selected during a previous class meeting.   Use Piazza to vote on which topics you would like to cover (please pick 3).
1. Mixing Options, first paragraph (including eqs. ES7-ES8)
2. Mixing Options - stable & neutral, second paragraph (including eqs. ES9 - ES14)
3. Mixing Options - unstable, second paragraph continued (including eqs. ES15 - ES17)
4. Mixing Options - TKE, third paragraph (including eqs. ES18-ES21)
5. Mixing Options - TKE, third paragraph continued (including eqs. ES22-ES25)
6. Stability & Mixing Depth Calculation Options
7. Dynamic Management...for Particles and Puffs.
8. Time-varying Emissions (first half) - Wind-Blown Dust
9. Time-varying Emissions (last portion, on Plume Rise)
10. Deep Convection / CAPE
11. In-Cloud Wet Scavenging
    

1. Be prepared to discuss how mean wind and turbulence are included in Lagrangian models for atmospheric dispersion.
2. Be prepared to discuss how Lagrangian models for dispersion in turbulent atmospheres differ from the Langevin model for molecular diffusion.
3. Be prepared to discuss the boundary layer and the dispersion parameterizations in the BAMS supplement.

 Stull and/or Tim covers Panoply installation and use. 

Topics for in-class exercises: Get NetCDF files and view them using NASA's Panoply program.  Convert WRF meteorology NetCDF files into ARL format that Hysplit can read.  Then use those ARL files to run your own simulated emissions from UBC.

1. In class, download to your laptop the zipped WRF output files , and unzip it to give a folder of NetCDF files.
2. Add the suffix .nc to each of those filenames.   (This identifies them as NetCDF files.)
3. From the NASA Panoply website ( https://www.giss.nasa.gov/tools/panoply/ ) Mac users must use the Safari browser (not Chrome) to:
1. first download to your laptop the updated version of Java runtime environment, follow the instructions to install Java, and to
   
2. use Safari/preferences/websites to activate (using the check box) the Java plug-in, and to allow the giss.nasa.gov website to use this plug-in. 
   
3. Finally download the panoply application to your applications folder.
   
4. Use Panoply to view the contents of any one of the NetCDF files and to view maps of different fields at different heights.  This is to get you familiar with the type of info that is in NetCDF files.   (Hint: to zoom any map, on the Mac use command drag to create an outline of the region you want to zoom.)
5. Convert the WRF meteorological files from NetCDF to ARL format.  If you are unable to do that conversion, then use an already converted ARL file (for a different date) provided in the next item below.
   
6. If you are unable to convert the NetCDF files to ARL files yourself, then instead you can download the large (1.5 GB)  ARL meteorological file we provide for Canada (at 12 km grid spacing) to your laptop, to calculate  trajectories and concentrations for emissions from a hypothetical source as specified in class.  (Note: it might take a few minutes to download the large data file.) 
7. Use the Hysplit gui tools to open and view the weather patterns in the ARL file.
8. Create emission from UBC at different heights to track using the trajectory option.  Also try emissions at one height to show concentrations.  Also try to view a cloud of particles emitted from one source.

Here is the Hysplit dispersion-model run assignment, which you can start in class, and finish at home to submit via Canvas.

(If there is extra time this period, Start can discussing how Calpuff works. Stull: Intro to Calpuff "slugs".  Review of sections 2.1 and 2.2 of Calpuff Users Guide v5. ) 

The following homework is due today. (But OK to turn in late.)

# HW4 - Homework, Write your own 1-D Lagrangian particle dispersion model
# R. Stull, 22 Feb 2016, 10 Feb 2024
# ES = equation number in BAMS supplementary material by Stein et al (Dec 2015) on Hysplit model
# https://journals.ametsoc.org/supplemental//journals/bams/96/12/bams-d-14-00110.1.xml/10.1175_bams-d-14-00110.2.pdf

# Givens:
# general
g_by_tv = 0.0333     # g/Tv approx constant (m s^-2 K^-1)
rho_cp = 1231.         # air density times specific heat at const press, approx as constant. (W/m2)/(K m/s)

# unstable pbl
Hsfc = 280.       # surface heat flux (W/m2)
FHsfc = Hsfc / rho_cp       # kinematic heat flux (K m/s) at the surface
zi = 1000.         # mixed-layer depth for unstable pbl (m)
tlvert = 200.      # lagrangian time scale for vertical motions for unstable pbl (s)

# plume info
nparticles = 20    # number of particles (emitted as a single burst)
zs = 100.             # source emission height (m) for all particles

# initial conditions
t = 0.       # initial time
it = 0.      # initial time index
zp = zs     # initial particle height for any one particle (same for all particles)
wp = 0.     # initial W' vertical velocity value for that same one particle (same for all particles)

# model integration
delt = 5.         # timestep increment (s) for each time step
itmax = 200     # max number of timesteps

# Hint 1: for unstable pbl, use w* = Deardorff velocity from Stull Practical Meteor, eq (18.19a)
# Hint 2: kinematic heat flux FHsfc = (traditional heat flux) / rho_cp .
# Hint 3: include particle reflection at the top and bottom boundaries

# Homework, create a Lagrangian model to track the vertical movement of a bunch particles released
# simultaneously from the same source height into a convective boundary layer. Consider only vertical
# particle movement, and plot the vertical position of all the particles as a function of timestep.
# Use a vertical domain of 0 ≤ z ≤ zi, where zi is the mixed layer depth.
# No mean wind, thus no horizontal domain. Instead, the horizontal axis in your plot is timestep index.
#
# Use the sigma_w expression given by Stein el al (2015) eq. ES15 in their supplement, to find the
# the value of sigma_w at the particle location.
# Use the position forecast eq (simplified Draxler eq 53). Z(t+?t) = Z(t) + W' * ?t
#
# Make two different experiments, with slightly different eqs for forecasting W' :
#
# 1) Here is a simplified version of Stein et al Supplement eq. (ES4), but without the
# sigma-gradient term (with notation the same as in that journal paper):
# wp = (R*wp) + (sigma_w * lambda * sqrt(1 - R^2) )
#
# 2) Same as (1), but with the additional sigma_w gradient term dsigma_w/dz included:
# wp = (R*wp) + (sigma_w * lambda * sqrt(1 - R^2) ) + (sigma_w * tlvert * (1-R) * dsigma_w/dz )
#
# 3) Discuss the difference between your plots from (1) and (2) above. Why is it important?
#
# where
# tlvert = Lagrangian time scale for vertical motion,
# R = autocorrelation coefficient,
# lambda is a random number from a Gaussian distribution with mean of zero and st.dev of 1
#
# Hint 4: Caution calculating the vertical sigma_w gradient at the top & bottom domain boundaries.
#
# Check your answers and discuss the significance. 
# Turn in your code and the contour plots.

If needed, Tim Chui will finish showing us how to install Calpuff, Calmet, Calpost, and other components.

Before class, please skim the Calpuff User Instructions.

Meanwhile, start working on the following homework, due 6 March:

=========

# HW5 - Homework:  2-D  (x-z) Lagrangian particle transport.  Use any language: R, Matlab, python, etc.
# Hint: You might be able to re-use parts of the code from your 1-D Lagrangian homework.

# A) First: create an artificial mean-wind environment where the mean U wind varies in x and z,
#       within which you can test your 1-D (x-direction) mean advection by the U component of wind. 
#      Assume V = 0 , W = 0.  And neglect turbulence in the x & y-directions (i.e., u' = 0, v' = 0 ).
#      For this wind field, let U obey the log wind profile with z, and that u* varies with x.
# Given:
ustar0 = 0.2     # first u* constant (m/s) in equation below
ustar1 = 0.5     # second u* constant (m/s) in equation below
zo = 0.2          # roughness length (m)
L = 50000.      # horizontal domain (m).   0 ≤ x ≤ 50 km
zi = 1000.        # (m) your veritical domain is 0 ≤ z ≤ zi
#
# Hint: create functions for the following equations (where L = domain size, NOT Obukhov length):
# Use u* as a function of x:  ustar = ustar0 + ustar1*( (x/L) * (1. - (x/L)) )
# Use U as a function of z and u*:  U = (ustar/0.4) * log( (z+zo)/zo )
#
# To see what you produced, calculate U(x,z) for fairly high resolution in x, and z,
# and then plot it as a contour plot within the domain given above. 

#

# ===========
# B) Implement the Hysplit mean-wind advection scheme in your own code, for just the U wind.
#
Implement the Hypslit mean advection scheme, and test it for a single particle released
  at (x, z) = (0, 100. m).   Namely the source height is zs = 100.
  This scheme is given by eq. (1) in the Stein et al paper (2015) in BAMS.  (See our course web page.)
  Namely, neglecting vertical turbulence for now, calculate the postion of the one particle at each
  time step (dt = 60 s) as it advects to the East.
  Although you do NOT need to submit this, you might want to plot a point at the position of the
  particle at each timestep as it advects across the (x, z) domain,
  to convince yourself that your horizontal advection is working.
  Your answer should be a straight line of dots that are spaced farther apart in regions of faster wind.


===========
C) Add the vertical-turbulent dispersion code from your previous 1-D Lagrangian homework
to the horizontal-advection code that you just created , so we can track 5 particles. 
Use your code that included reflection from z=0 and z=zi, and which prevented accumulation of
particles near the top and bottom edges.

# Givens:
# general
g_by_tv = 0.0333    # g/Tv approx constant  (m s^-2 K^-1)
rho_cp = 1231.      # air density x specific heat at const press, is constant. (W/m2)/(K m/s)

# unstable pbl
Hsfc = 50.                # surface heat flux (W/m2)
FHsfc = Hsfc / rho_cp     # kinematic heat flux (K m/s) at the surface
wstar = (g_by_tv * FHsfc * zi)^(1/3)          # Deardorff convective velocity
tlvert = 200.             # lagrangian time scale for vertical motions for unstable pbl (s)

dt = 60.                  # timestep increment (s) for tracking any one particle

# plume info
nparticles = 5           # number of particles  (emitted as a single burst at t = 0)
zs = 100.                 # source emission height (m) for all particles

Plot the tracks of these particles in your (x, z) domain, until they exit the domain.
Use the following special colour scheme.  First 10 positions use gray.
Next 2 positions, use red.  Next 10 positions gray.  Next 2 positions, use orange.
Next 10 gray.  Next 2 gold.  next 10 gray.  next 2 green.  next 10 gray.  next 2 blue.
Next 10 gray.  Next 2 violet.
Repeat the colour sequence (grays and rainbow colours)
as needed to cover the whole track of the particle until it exits the domain. Except plot two 
points as BLACK when the first of these 2 points is t = 120 time steps from their release time. 


========
D) Synthesis and Discussion

1) Even though all the particles were released at the same time, they obviously dispersed in the vertical (z).
Did they also disperse in the horizontal (x)? Why?

2) How could you covert positions of many particles into concentration values at any fixed location?

=== end ===

1. Start by viewing the first item linked in the Resources section of the Hysplit web page for this course.
2. Next, Read the two items listed in the Theory section of the Hysplit web page for this course.
   
3. Then read the BAMS journal paper about Hysplit, which is the second link in the Resources section of the Hysplit web page for this course.  Also, read the Supplement to the BAMS paper.
   

1. Register with NOAA ARL to gain access to the Hysplit model.
2. See the instructions on Hysplit web page  to download the tcl software, used as a graphic user interface for Hysplit on the Mac.
   
3. Follow our TA's instructions on our Hysplit web page to download a pre-compiled version of Hysplit for the Mac.  If you have problems, perhaps you classmates or our TA can help.
4. Run the Hysplit4 , Trajectory , QuickStart , Run Example, as outlined in the Tips by our TA.  If you did not get the same output as in the Tips, consult our TA.
   
5. Run the Hysplit4 , Concentration, QuickStart, Run Example, as outlined in the Tips by our TA.  If you did not get the same output as in the Tips, consult our TA.
6. From the online Hysplit tutorial ( https://ready.arl.noaa.gov/documents/Tutorial/html/index.html ), read sections 2.1 thru 2.4 (and view the "animation" linked from section 2.4) .
   

1. Be prepared to discuss how mean wind and turbulence are included in Lagrangian models for atmospheric dispersion.
2. Be prepared to discuss how Lagrangian models for dispersion in turbulent atmospheres differ from the Langevin model for molecular diffusion.
3. Be prepared to discuss the boundary layer and the dispersion parameterizations in the BAMS supplement.
   

1. On your own computer, navigate to your hysplit4 directory, and create a new subdirectory called "zcaptex". 
   
2. Into zcaptex, copy the following meteorology binary files from the Hysplit tutorial lesson 3.1 .
• North American Regional Reanalysis captex2_narr.bin 
  
• Weather Research Forecast Model (27 km) captex2_wrf27.bin
  
• Weather Research Forecast Model (9 km) captex2_wrf09.bin
  
• Weather Research Forecast Model (3 km) captex2_wrf03.bin
  
• NCAR Global Reanalysis Extract-1 captex2_gblr.bin
  
• NCAR Global Reanalysis Extract-2 RP198309.gbl
  
• ECMWF ERA-40 Reanalysis captex2_era40.bin 
  
3. In Lesson 3.1, item 1, follow the tutorial instructions click on the Set File Name button to the first file above, then click on the Create Map button to view the domain.
   Repeat for all the files above to view their domains.
4. In Lesson 3.1, item 2, follow the tutorial instructions to see details about each of those files.  If the resulting Simulation Log window is blank, just click on it or scroll in it to make the info appear.  Confirm that the first dozen or so lines in this file give info that agrees with the summary in the tutorial.
5. In Hysplit4, create a new subdirectory "zreanalysis".  Do tutorial  Lesson 3.2, item 1, by accessing Meteorology / ARL Data ftp / Reanalysis, to download that file.  It will probably end up in your Hysplit4/working directory, after 20 minutes or so for the download.   (The Download progress bar is flakey on the Mac, but if you click on it and then click off of it, the progress temporarily appears.  Repeat on and off clicks to make it re-appear.)
   When the download is finished, you can move it into the zreanalysis subdirectory for future use. 
   
6. In Hysplit4, create a new subdirectory "zwfrt".  Follow our TA's Tips to Download WRF from our own ensemble forecast system.   There will be 2 netCDF files to download, both for our WRF(ARW) GFS runs:  36 km (7.5 days), and 1.3 km (3.5 days).   This relates to tutorial lesson 3.3 on converting data.  Pick a day for which the initial weather has boundary-layer and mid-level winds predominantly from the west (anywhere between southwest to northwest), so that we can track emissions from UBC as they move through the domain.
   After acquiring these netCDF files, follow the directions in tutorial lesson 3.3 to convert them to the ARL format that Hysplit can read as input.  We will use these later.
7. Read only:  tutorial lessons 3.4 - 3.6. 
   We will come back to hysplit later.  Now it is time to write your own code to explore some of the parameterizations used in hysplit (see next section).