Example on retrieving and plotting windsΒΆ

This is a simple example for how to retrieve and plot winds from 2 radars using PyDDA.

Author: Robert C. Jackson

  • ../../_images/sphx_glr_plot_examples_001.png
  • ../../_images/sphx_glr_plot_examples_002.png
  • ../../_images/sphx_glr_plot_examples_003.png

Out:

Calculating weights for radars 0 and 1
Calculating weights for models...
Starting solver
rmsVR = 6.776753068494562
Total points:92092.0
| Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel | Max w
|  22.8056| 179.8659|   0.0000|   0.0000|   0.0000|   0.0000|  14.2882
Norm of gradient: 0.060347431274420944
Iterations before filter: 10
| Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel | Max w
|   3.8713|  66.1170|   0.0000|   0.0000|   0.0000|   0.0000|  23.7745
Norm of gradient: 0.02372694747570922
Iterations before filter: 20
| Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel | Max w
|   4.0472|  37.8079|   0.0000|   0.0000|   0.0000|   0.0000|  24.6297
Norm of gradient: 0.02779702290042135
Iterations before filter: 30
| Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel | Max w
|   0.8673|  28.6315|   0.0000|   0.0000|   0.0000|   0.0000|  26.5968
Norm of gradient: 0.006835527283244816
Iterations before filter: 40
| Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel | Max w
|   1.6686|  20.8826|   0.0000|   0.0000|   0.0000|   0.0000|  28.4437
Norm of gradient: 0.023224596559256867
Iterations before filter: 50
| Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel | Max w
|   0.5065|  18.2438|   0.0000|   0.0000|   0.0000|   0.0000|  28.6566
Norm of gradient: 0.006565409425024461
Iterations before filter: 60
| Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel | Max w
|   0.5679|  15.6097|   0.0000|   0.0000|   0.0000|   0.0000|  29.2501
Norm of gradient: 0.006177474586057
Iterations before filter: 70
| Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel | Max w
|   0.3157|  14.5197|   0.0000|   0.0000|   0.0000|   0.0000|  30.2715
Norm of gradient: 0.0035443644107431118
Iterations before filter: 80
| Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel | Max w
|   0.2528|  13.5275|   0.0000|   0.0000|   0.0000|   0.0000|  30.3527
Norm of gradient: 0.002892126199622085
Iterations before filter: 90
| Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel | Max w
|   0.2858|  13.0896|   0.0000|   0.0000|   0.0000|   0.0000|  30.6199
Norm of gradient: 0.0033178819236652727
Iterations before filter: 100
| Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel | Max w
|   0.2688|  13.0232|   0.0000|   0.0000|   0.0000|   0.0000|  30.6529
Norm of gradient: 0.003424754964509097
Iterations before filter: 110
| Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel | Max w
|   0.2203|  13.0222|   0.0000|   0.0000|   0.0000|   0.0000|  30.6752
Norm of gradient: 0.003483670630141192
Iterations before filter: 120
| Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel | Max w
|   0.2194|  13.0228|   0.0000|   0.0000|   0.0000|   0.0000|  30.6785
Norm of gradient: 0.0034925009878411786
Iterations before filter: 130
Applying low pass filter to wind field...
Iterations after filter: 1
Iterations after filter: 2
Done! Time = 540.5

import pyart
import pydda
from matplotlib import pyplot as plt


berr_grid = pyart.io.read_grid(pydda.tests.EXAMPLE_RADAR0)
cpol_grid = pyart.io.read_grid(pydda.tests.EXAMPLE_RADAR1)

sounding = pyart.io.read_arm_sonde(
    pydda.tests.SOUNDING_PATH)


# Load sounding data and insert as an intialization
u_init, v_init, w_init = pydda.initialization.make_wind_field_from_profile(
        cpol_grid, sounding[1], vel_field='corrected_velocity')

# Start the wind retrieval. This example only uses the mass continuity
# and data weighting constraints.
Grids = pydda.retrieval.get_dd_wind_field([berr_grid, cpol_grid], u_init,
                                          v_init, w_init, Co=1.0, Cm=1500.0,
                                          Cz=0,
                                          frz=5000.0, filt_iterations=2,
                                          mask_outside_opt=True, upper_bc=1)
# Plot a horizontal cross section
plt.figure(figsize=(9, 9))
pydda.vis.plot_horiz_xsection_barbs(Grids, background_field='reflectivity',
                                    level=6,
                                    w_vel_contours=[3, 6, 9, 12, 15],
                                    barb_spacing_x_km=5.0,
                                    barb_spacing_y_km=15.0)
plt.show()

# Plot a vertical X-Z cross section
plt.figure(figsize=(9, 9))
pydda.vis.plot_xz_xsection_barbs(Grids, background_field='reflectivity',
                                 level=40,
                                 w_vel_contours=[3, 6, 9, 12, 15],
                                 barb_spacing_x_km=10.0,
                                 barb_spacing_z_km=2.0)
plt.show()

# Plot a vertical Y-Z cross section
plt.figure(figsize=(9, 9))
pydda.vis.plot_yz_xsection_barbs(Grids, background_field='reflectivity',
                                 level=40,
                                 w_vel_contours=[3, 6, 9, 12, 15],
                                 barb_spacing_y_km=10.0,
                                 barb_spacing_z_km=2.0)
plt.show()

Total running time of the script: ( 9 minutes 2.931 seconds)

Gallery generated by Sphinx-Gallery