search

ARM-DOE / pyart

The Python-ARM Radar Toolkit. A data model driven interactive toolkit for working with weather radar data.
https://arm-doe.github.io/pyart/
Other
516 stars 268 forks source link

plotting rhi's created with over the top scanning mode #1165

Closed srbrodzik closed 2 years ago

srbrodzik commented 2 years ago

I'm running python3 on a debian 11 system to create plots from SPOL cfradial RHI sector files at various azimuths. They all plot nicely except for the one scan sequence of 2 rhi's at 101 and 281 deg. They are created by scanning the radar up at one azimuth and then over the top into the other azimuth. When I plot the data, only one of the azimuth angles plots correctly. The other one produces a bunch of blank plots. I'm thinking it must be the cfradial files themselves. Would love any advice out there. Here's the plotting portion of my code. I don't see a way to submit a sample cfradial file but that would probably be helpful.

` # create display display = pyart.graph.RadarDisplay(radar) fig = plt.figure(figsize=(12, 10))

fig.suptitle('SPOL '+str(angle)+' deg '+start_str, fontsize=20)

        fig.suptitle('SPOL '+start_str+' UTC RHI '+str(angle)+' deg', fontsize=20)

        ax = fig.add_subplot(321)
        display.plot('DBZ', sweep=index, ax=ax, title='',
                     colorbar_label='Reflectivity (dBZ)',
                     axislabels=('', 'Distance above radar (km)'))
        display.set_limits((0, 150), (0, 18), ax=ax)

        ax = fig.add_subplot(322)
        display.plot('VEL', sweep=index, ax=ax,
                     title='', colorbar_label='Radial Velocity (m/s)',
                     cmap='pyart_Carbone42',
                     axislabels=('', ''))
        display.set_limits((0, 150), (0, 18), ax=ax)

        ax = fig.add_subplot(323)
        display.plot('ZDR', sweep=index, ax=ax,vmin=-1,vmax=4,
                     cmap='nipy_spectral',
                     axislabels=('', 'Distance above radar (km)'),
                     title='', colorbar_label='ZDR (dB)')
        display.set_limits((0, 150), (0, 18), ax=ax)

        ax = fig.add_subplot(324)
        display.plot('RHOHV', sweep=index, ax=ax,vmin=0,vmax=1,
                     cmap='nipy_spectral',
                     title='', colorbar_label='RHOHV',
                     axislabels=('', ''))
        display.set_limits((0, 150), (0, 18), ax=ax)

        ax = fig.add_subplot(325)
        display.plot('KDP', sweep=index, ax=ax,vmin=-1,vmax=4,
                     cmap='nipy_spectral',
                     title='', colorbar_label='KDP (deg/km)',
                     axislabels=('Distance from radar (km)', 'Distance above radar (km)'))
        display.set_limits((0, 150), (0, 18), ax=ax)

        ax = fig.add_subplot(326)
        display.plot('PHIDP', sweep=index, ax=ax,vmin=-180,vmax=180,
                     cmap='nipy_spectral',
                     title='', colorbar_label='PHIDP (deg)',
                     axislabels=('Distance from radar (km)', ''))
        display.set_limits((0, 150), (0, 18), ax=ax)

        # save image
        #plt.show()
        plt.savefig(outdir+'/'+file_out)

`

mgrover1 commented 2 years ago

Can you share a sample data file? You can send it to mgrover@anl.gov!

mgrover1 commented 2 years ago

@srbrodzik - I had a chance to look at this today. Here are some thoughts!

I'm running python3 on a debian 11 system to create plots from SPOL cfradial RHI sector files at various azimuths. They all plot nicely except for the one scan sequence of 2 rhi's at 101 and 281 deg.

First, I noticed that we need to manually set the scan mode (assuming these are RHIs), although it sounds like there is a bit of a custom scan strategy here?

radar = pyart.io.read("cfrad.20220531_080000.424_to_20220531_080047.680_SPOL_PrecipRhi1_RHI.nc")
radar.scan_type = 'rhi'

I was able to reproduce the issue here, using the data sent over. Here is a plot at each of the azimuths, using this block of code:

cmap_dict = {"DBZ": "pyart_HomeyerRainbow",
             "VEL": "pyart_balance",
             "ZDR": "pyart_Carbone42",
             "RHOHV": "pyart_Carbone42",
             "KDP": "pyart_Carbone42",
             "PHIDP": "pyart_Carbone42",}
display = pyart.graph.RadarDisplay(radar)
for sweep in range(radar.nsweeps):
    for field in cmap_dict.keys():
        display.plot(field, sweep, cmap=cmap_dict[field])
        plt.ylim(0, 20)
        plt.show()
        plt.close()

Focusing on just the reflectivity field (DBZ), these are the two plots created

Screen Shot 2022-06-02 at 9 46 18 AM

It looks like it is having an issue with the second one (281 deg)...

They are created by scanning the radar up at one azimuth and then over the top into the other azimuth.

Could you expand on this? What do you mean by "over the top into the other azimuth"? Is this a typical RHI scan? A vertically pointing scan?

I'm thinking it must be the cfradial files themselves.

I suspect this is the case... maybe it would be helpful to setup a call to chat about this? I know that ya'll are in a very different timezone, so perhaps we could do earlier tomorrow (June 3) early in the morning my time (7:30 AM CDT), which would be 8:30 PM in Taiwan?

mgrover1 commented 2 years ago

Also, here is a plot of the second RHI zoomed out a bit - there appears to be something picked up by the radar, but I am not exactly sure what it is... Screen Shot 2022-06-02 at 10 05 30 AM

srbrodzik commented 2 years ago

Thanks Max.

Let me try your code tomorrow (it's Thursday night around 11pm here right now) and I'll let you know what I get. I'm not sure why your code shows data at 101 deg and mine doesn't. Would love to talk it over Friday night my time, Friday morning your time. We'd have to talk via Signal or WhatsApp since I currently have a Taiwan phone number. Or we could always use Zoom. I've set up Signal and WhatsApp to use my US number so you'd just have to dial that (206-660-0854) to reach me. Or I can set up a Zoom call and share the link with you.

FYI, this is kind of a weird RHI 'volume'. It goes up at 101 deg (I think that angle is first but it could be the other way around), stops at the top, spins the dish around and then comes down at 281 deg. We're doing this because there's another radar to the east and this is the line from one radar to the other. The second radar will do some similar scans so we can compare data.

Stacy

On Thu, Jun 2, 2022 at 11:08 PM Max Grover @.***> wrote:

Also, here is a plot of the second RHI zoomed out a bit - there appears to be something picked up by the radar, but I am not exactly sure what it is... [image: Screen Shot 2022-06-02 at 10 05 30 AM] https://user-images.githubusercontent.com/26660300/171660751-d5ff9cb1-1c0e-487d-8a26-d6f4d27a0d77.png

— Reply to this email directly, view it on GitHub https://github.com/ARM-DOE/pyart/issues/1165#issuecomment-1144976329, or unsubscribe https://github.com/notifications/unsubscribe-auth/ABHT7FU6IOCEFRMTT4AULXLVNDE5VANCNFSM5XQN7L7A . You are receiving this because you were mentioned.Message ID: @.***>

mgrover1 commented 2 years ago

Great! Sounds good - I will send you a calendar invite for Zoom to the email you used to send the data file. I looped through the different sweeps, so perhaps the sweep number was set to 1 when you tried it? Here is the full script if you want to copy/paste and give it a try. Looking forward to chatting.

import pyart
import matplotlib.pyplot as plt

radar = pyart.io.read("cfrad.20220531_080000.424_to_20220531_080047.680_SPOL_PrecipRhi1_RHI.nc")
radar.scan_type = 'rhi'

cmap_dict = {"DBZ": "pyart_HomeyerRainbow",
             "VEL": "pyart_balance",
             "ZDR": "pyart_Carbone42",
             "RHOHV": "pyart_Carbone42",
             "KDP": "pyart_Carbone42",
             "PHIDP": "pyart_Carbone42",}
display = pyart.graph.RadarDisplay(radar)
for sweep in range(radar.nsweeps):
    for field in cmap_dict.keys():
        display.plot(field, sweep, cmap=cmap_dict[field])
        plt.ylim(0, 20)
        plt.show()
        plt.close()
srbrodzik commented 2 years ago

Thanks for the help. I'll check the sweep index in my code when I get to the Ops Center. I hope it's not something that stupid:-)

On Thu, Jun 2, 2022 at 11:24 PM Max Grover @.***> wrote:

Great! Sounds good - I will send you a calendar invite for Zoom to the email you used to send the data file. I looped through the different sweeps, so perhaps the sweep number was set to 1 when you tried it? Here is the full script if you want to copy/paste and give it a try. Looking forward to chatting.

import pyartimport matplotlib.pyplot as plt radar = pyart.io.read("cfrad.20220531_080000.424_to_20220531_080047.680_SPOL_PrecipRhi1_RHI.nc")radar.scan_type = 'rhi' cmap_dict = {"DBZ": "pyart_HomeyerRainbow", "VEL": "pyart_balance", "ZDR": "pyart_Carbone42", "RHOHV": "pyart_Carbone42", "KDP": "pyart_Carbone42", "PHIDP": "pyart_Carbone42",}display = pyart.graph.RadarDisplay(radar)for sweep in range(radar.nsweeps): for field in cmap_dict.keys(): display.plot(field, sweep, cmap=cmap_dict[field]) plt.ylim(0, 20) plt.show() plt.close()

— Reply to this email directly, view it on GitHub https://github.com/ARM-DOE/pyart/issues/1165#issuecomment-1144994258, or unsubscribe https://github.com/notifications/unsubscribe-auth/ABHT7FQN7KB6KTINMNUKYNDVNDG3FANCNFSM5XQN7L7A . You are receiving this because you were mentioned.Message ID: @.***>

srbrodzik commented 2 years ago

Hi Max,

OK, your snippet of code worked for me so I added things from my code to it until I discovered the problem. I was putting in xlim = (0,150) but for this particular angle in this particular scanning sequence, I needed to set xlim = (-150,0).

Guess we don't need to chat about this tonight. Thanks for your willingness to help me figure this out. I really appreciate it.

I'm sure I'll be in touch down the road.

Stacy

On Thu, Jun 2, 2022 at 3:26 PM Stacy Brodzik @.***> wrote:

Thanks for the help. I'll check the sweep index in my code when I get to the Ops Center. I hope it's not something that stupid:-)

On Thu, Jun 2, 2022 at 11:24 PM Max Grover @.***> wrote:

Great! Sounds good - I will send you a calendar invite for Zoom to the email you used to send the data file. I looped through the different sweeps, so perhaps the sweep number was set to 1 when you tried it? Here is the full script if you want to copy/paste and give it a try. Looking forward to chatting.

import pyartimport matplotlib.pyplot as plt radar = pyart.io.read("cfrad.20220531_080000.424_to_20220531_080047.680_SPOL_PrecipRhi1_RHI.nc")radar.scan_type = 'rhi' cmap_dict = {"DBZ": "pyart_HomeyerRainbow", "VEL": "pyart_balance", "ZDR": "pyart_Carbone42", "RHOHV": "pyart_Carbone42", "KDP": "pyart_Carbone42", "PHIDP": "pyart_Carbone42",}display = pyart.graph.RadarDisplay(radar)for sweep in range(radar.nsweeps): for field in cmap_dict.keys(): display.plot(field, sweep, cmap=cmap_dict[field]) plt.ylim(0, 20) plt.show() plt.close()

— Reply to this email directly, view it on GitHub https://github.com/ARM-DOE/pyart/issues/1165#issuecomment-1144994258, or unsubscribe https://github.com/notifications/unsubscribe-auth/ABHT7FQN7KB6KTINMNUKYNDVNDG3FANCNFSM5XQN7L7A . You are receiving this because you were mentioned.Message ID: @.***>

srbrodzik commented 2 years ago

I found another data oddness. If you need me to submit this formally via github, I'll do that. If not, here's the story.

We're collecting a different kind of rhi volume that has 12 widely-spaced azimuths that range from 0 to 360. They all plot up nicely except the 270 deg azimuth scan which has some artifacts above the radar. I used your code snippet to do the plotting and added a couple things so I could tell one plot from the next.

I thought it might be the cfradial file but the NCAR CIDD software displays these 270 deg angle rhi's without any artifacts. See the image called 'spol.202206030411.rhi_270_cidd.png'.

Everything is attached. Let me know if you have any thoughts. You'll probably need to take a look through this so I don't expect to talk about it tonight. Just take your time and let me know what you find. If we need to chat at that point, we can set something up.

Thanks for your help. I really appreciate it especially during this field work.

Stacy

cfrad.20220603_040943.316_to_20220603_041149.93... https://drive.google.com/file/d/1Ij1gccjsSCk6m2ExkdKE5G72KROf0-n5/view?usp=drive_web

On Fri, Jun 3, 2022 at 6:45 AM Stacy Brodzik @.***> wrote:

Hi Max,

OK, your snippet of code worked for me so I added things from my code to it until I discovered the problem. I was putting in xlim = (0,150) but for this particular angle in this particular scanning sequence, I needed to set xlim = (-150,0).

Guess we don't need to chat about this tonight. Thanks for your willingness to help me figure this out. I really appreciate it.

I'm sure I'll be in touch down the road.

Stacy

On Thu, Jun 2, 2022 at 3:26 PM Stacy Brodzik @.***> wrote:

Thanks for the help. I'll check the sweep index in my code when I get to the Ops Center. I hope it's not something that stupid:-)

On Thu, Jun 2, 2022 at 11:24 PM Max Grover @.***> wrote:

Great! Sounds good - I will send you a calendar invite for Zoom to the email you used to send the data file. I looped through the different sweeps, so perhaps the sweep number was set to 1 when you tried it? Here is the full script if you want to copy/paste and give it a try. Looking forward to chatting.

import pyartimport matplotlib.pyplot as plt radar = pyart.io.read("cfrad.20220531_080000.424_to_20220531_080047.680_SPOL_PrecipRhi1_RHI.nc")radar.scan_type = 'rhi' cmap_dict = {"DBZ": "pyart_HomeyerRainbow", "VEL": "pyart_balance", "ZDR": "pyart_Carbone42", "RHOHV": "pyart_Carbone42", "KDP": "pyart_Carbone42", "PHIDP": "pyart_Carbone42",}display = pyart.graph.RadarDisplay(radar)for sweep in range(radar.nsweeps): for field in cmap_dict.keys(): display.plot(field, sweep, cmap=cmap_dict[field]) plt.ylim(0, 20) plt.show() plt.close()

— Reply to this email directly, view it on GitHub https://github.com/ARM-DOE/pyart/issues/1165#issuecomment-1144994258, or unsubscribe https://github.com/notifications/unsubscribe-auth/ABHT7FQN7KB6KTINMNUKYNDVNDG3FANCNFSM5XQN7L7A . You are receiving this because you were mentioned.Message ID: @.***>