# Over aggressive filtration **URL:** \ **Category:** Corrections/Data Quality\ **Tags:** py-art\ **Created:** [February 8, 2023, 9:14pm UTC](https://openradar.discourse.group/t/over-aggressive-filtration/160 "2023-02-08T21:14:54Z")\ **Posts on this page:** 5\ **Page:** 1
**Author:** ![syedhamidali](https://yyz2.discourse-cdn.com/free1/user_avatar/openradar.discourse.group/syedhamidali/32/377_2.png) [@syedhamidali](https://openradar.discourse.group/u/syedhamidali)\ **Post date:** [February 8, 2023, 9:14pm UTC](https://openradar.discourse.group/t/over-aggressive-filtration/160/1 "2023-02-08T21:14:54Z")
Hi All, I want your help with noise removal and dealiasing. I am conducting quality checks and filtration on my data. However, the process also removes good data. I attempted to adjust the values based on the characteristics of the field, but it has yet to be effective. The current method of keeping the VT above a certain threshold still leaves some noise present, while if the threshold is set too low, the filtration becomes too aggressive. I would greatly appreciate any advice or guidance you can give me. Please look at the plot; if you think “it is a tiny patch”, it is possibly a tornado, and I want to perform some dual doppler analysis, and neither can’t afford noise nor lose even that much of the data. Thank you! ![download](https://global.discourse-cdn.com/free1/uploads/openradar/original/1X/ce50c1787d3c93e21cb0a9d0bd8aa1fe610a3f8c.png) ```auto import pyart import numpy as np from pyart.retrieve import get_freq_band from matplotlib import pyplot as plt import glob, os def align_field(field): values, counts = np.unique(field['data'], return_counts=True) c_value = values[np.argmax(counts)] field['data'] = np.array([c_value]) return field def align_radar_coords(radar): radar.longitude = align_field(radar.longitude) radar.latitude = align_field(radar.latitude) radar.altitude = align_field(radar.altitude) radar.altitude_agl = align_field(radar.altitude_agl) return radar def filter_radar(radar, vel_field="VEL_F", refl_field="DBZHCC_F", ncp_field="NCP", rhv_field="RHOHV", phi_field="PHIDP"): '''Remove noise based on velocity texture and mask all the fields''' ##Drop some fields fields_to_drop = ["DBMHC", "DBMVC", "VS", "VS_F", "VL", "VL_F", "DBZHCC"] for field in fields_to_drop: if field in radar.fields: del radar.fields[field] texture = pyart.retrieve.calculate_velocity_texture(radar, vel_field=vel_field, wind_size=5, check_nyq_uniform=False) radar.add_field('VT',texture,replace_existing=True) ##create gatefilter gf = pyart.filters.GateFilter(radar) BAND = get_freq_band(radar.instrument_parameters['frequency']['data']) if BAND == "X": vt = 5 if BAND == "C": vt = 25 gf.exclude_above('VT', vt) gf.exclude_outside(refl_field, -30, 100) gf.exclude_below("SNRHC", -3) gf_despeckeld = pyart.correct.despeckle_field(radar, refl_field, gatefilter=gf) corr_ZH = radar.fields[refl_field].copy() ZH_array = np.ma.masked_where(gf_despeckeld.gate_included == False, radar.fields[refl_field]['data']) corr_ZH['data'] = ZH_array radar.add_field('DBZH',corr_ZH,replace_existing=True) ##get the mask mask = np.ma.getmask(radar.fields['DBZH']['data']) radar.scan_type = b'ppi' ##iterate through remaining fields skip_fields = ["DBZHC", "VEL"] for field in radar.fields.keys(): if any(field == skip_field for skip_field in skip_fields): continue ##mask the field radar.fields[field]['data'] = np.ma.masked_where(mask, radar.fields[field]['data']) return radar def dealiase(radar, vel_name = "VEL_F"): try: gatefilter = pyart.correct.GateFilter(radar) corr_vel = pyart.correct.dealias_region_based( radar, vel_field=vel_name, keep_original=False, gatefilter = gatefilter) radar.add_field(vel_name, corr_vel, replace_existing=True) except: None def qc(file=None, savedir = "COW1/QC/"): radar = pyart.io.read(file) radar = align_radar_coords(radar) print("") print(radar.fields.keys()) print("Coords: ", radar.longitude['data'], radar.latitude['data'], radar.altitude['data'],) radar = filter_radar(radar) dealiase(radar) ##pyart.io.write(os.path.join(basedir, savedir, file)) return radar def plot_radar(radar): # Let us view the vleocity with the filter applied. fig = plt.figure(figsize=[10, 8]) display = pyart.graph.RadarMapDisplay(radar) ax = plt.subplot(221) display.plot_ppi("DBZHC", vmin=-10., vmax=70., cmap='pyart_HomeyerRainbow') ax = plt.subplot(222) display.plot_ppi("DBZH", vmin=-10., vmax=70., cmap='pyart_HomeyerRainbow', title='filtered') ax = plt.subplot(223) display.plot_ppi("VEL", sweep=0, cmap="pyart_NWSVel", vmin = -30, vmax = 30,) ax = plt.subplot(224) display.plot_ppi("VEL_F", sweep=0, cmap="pyart_NWSVel", vmin = -30, vmax = 30, title = "filtered") plt.show() basedir = "/depot/dawson29/data/Projects/PERiLS/obsdata/2022/Illinois_Mobile_Radar/IOP1/" cow = os.path.join(basedir, "COW1/cmerged/20220322/") files = glob.glob(os.path.join(cow, "*")) files.sort() len(files) radar = qc(files[112]) plot_radar(radar) ``` In case you need a data file, here is a [link](https://drive.google.com/file/d/1X2qcyosVPA6eubH7GAfwiN6GIfosXKFt/view?usp=share_link). ---
**Author:** ![scollis](https://yyz2.discourse-cdn.com/free1/user_avatar/openradar.discourse.group/scollis/32/21_2.png) [@scollis](https://openradar.discourse.group/u/scollis)\ **Post date:** [February 8, 2023, 9:49pm UTC](https://openradar.discourse.group/t/over-aggressive-filtration/160/2 "2023-02-08T21:49:46Z")
Yes we have had challenges with this as well. The TVS is creating texture… Have you tried adding an include\_inside to force Z \> some threshold back in unless the TVS is right on the edge of a BWER or suchlike it should have a decent Z signature as well… Or you can use a fuzzy logic approach… This code really needs some refresh but this is how we do it on ARM data [cmac2.0/cmac\_processing.py at main · EVS-ATMOS/cmac2.0 · GitHub](https://github.com/EVS-ATMOS/cmac2.0/blob/main/cmac/cmac_processing.py#L108) ---
**Author:** ![rcjackson](https://yyz2.discourse-cdn.com/free1/user_avatar/openradar.discourse.group/rcjackson/32/14_2.png) [@rcjackson](https://openradar.discourse.group/u/rcjackson)\ **Post date:** [February 9, 2023, 4:06pm UTC](https://openradar.discourse.group/t/over-aggressive-filtration/160/3 "2023-02-09T16:06:00Z")
In addition to what Scott suggests after putting the reflectivity greater than ~10 dBZ back in the field I would apply despeckling to get rid of the isolated speckles that can create errors in dealiasing. This should not affect your tornado since it’s embedded in an MCS. ---
**Author:** ![syedhamidali](https://yyz2.discourse-cdn.com/free1/user_avatar/openradar.discourse.group/syedhamidali/32/377_2.png) [@syedhamidali](https://openradar.discourse.group/u/syedhamidali)\ **Post date:** [February 18, 2023, 5:58am UTC](https://openradar.discourse.group/t/over-aggressive-filtration/160/4 "2023-02-18T05:58:27Z")
Thank you, @scollis and @rcjackson. Your recommendations worked; However, I’d like to mention that I had to fine-tune the wind size. ---
**Author:** ![scollis](https://yyz2.discourse-cdn.com/free1/user_avatar/openradar.discourse.group/scollis/32/21_2.png) [@scollis](https://openradar.discourse.group/u/scollis)\ **Post date:** [February 18, 2023, 2:42pm UTC](https://openradar.discourse.group/t/over-aggressive-filtration/160/5 "2023-02-18T14:42:41Z")
Great to know! Would be great if you (eventually) share some code… What would be amazing is a cookbook on this cc @mgrover1