.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "examples/mapping/plot_compare_two_radars_gatemapper.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code .. rst-class:: sphx-glr-example-title .. _sphx_glr_examples_mapping_plot_compare_two_radars_gatemapper.py: ====================================== Compare Two Radars Using Gatemapper ====================================== Map the reflectivity field of a single radar in Antenna coordinates to another radar in Antenna coordinates and compare the fields. .. GENERATED FROM PYTHON SOURCE LINES 10-26 .. code-block:: Python print(__doc__) # Author: Max Grover (mgrover@anl.gov) and Bobby Jackson (rjackson@anl.gov) # License: BSD 3 clause import warnings import cartopy.crs as ccrs import matplotlib.pyplot as plt import numpy as np import pyart from pyart.testing import get_test_data warnings.filterwarnings("ignore") .. GENERATED FROM PYTHON SOURCE LINES 27-30 **Read in the Data** For this example, we use two XSAPR radars from our test data. .. GENERATED FROM PYTHON SOURCE LINES 30-37 .. code-block:: Python # read in the data from both XSAPR radars xsapr_sw_file = get_test_data("swx_20120520_0641.nc") xsapr_se_file = get_test_data("sex_20120520_0641.nc") radar_sw = pyart.io.read_cfradial(xsapr_sw_file) radar_se = pyart.io.read_cfradial(xsapr_se_file) .. GENERATED FROM PYTHON SOURCE LINES 38-46 **Filter and Configure the GateMapper** We are interested in mapping the southwestern radar to the southeastern radar. Before running our gatemapper, we add a filter for only positive reflectivity values. We also need to set a distance (meters) and time (seconds) between the source and destination gate allowed for an adequate match), using the distance_tolerance/time_tolerance variables. .. GENERATED FROM PYTHON SOURCE LINES 46-58 .. code-block:: Python gatefilter = pyart.filters.GateFilter(radar_sw) gatefilter.exclude_below("reflectivity_horizontal", 20) gmapper = pyart.map.GateMapper( radar_sw, radar_se, distance_tolerance=500.0, time_tolerance=60, gatefilter_src=gatefilter, ) radar_sw_mapped_to_radar_se = gmapper.mapped_radar(["reflectivity_horizontal"]) .. GENERATED FROM PYTHON SOURCE LINES 59-63 **Plot the Original Data** Let's take a look at our original fields - notice the difference in reflectivity values! .. GENERATED FROM PYTHON SOURCE LINES 63-99 .. code-block:: Python fig = plt.figure(figsize=(16, 6)) ax = plt.subplot(121, projection=ccrs.PlateCarree()) # Plot the southwestern radar disp1 = pyart.graph.RadarMapDisplay(radar_sw) disp1.plot_ppi_map( "reflectivity_horizontal", sweep=1, ax=ax, vmin=-20, vmax=70, min_lat=36, max_lat=37, min_lon=-98, max_lon=-97, lat_lines=np.arange(36, 37.25, 0.25), lon_lines=np.arange(-98, -96.75, 0.25), ) ax2 = plt.subplot(122, projection=ccrs.PlateCarree()) disp2 = pyart.graph.RadarMapDisplay(radar_se) disp2.plot_ppi_map( "reflectivity_horizontal", sweep=1, ax=ax2, vmin=-20, vmax=70, min_lat=36, max_lat=37, min_lon=-98, max_lon=-97, lat_lines=np.arange(36, 37.25, 0.25), lon_lines=np.arange(-98, -96.75, 0.25), ) .. image-sg:: /examples/mapping/images/sphx_glr_plot_compare_two_radars_gatemapper_001.png :alt: xsapr-sg 1.5 Deg. 2011-05-20T06:42:26Z Equivalent reflectivity factor, xsapr-sg 1.5 Deg. 2011-05-20T06:41:37Z Equivalent reflectivity factor :srcset: /examples/mapping/images/sphx_glr_plot_compare_two_radars_gatemapper_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 100-102 Now, we can compare our original field from the southwestern radar, to the new remapped field - there are similarities... .. GENERATED FROM PYTHON SOURCE LINES 102-139 .. code-block:: Python fig = plt.figure(figsize=(16, 6)) ax = plt.subplot(121, projection=ccrs.PlateCarree()) # Plot the southeastern radar disp1 = pyart.graph.RadarMapDisplay(radar_se) disp1.plot_ppi_map( "reflectivity_horizontal", sweep=1, ax=ax, vmin=-20, vmax=70, min_lat=36, max_lat=37, min_lon=-98, max_lon=-97, lat_lines=np.arange(36, 37.25, 0.25), lon_lines=np.arange(-98, -96.75, 0.25), ) # Plot the southwestern radar mapped to the southeastern radar ax2 = plt.subplot(122, projection=ccrs.PlateCarree()) disp2 = pyart.graph.RadarMapDisplay(radar_sw_mapped_to_radar_se) disp2.plot_ppi_map( "reflectivity_horizontal", sweep=1, ax=ax2, vmin=-20, vmax=70, min_lat=36, max_lat=37, min_lon=-98, max_lon=-97, lat_lines=np.arange(36, 37.25, 0.25), lon_lines=np.arange(-98, -96.75, 0.25), ) .. image-sg:: /examples/mapping/images/sphx_glr_plot_compare_two_radars_gatemapper_002.png :alt: xsapr-sg 1.5 Deg. 2011-05-20T06:41:37Z Equivalent reflectivity factor, xsapr-sg 1.5 Deg. 2011-05-20T06:41:37Z Equivalent reflectivity factor :srcset: /examples/mapping/images/sphx_glr_plot_compare_two_radars_gatemapper_002.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 140-147 **Calculate and Plot the Difference** It can be difficult to "eyeball" the difference between these two fields. Fortunately, now that our radars match coordinates, we can plot a difference. Keep in mind there is a time difference of ~ 1 minute between these plots, leading to small difference due to the precipitation moving through the domain over the course of that minute. .. GENERATED FROM PYTHON SOURCE LINES 147-187 .. code-block:: Python # Extract the numpy arrays for our reflectivity fields reflectivity_se_radar = radar_se.fields["reflectivity_horizontal"]["data"] reflectivity_sw_radar = radar_sw_mapped_to_radar_se.fields["reflectivity_horizontal"][ "data" ] # Calculate the difference between the southeastern and southwestern radar reflectivity_difference = reflectivity_se_radar - reflectivity_sw_radar # Add a field like this to the radar_se radar object radar_se.add_field_like( "reflectivity_horizontal", field_name="reflectivity_bias", data=reflectivity_difference, ) # Setup our figure fig = plt.figure(figsize=(8, 6)) ax = plt.subplot(111, projection=ccrs.PlateCarree()) # Plot the difference field disp1 = pyart.graph.RadarMapDisplay(radar_se) disp1.plot_ppi_map( "reflectivity_bias", cmap="pyart_balance", title="Reflectivity Difference \n XSAPR Southwest - XSPAR Southeast", sweep=1, ax=ax, vmin=-30, vmax=30, min_lat=36, max_lat=37, min_lon=-98, max_lon=-97, lat_lines=np.arange(36, 37.25, 0.25), lon_lines=np.arange(-98, -96.75, 0.25), ) .. image-sg:: /examples/mapping/images/sphx_glr_plot_compare_two_radars_gatemapper_003.png :alt: Reflectivity Difference XSAPR Southwest - XSPAR Southeast :srcset: /examples/mapping/images/sphx_glr_plot_compare_two_radars_gatemapper_003.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 188-193 **Plot a Histogram for Comparison** Another way of plotting the comparison here is using a 2-dimensional histogram,which is more helpful in this case where our scans don't neccessarily match exactly in time. .. GENERATED FROM PYTHON SOURCE LINES 193-226 .. code-block:: Python # Include elevations above the lowest one incl_gates = np.argwhere(radar_sw_mapped_to_radar_se.elevation["data"] > 1.0) # Filter the reflectivity fields using the filter created above refl_se = reflectivity_se_radar[incl_gates, :] refl_sw = reflectivity_sw_radar[incl_gates, :] # Make sure not include masked values values_without_mask = np.logical_and(~refl_se.mask, ~refl_sw.mask) refl_se = refl_se[values_without_mask] refl_sw = refl_sw[values_without_mask] # Set the bins for our histogram bins = np.arange(-10, 60, 1) # Create the 2D histogram using the flattened numpy arrays hist = np.histogram2d(refl_se.flatten(), refl_sw.flatten(), bins=bins)[0] hist = np.ma.masked_where(hist == 0, hist) # Setup our figure fig = plt.figure(figsize=(8, 6)) # Create a 1-1 comparison x, y = np.meshgrid((bins[:-1] + bins[1:]) / 2.0, (bins[:-1] + bins[1:]) / 2.0) c = plt.pcolormesh(x, y, np.log10(hist.T), cmap="pyart_HomeyerRainbow") # Add a colorbar and labels plt.colorbar(c, label="$log_{10}$ counts") plt.xlabel("XSAPR Southeast $Z_{H}$ [dBZ]") plt.ylabel("XSAPR Southwest $Z_{H}$ [dBZ]") plt.show() .. image-sg:: /examples/mapping/images/sphx_glr_plot_compare_two_radars_gatemapper_004.png :alt: plot compare two radars gatemapper :srcset: /examples/mapping/images/sphx_glr_plot_compare_two_radars_gatemapper_004.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 37.317 seconds) .. _sphx_glr_download_examples_mapping_plot_compare_two_radars_gatemapper.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_compare_two_radars_gatemapper.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_compare_two_radars_gatemapper.py ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_