.. only:: html

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

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

    .. _sphx_glr_source_auto_examples_plot_fun_with_constraints.py:


Example on geographic plotting and constraint variation
-------------------------------------------------------

In this example we show how to plot wind fields on a map and change
the default constraint coefficients using PyDDA.

This shows how important it is to have the proper intitial state and
constraints when you derive your wind fields. In the first figure,
the sounding was used as the initial state, but for the latter
two examples we use a zero initial state which provides for more 
questionable winds at the edges of the Dual Doppler Lobes.



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


    *

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

    *

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

    *

      .. image:: /source/auto_examples/images/sphx_glr_plot_fun_with_constraints_003.png
          :alt: PyDDA retreived winds @4.54 km
          :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[:]
    /home/rjackson/anaconda3/envs/pyart-2020/lib/python3.7/site-packages/pydda/vis/streamline_plot.py:396: RuntimeWarning: invalid value encountered in less
      w_filled = np.ma.masked_where(w[level, :, :] < np.min(w_vel_contours),
    /home/rjackson/anaconda3/envs/pyart-2020/lib/python3.7/site-packages/cartopy/mpl/geoaxes.py:1400: UserWarning: linewidths is ignored by contourf
      result = matplotlib.axes.Axes.contourf(self, *args, **kwargs)
    /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/cartopy/mpl/geoaxes.py:1366: UserWarning: The following kwargs were not used by contour: 'color'
      result = matplotlib.axes.Axes.contour(self, *args, **kwargs)
    /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/cartopy/mpl/geoaxes.py:1366: UserWarning: The following kwargs were not used by contour: 'color'
      result = matplotlib.axes.Axes.contour(self, *args, **kwargs)
    /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.776753068494561
    Total points: 92092
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |  17.5234| 117.2994|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  11.8763
    Norm of gradient: 0.03363053426960669
    Iterations before filter: 10
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   3.3593|  44.3660|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  18.9330
    Norm of gradient: 0.03594766688047328
    Iterations before filter: 20
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   1.6825|  26.6090|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  20.4340
    Norm of gradient: 0.019569191810207776
    Iterations before filter: 30
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   0.9598|  18.2610|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  23.1350
    Norm of gradient: 0.011461364156613545
    Iterations before filter: 40
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   0.5382|  14.2669|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  25.8571
    Norm of gradient: 0.005058645720581414
    Iterations before filter: 50
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   0.5519|  11.2071|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  27.7493
    Norm of gradient: 0.00913384867304721
    Iterations before filter: 60
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   0.5063|   9.0763|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  29.5937
    Norm of gradient: 0.004085182612931341
    Iterations before filter: 70
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   0.3574|   7.4804|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  30.4240
    Norm of gradient: 0.010016895672519415
    Iterations before filter: 80
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   0.2502|   7.4561|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  30.4265
    Norm of gradient: 0.0009764931617796537
    Iterations before filter: 90
    Applying low pass filter to wind field...
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |6575.6187|   6.1340|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|   8.7456
    Norm of gradient: 0.7462737969222745
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |5966.0414|   6.4019|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|   8.6014
    Norm of gradient: 0.6983294931210702
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    | 396.1297|  84.0446|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  11.6213
    Norm of gradient: 0.20222520819243076
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    | 210.3138|  87.8755|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  11.6785
    Norm of gradient: 0.12040101890762513
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |  46.4198|  87.9840|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  11.2879
    Norm of gradient: 0.08702823439023104
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |  47.4066|  74.8731|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  11.3296
    Norm of gradient: 0.1633964170150081
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |  13.2350|  73.1213|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  11.3530
    Norm of gradient: 0.03300072469301934
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   8.7393|  68.0998|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  11.3983
    Norm of gradient: 0.01919300030010946
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   6.1531|  54.8219|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  11.6947
    Norm of gradient: 0.028416505852986308
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   5.4481|  44.6685|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  12.4289
    Norm of gradient: 0.019417630225052232
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |  16.3036|  33.2967|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  14.2565
    Norm of gradient: 0.07539028471955857
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   6.9822|  36.7060|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  13.3551
    Norm of gradient: 0.04778436338484203
    Iterations after filter: 1
    Iterations after filter: 2
    Done! Time = 117.2
    /home/rjackson/anaconda3/envs/pyart-2020/lib/python3.7/site-packages/pydda/vis/streamline_plot.py:396: RuntimeWarning: invalid value encountered in less
      w_filled = np.ma.masked_where(w[level, :, :] < np.min(w_vel_contours),
    /home/rjackson/anaconda3/envs/pyart-2020/lib/python3.7/site-packages/cartopy/mpl/geoaxes.py:1400: UserWarning: linewidths is ignored by contourf
      result = matplotlib.axes.Axes.contourf(self, *args, **kwargs)
    /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/cartopy/mpl/geoaxes.py:1366: UserWarning: The following kwargs were not used by contour: 'color'
      result = matplotlib.axes.Axes.contour(self, *args, **kwargs)
    /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/cartopy/mpl/geoaxes.py:1366: UserWarning: The following kwargs were not used by contour: 'color'
      result = matplotlib.axes.Axes.contour(self, *args, **kwargs)
    /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.776753068494561
    Total points: 92092
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |  17.5234| 117.2994|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  11.8763
    Norm of gradient: 0.03363053426960669
    Iterations before filter: 10
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   3.3593|  44.3660|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  18.9330
    Norm of gradient: 0.03594766688047328
    Iterations before filter: 20
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   1.6825|  26.6090|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  20.4340
    Norm of gradient: 0.019569191810207776
    Iterations before filter: 30
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   0.9598|  18.2610|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  23.1350
    Norm of gradient: 0.011461364156613545
    Iterations before filter: 40
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   0.5382|  14.2669|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  25.8571
    Norm of gradient: 0.005058645720581414
    Iterations before filter: 50
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   0.5519|  11.2071|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  27.7493
    Norm of gradient: 0.00913384867304721
    Iterations before filter: 60
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   0.5063|   9.0763|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  29.5937
    Norm of gradient: 0.004085182612931341
    Iterations before filter: 70
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   0.3574|   7.4804|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  30.4240
    Norm of gradient: 0.010016895672519415
    Iterations before filter: 80
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   0.2502|   7.4561|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  30.4265
    Norm of gradient: 0.0009764931617796537
    Iterations before filter: 90
    Applying low pass filter to wind field...
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |6575.6187|   6.1340|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|   8.7456
    Norm of gradient: 0.7462737969222745
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |5966.0414|   6.4019|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|   8.6014
    Norm of gradient: 0.6983294931210702
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    | 396.1297|  84.0446|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  11.6213
    Norm of gradient: 0.20222520819243076
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    | 210.3138|  87.8755|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  11.6785
    Norm of gradient: 0.12040101890762513
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |  46.4198|  87.9840|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  11.2879
    Norm of gradient: 0.08702823439023104
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |  47.4066|  74.8731|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  11.3296
    Norm of gradient: 0.1633964170150081
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |  13.2350|  73.1213|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  11.3530
    Norm of gradient: 0.03300072469301934
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   8.7393|  68.0998|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  11.3983
    Norm of gradient: 0.01919300030010946
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   6.1531|  54.8219|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  11.6947
    Norm of gradient: 0.028416505852986308
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   5.4481|  44.6685|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  12.4289
    Norm of gradient: 0.019417630225052232
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |  16.3036|  33.2967|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  14.2565
    Norm of gradient: 0.07539028471955857
    | Jvel    | Jmass   | Jsmooth |   Jbg   | Jvort   | Jmodel  | Jpoint  | Max w  
    |   6.9822|  36.7060|   0.0000|   0.0000|   0.0000|   0.0000|   0.0000|  13.3551
    Norm of gradient: 0.04778436338484203
    Iterations after filter: 1
    Iterations after filter: 2
    Done! Time = 116.8
    /home/rjackson/anaconda3/envs/pyart-2020/lib/python3.7/site-packages/pydda/vis/streamline_plot.py:396: RuntimeWarning: invalid value encountered in less
      w_filled = np.ma.masked_where(w[level, :, :] < np.min(w_vel_contours),
    /home/rjackson/anaconda3/envs/pyart-2020/lib/python3.7/site-packages/cartopy/mpl/geoaxes.py:1400: UserWarning: linewidths is ignored by contourf
      result = matplotlib.axes.Axes.contourf(self, *args, **kwargs)
    /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/cartopy/mpl/geoaxes.py:1366: UserWarning: The following kwargs were not used by contour: 'color'
      result = matplotlib.axes.Axes.contour(self, *args, **kwargs)
    /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/cartopy/mpl/geoaxes.py:1366: UserWarning: The following kwargs were not used by contour: 'color'
      result = matplotlib.axes.Axes.contour(self, *args, **kwargs)






|


.. code-block:: default


    import pydda
    import pyart
    import cartopy.crs as ccrs
    import matplotlib.pyplot as plt


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

    # Load our radar data
    sounding = pyart.io.read_arm_sonde(
        pydda.tests.SOUNDING_PATH)
    u_init, v_init, w_init = pydda.initialization.make_constant_wind_field(
        berr_grid, (0.0, 0.0, 0.0))

    # Let's make a plot on a map
    fig = plt.figure(figsize=(7, 7))
    ax = plt.axes(projection=ccrs.PlateCarree())

    pydda.vis.plot_horiz_xsection_streamlines_map(
        [cpol_grid, berr_grid], ax=ax, bg_grid_no=-1, level=7, w_vel_contours=[3, 5, 8])
    plt.show()

    # Let's see what happens when we use a zero initialization
    new_grids = pydda.retrieval.get_dd_wind_field([cpol_grid, berr_grid],
                                        u_init, v_init, w_init,
                                        Co=1.0, Cm=1500.0, frz=5000.0,
                                        mask_outside_opt=False)

    fig = plt.figure(figsize=(7, 7))
    ax = plt.axes(projection=ccrs.PlateCarree())

    pydda.vis.plot_horiz_xsection_streamlines_map(
        new_grids, ax=ax, bg_grid_no=-1, level=7, w_vel_contours=[3, 5, 8])
    plt.show()

    # Or, let's make the radar data more important!
    new_grids = pydda.retrieval.get_dd_wind_field([cpol_grid, berr_grid],
                                        u_init, v_init, w_init,
                                        Co=1.0, Cm=1500.0, frz=5000.0,
                                        mask_outside_opt=False)
    fig = plt.figure(figsize=(7, 7))
    ax = plt.axes(projection=ccrs.PlateCarree())

    pydda.vis.plot_horiz_xsection_streamlines_map(
        new_grids, ax=ax, bg_grid_no=-1, level=7, w_vel_contours=[3, 5, 8])
    plt.show()


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

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


.. _sphx_glr_download_source_auto_examples_plot_fun_with_constraints.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_fun_with_constraints.py <plot_fun_with_constraints.py>`



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

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


.. only:: html

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

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