.. only:: html

    .. note::
        :class: sphx-glr-download-link-note

        Click :ref:`here <sphx_glr_download_source_auto_examples_plot_examples_nested.py>`     to download the full example code
    .. rst-class:: sphx-glr-example-title

    .. _sphx_glr_source_auto_examples_plot_examples_nested.py:


Example on retrieving and plotting winds on a distributed cluster
-----------------------------------------------------------------

This is a simple example for how to retrieve winds using the
nested grid features of PyDDA.

Author: Robert C. Jackson



.. rst-class:: sphx-glr-horizontal


    *

      .. image:: /source/auto_examples/images/sphx_glr_plot_examples_nested_001.png
          :alt: PyDDA retreived winds @4.05 km
          :class: sphx-glr-multi-img

    *

      .. image:: /source/auto_examples/images/sphx_glr_plot_examples_nested_002.png
          :alt: PyDDA retreived winds @20.0 km south of origin.
          :class: sphx-glr-multi-img

    *

      .. image:: /source/auto_examples/images/sphx_glr_plot_examples_nested_003.png
          :alt: PyDDA retreived winds @20.0 km west of origin.
          :class: sphx-glr-multi-img


.. rst-class:: sphx-glr-script-out

 Out:

 .. code-block:: none

    /blues/gpfs/home/rjackson/pyart/pyart/io/cfradial.py:365: UserWarning: WARNING: valid_min not used since it
    cannot be safely cast to variable data type
      data = self.ncvar[:]
    /blues/gpfs/home/rjackson/pyart/pyart/io/cfradial.py:365: UserWarning: WARNING: valid_max not used since it
    cannot be safely cast to variable data type
      data = self.ncvar[:]
    LocalCluster('tcp://127.0.0.1:34918', workers=2, threads=36, memory=135.13 GB)
    <Client: 'tcp://127.0.0.1:34918' processes=2 threads=36, memory=135.13 GB>
    /home/rjackson/anaconda3/envs/pyart-2020/lib/python3.7/site-packages/pydda/retrieval/angles.py:24: RuntimeWarning: invalid value encountered in arccos
      elev = np.arccos((Re**2 + slantrsq - rh**2)/(2 * Re * slantr))
    Calculating weights for radars 0 and 1
    /home/rjackson/anaconda3/envs/pyart-2020/lib/python3.7/site-packages/pydda/retrieval/wind_retrieve.py:653: RuntimeWarning: invalid value encountered in arccos
      theta_2 = np.arccos((x-rad2[1])/b)
    Calculating weights for models...
    Starting solver 
    rmsVR = 6.763780749907824
    Total points: 23059
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   7.3871|  31.7975|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  11.4149
    Norm of gradient: 0.06526775209383506
    Iterations before filter: 10
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   0.7160|  14.2379|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  25.1676
    Norm of gradient: 0.022656970042356808
    Iterations before filter: 20
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   0.5321|   9.8464|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  47.6630
    Norm of gradient: 0.019575025191868944
    Iterations before filter: 30
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   0.1613|   7.4489|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  54.7406
    Norm of gradient: 0.024242219942816612
    Iterations before filter: 40
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   0.0950|   6.6185|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  55.4003
    Norm of gradient: 0.006852171865775456
    Iterations before filter: 50
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   0.0558|   5.5634|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  48.4985
    Norm of gradient: 0.007062978135599922
    Iterations before filter: 60
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   0.0635|   5.1455|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  47.1794
    Norm of gradient: 0.007387512143371025
    Iterations before filter: 70
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   0.0743|   4.7404|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  41.4687
    Norm of gradient: 0.012528253916469778
    Iterations before filter: 80
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   0.0963|   4.4825|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  38.0718
    Norm of gradient: 0.007383190140180377
    Iterations before filter: 90
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   0.0248|   4.2487|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  35.5239
    Norm of gradient: 0.0063945697721171195
    Iterations before filter: 100
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   0.0282|   4.2277|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  35.1107
    Norm of gradient: 0.006355652235950298
    Iterations before filter: 110
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   0.0161|   4.2245|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  35.0654
    Norm of gradient: 0.006053121760374165
    Iterations before filter: 120
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   0.0212|   4.1906|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  34.8727
    Norm of gradient: 0.011452701382835957
    Iterations before filter: 130
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   0.0146|   4.1923|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  34.8769
    Norm of gradient: 0.0030135879499404877
    Iterations before filter: 140
    Applying low pass filter to wind field...
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |3746.0950|   2.8641|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  13.5970
    Norm of gradient: 0.998609813914212
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |3378.9225|   2.9509|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  13.5978
    Norm of gradient: 0.9547921379803676
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    | 221.0979|  23.0871|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  13.6217
    Norm of gradient: 0.23792226768098468
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    | 114.3039|  24.1503|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  13.8287
    Norm of gradient: 0.15061385334713479
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |  21.2646|  23.0345|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  13.6886
    Norm of gradient: 0.06481576024485435
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |  24.7439|  19.7879|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  13.7359
    Norm of gradient: 0.140284224757774
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |  10.5934|  21.2488|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  13.7121
    Norm of gradient: 0.045492655928316694
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   6.0063|  19.7451|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  13.7295
    Norm of gradient: 0.03183975927355813
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   1.9115|  16.4765|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  13.7880
    Norm of gradient: 0.017344264451652755
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   1.3409|  14.3912|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  13.8497
    Norm of gradient: 0.02718263559465569
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   1.4336|  12.4225|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  13.9390
    Norm of gradient: 0.015722180537225953
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   0.8802|  11.7888|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  13.9765
    Norm of gradient: 0.010982739882723191
    Iterations after filter: 1
    Iterations after filter: 2
    Done! Time = 138.8
    Waiting for nested grid to be retrieved...
    /home/rjackson/anaconda3/envs/pyart-2020/lib/python3.7/site-packages/pydda/vis/barb_plot.py:175: UserWarning: linewidths is ignored by contourf
      alpha=contour_alpha)
    /home/rjackson/anaconda3/envs/pyart-2020/lib/python3.7/site-packages/pydda/retrieval/wind_retrieve.py:653: RuntimeWarning: invalid value encountered in arccos
      theta_2 = np.arccos((x-rad2[1])/b)
    /home/rjackson/anaconda3/envs/pyart-2020/lib/python3.7/site-packages/pydda/vis/barb_plot.py:214: UserWarning: The following kwargs were not used by contour: 'color'
      levels=[bca_min, bca_max], color='k')
    /home/rjackson/anaconda3/envs/pyart-2020/lib/python3.7/site-packages/pydda/retrieval/wind_retrieve.py:653: RuntimeWarning: invalid value encountered in arccos
      theta_2 = np.arccos((x-rad2[1])/b)
    /home/rjackson/anaconda3/envs/pyart-2020/lib/python3.7/site-packages/pydda/vis/barb_plot.py:214: UserWarning: The following kwargs were not used by contour: 'color'
      levels=[bca_min, bca_max], color='k')
    /home/rjackson/anaconda3/envs/pyart-2020/lib/python3.7/site-packages/pydda/vis/barb_plot.py:637: UserWarning: linewidths is ignored by contourf
      alpha=contour_alpha)
    /home/rjackson/anaconda3/envs/pyart-2020/lib/python3.7/site-packages/pydda/vis/barb_plot.py:825: UserWarning: linewidths is ignored by contourf
      alpha=contour_alpha)






|


.. code-block:: default


    import pyart
    import pydda
    from matplotlib import pyplot as plt
    from distributed import LocalCluster, Client

    # Needed so that distributed doesn't run all of your code when the worker 
    # starts!
    if __name__ == '__main__':

        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 our dask distributed cluster. You can use any distributed cluster
        # for this...a LocalCluster is used here for the sake of being able to run
        # this example locally.
        cluster = LocalCluster(n_workers=2)
        print(cluster)
        client = Client(cluster)
        print(client)

        # Start the wind retrieval. This example only uses the mass continuity
        # and data weighting constraints.
        Grids = pydda.retrieval.get_dd_wind_field_nested(
            [berr_grid, cpol_grid], u_init,  v_init, w_init, client, 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()


.. rst-class:: sphx-glr-timing

   **Total running time of the script:** ( 4 minutes  47.999 seconds)


.. _sphx_glr_download_source_auto_examples_plot_examples_nested.py:


.. only :: html

 .. container:: sphx-glr-footer
    :class: sphx-glr-footer-example



  .. container:: sphx-glr-download sphx-glr-download-python

     :download:`Download Python source code: plot_examples_nested.py <plot_examples_nested.py>`



  .. container:: sphx-glr-download sphx-glr-download-jupyter

     :download:`Download Jupyter notebook: plot_examples_nested.ipynb <plot_examples_nested.ipynb>`


.. only:: html

 .. rst-class:: sphx-glr-signature

    `Gallery generated by Sphinx-Gallery <https://sphinx-gallery.github.io>`_
