Wolfram Data Repository

Immediate Computable Access to Curated Contributed Data

Nocturnal Jet Wind Profiles Measured with WiPPR

Source Notebook

Measurements of vertical wind profiles measured at night in the Arizona desert using a high dynamic range, Ka band radar

Details

The Wind Profiling Portable Radar (WiPPR) was developed to detect clear air scatter in the convective boundary layer. The radar operates at a carrier frequency of 33.4 GHz in the Ka band and uses a fixed 48 MHz sweep width. The radar can detect clear air scatter targets at altitudes up to 1500-2000 m. This altitude range represents the upper limits of the convective boundary layer (also called atmospheric boundary layer). These clear air scatterers are turbulent motions of the air associated with ever-present hydrodynamic-thermodynamic instabilities in the atmosphere. These phenomena have been observed all over the continental US. Their prevalence is most pronounced during time periods when solar illumination is high and the atmosphere is unstable. This turbulence is also present during stable atmospheric conditions but usually at lower altitudes. The radar can also detect trace precipitation, rainfall, and certain types of clouds. WiPPR radar data can also be processed to produce turbulent intensity as a function of altitude. This quantity is important in wind pressure loading problems. The radar has a dynamic range in excess of 70 dB and a very low system noise temperature. This allows for the detection of faint echoes in the presence of much larger echoes.
A fundamental design consideration for a radar is the choice of carrier frequency. The WiPPR signal needs to propagate up to the top of the convective boundary layer (CBL) and back without suffering large absorption losses. This is nominally a round trip distance of 3 km. Absorption of propagating radar waves in the atmosphere in the range 10-100 GHz is greatest at 22 GHz and 60 GHz due to absorption by water vapor molecules (22 GHz) and by oxygen molecules (60 GHz). For a fixed frequency, absorption is less at higher elevations than at lower elevations. Two-way absorption losses to the top of the CBL and back at 33.4 GHz, 60 GHz and 80 GHz are respectively 0.21 dB/km, 44 dB/km and 0.97 dB/km. It is obviously vital to avoid the oxygen absorption at 60 GHz. A carrier frequency of 33.4 GHz was chosen for WiPPR because of its low absorption and the slight advantage it presented over frequencies in the 80-100 GHz range during periods of rainfall. Falling hydrometers introduce additional absorption which is lower at 33.4 GHz than at frequencies in the 80-100 GHz band.
The fundamental data building block produced by WiPPR is the range velocity matrix (RVM). This data structure localizes echo producing objects in slant range-Doppler velocity space and enables the computation of wind velocity profiles. The size of each range velocity matrix (RVM) in this data set is 512 by 256. The extent in slant range is rmax=512d z where dz=3.125 meter. The horizontal span of the RVM is approximately (-c/4fctp,c/4fctp) where c is the speed of light, fc is the radar carrier frequency and tp is the radar pulse length. This works out to be approximately (-11.4,11.4) m/s. The relationship between d z and the bandwidth B of the radar pulse is d z=c/2B.
The data presented here were recorded in early morning in the Arizona desert in the summer of 2017 at a location just south of Quartzsite, AZ. During the timespan of these data, nocturnal jet winds had formed centered at about 500 m in altitude. Nocturnal jets generally form above the land at night under clear sky conditions and can produce wind speeds in excess of the geostrophic winds at higher altitudes. Nocturnal jets can produce very large wind shears in the lower atmospheric boundary layer.

Examples

Basic Examples (3) 

View the system as used in the data collection.:
In[1]:=
Thumbnail[ ResourceData["Nocturnal Jet Wind Profiles Measured with WiPPR", "systemPicture"], Large]
Out[1]=
Look at the first 5 data entries. The right hand number is a pointer to the beam direction associated with the RVM:
In[2]:=
rd = ResourceData[\!\(\* TagBox["\"\<Nocturnal Jet Wind Profiles Measured with WiPPR\>\"", #& , BoxID -> "ResourceTag-Nocturnal Jet Wind Profiles Measured with WiPPR-Input", AutoDelete->True]\)]; Take[rd, {1, 5}]
Out[3]=
The system operating geometry is best understood with the aid of a figure. Beams 1-4 respectively correspond to the directions north, east, south and west. The radar operates in a near vertical configuration and revolves clockwise in a 'stop and go' fashion sequentially. The radar is not moving during time periods when data is being recorded:
In[4]:=
ResourceData["Nocturnal Jet Wind Profiles Measured with WiPPR", "beamPicture"]
Out[4]=

Visualizations (5) 

Parameter values useful in the visualizations and analysis that follows are defined here:
In[5]:=
rd = ResourceData["Nocturnal Jet Wind Profiles Measured with WiPPR"]; a = ResourceData["Nocturnal Jet Wind Profiles Measured with WiPPR", "LogEchoMin"]; b = ResourceData["Nocturnal Jet Wind Profiles Measured with WiPPR", "LogEchoMax"]; VMin = ResourceData["Nocturnal Jet Wind Profiles Measured with WiPPR", "VDopplerMin"]; VMax = ResourceData["Nocturnal Jet Wind Profiles Measured with WiPPR", "VDopplerMax"]; dVDoppler = ResourceData["Nocturnal Jet Wind Profiles Measured with WiPPR", "dVDoppler"]; dz = ResourceData["Nocturnal Jet Wind Profiles Measured with WiPPR", "dz"];
For visualization purposes the range velocity matrices(RVMs) in the data set were rescaled prior to creation of resource data. The original RVM was recorded on a power scale of positive values. For display purpose the data were rescaled and displayed in a human viewable data format. This rescaling is reversed here. The resulting RVM is displayed on a decibel scale over a range of 40 dB. The roll-off color effect in the figure is caused by the shape of intermediate stage low pass filter in the radar electronics:
In[6]:=
nf = 30; image = First@rd[[nf]]; rvm = Exp[ImageData[ImageReflect[ImageApply[(b*# + a) &, image]]]]; rvmdb = 10 Log[10, rvm]; dBplotmin = 39; dBplotmax = 79; ReliefPlot[rvmdb, Sequence[ LightingAngle -> None, PlotRange -> {dBplotmin, dBplotmax}, ColorFunction -> ColorData[{"DarkBands", "Reverse"}], ClippingStyle -> {LightGray, Black}, FrameTicks -> Automatic, DataRange -> {{VMin, VMax}, {0, dz Length[rvmdb]}}, AspectRatio -> 1.2, FrameLabel -> {"Doppler velocity (m/s)", "altitude (m)"}, ImageSize -> 400, PlotLegends -> BarLegend[Automatic, LegendLabel -> "dB"]]]
Out[10]=
For the purpose of Doppler velocity estimation the range velocity matrix (RVM) data are converted to signal-to-noise (SNR) and binarized at a prescribed SNR threshold of 2 dB:
In[11]:=
snr = Map[#/Median[#] &, rvm]; SNRThresholdDB = 2; occupancy = ImageData[Binarize[Image[snr], 10^(SNRThresholdDB/10)]]; snrdb = 10 Log[10, snr]; snrdbmax = Max[snrdb]; gsnr = ReliefPlot[snr, Sequence[ LightingAngle -> None, PlotRange -> {snrdbmax - 40, snrdbmax}, ColorFunction -> ColorData[{"DarkBands", "Reverse"}], ClippingStyle -> {Gray, Black}, FrameTicks -> Automatic, DataRange -> {{VMin, VMax}, {0, dz Length[rvmdb]}}, AspectRatio -> 1.2, FrameLabel -> {"Doppler velocity (m/s)", "altitude (m)"}, ImageSize -> 400, PlotLabel -> "snr(dB)"]]; gocc = ReliefPlot[occupancy, Sequence[ LightingAngle -> None, PlotRange -> {0.25, 0.5}, ClippingStyle -> {White, Black}, FrameTicks -> Automatic, DataRange -> {{VMin, VMax}, {0, dz Length[rvmdb]}}, AspectRatio -> 1.2, FrameLabel -> {"Doppler velocity (m/s)", "altitude (m)"}, ImageSize -> 400, PlotLabel -> "occupancy"]]; GraphicsRow[{gsnr, gocc}, ImageSize -> 600]
Out[17]=
The information in the SNR and occupancy matrices can be used to estimate the Doppler velocity at a particular slant range. The product of SNR and occupancy is thresholded SNR. If we normalize a slice of thresholded SNR we have probability density function for SNR at that range. This can be used to estimate the Doppler velocity of the echo. In the following we find the Doppler velocity, slant range and SNR (power) of a clear air scatter echo:
In[18]:=
nr = 152; slice = snr[[nr]]*occupancy[[nr]]; If[Max[occupancy[[nr]]] > 0, pdf = slice/Total[slice]; Vecho = Table[ivv, {ivv, VMin, VMax, dVDoppler}]; Vest = pdf . Vecho; position = Ceiling[(Vest - VMin)/dVDoppler]; If[position > 256, position = 256]; {Vest, nr*dz, snr[[nr, position]]}, "No echo greater than threshold"]
Out[20]=
The SNR at this range as a function of Doppler velocity looks like this:
In[21]:=
ListPlot[slice, Sequence[ Joined -> True, Axes -> False, Frame -> True, DataRange -> {VMin, VMax}, FrameLabel -> {"Doppler velocity (m/s)", "SNR (power)"}]]
Out[21]=

Analysis (12) 

Wind velocity inversion is the process through which the wind velocity profile (vx(z),vy(z),vz(z)) at altitude z is estimated from scalar Doppler velocity measurements made on WiPPR radar beams pointed in the cardinal directions (N, E, S, W). If all scattering objects were moving horizontally then sophisticated procedures would not be required to determine wind velocity from Doppler velocity, but this is not the case. Falling hydrometers and updrafts are common in the convective boundary layer. The first step in the inversion process is to form a data catalog. The following code extracts all echo contacts with more than a specified amount of signal-to-noise(SNR) ratio and organizes the contacts in a catalog(table) with the structure: (altitude, file number, Doppler velocity, cosX, cosY, cosZ, beam number). The radar beam pointing direction of the echo in a catalog entry is (cosX, cosY, cosZ).

The following code constructs the WiPPR Doppler contact catalog:
In[22]:=
WiPPRDataCatalog[rd_] := Module[{a, b, SNRThresholdDB = 2.0, dVDoppler, VDopplerlocal, image, thisbeam, \[Phi]Correct, \[Theta], \[Theta]Vertical, \[Eta], dzLcoal, cosX, cosY, cosZ, occupancy, rvm, snr, snrFiltered, goodRanges, zLocal, slice, pdf, Vpdf, catalog, catalogHeadings, VDopMin, VDopMax}, a = ResourceData["Nocturnal Jet Wind Profiles Measured with WiPPR", "LogEchoMin"]; b = ResourceData["Nocturnal Jet Wind Profiles Measured with WiPPR", "LogEchoMax"]; VDopMin = -11.3094; VDopMax = 11.3985; dVDoppler = (VDopMax - VDopMin)/255; VDopplerlocal = Table[V, {V, VDopMin, VDopMax, dVDoppler}]; catalog = Reap[Do[image = First@rd[[nf]]; thisbeam = Last@rd[[nf]]; \[Phi]Correct = (thisbeam - 1) \[Pi]/ 2.0; \[Theta] = 80.0*\[Pi]/180.0; \[Theta]Vertical = \[Pi]/ 2 - \[Theta]; \[Eta] = {Cos[\[Theta]] Sin[\[Phi]Correct], Cos[\[Theta]] Cos[\[Phi]Correct], Sin[\[Theta]]}; {cosX, cosY, cosZ} = Chop[\[Eta]]; dzLcoal = 3.125*Cos[\[Theta]Vertical]; rvm = Exp[ ImageData[ImageReflect[ImageApply[(b*# + a) &, image]]]]; snr = Map[#/Median[#] &, rvm]; occupancy = ImageData[Binarize[Image[snr], 10^(SNRThresholdDB/10)]]; snrFiltered = occupancy*snr; goodRanges = Map[First, Position[Map[Max, occupancy], 1]]; Do[zLocal = nRange*dzLcoal; slice = snrFiltered[[nRange]]; pdf = slice/Total[slice]; Vpdf = pdf . VDopplerlocal; Sow[Join[{zLocal, nf, Vpdf, cosX, cosY, cosZ, thisbeam}]], {nRange, goodRanges}], {nf, 1, Length[rd]}]]; catalog[[2]][[1]]]
The catalog is calculated here:
In[23]:=
rd = ResourceData["Nocturnal Jet Wind Profiles Measured with WiPPR"]; catalog = WiPPRDataCatalog[rd]; catalogHeadings = {"z(m)", "nfile", "Vpdf(m/s)", "cosX", "cosY", "cosZ", "beam"};
A small portion of the contact catalog is:
In[24]:=
TableForm[Take[catalog, {1, 10}], TableHeadings -> {None, catalogHeadings}]
Out[24]=

Beams 1-4 are respectively pointed in the directions 1) north (0 deg), 2) east (90 deg), 3) south (180 deg) and 4) west (270 deg). The beams are nearly vertical with each beam elevated 80 deg from the horizontal, 10 deg off vertical. Provided that there is no vertical wind motion beam 3 directly measures the component of wind velocity component vy within the scale factor 1/cosθ where θ=80 deg. Beam 1 measures -vy within the scale factor. Beam 4 directly measures the vxcomponent of wind velocity within this scale factor. Beam 2 measures -vx within the scale factor.

The following block of code groups the catalog contents by beam and makes a quad plot of the beam Doppler data:
In[25]:=
ViewBeamDopplerData[catalog_, zmaxplot_] := Module[{cat, gD}, gD = Table[cat = Transpose@Select[catalog, #[[7]] == nbeam &]; ListPlot[Transpose[{cat[[3]], cat[[1]]}], Sequence[ PlotRange -> {{-11.3, 11.3}, {0, zmaxplot}}, Axes -> False, Frame -> True, AspectRatio -> 1., FrameLabel -> {"Doppler Velocity (m/s)", "Altitude (m)"}, PlotLabel -> Style[ ToString[nbeam], FontSize -> 9]], ImageSize -> 200], {nbeam, 1, 4}]; GraphicsGrid[{{gD[[1]], gD[[2]]}, {gD[[3]], gD[[4]]}}]]
When executed the result is:
In[26]:=
ViewBeamDopplerData[catalog, 1600]
Out[26]=

Note the 180 deg left-right symmetry in the data between beam 1 and 3 and also between beams 2 and 4. This provides a strong indication that the WiPPR system was correctly functioning.

The following code computes the data necessary for a visual estimation of (vx,vy) that exploits the symmetry in the WiPPR beam geometry. The possible effect of falling hydrometers or a vertical component of the wind is ignored in this computation but is considered later:
In[27]:=
ApparentVelocityBeamGrouped[catalog_] := Module[{\[Theta]Vertical, catT, catx4, catx2, caty1, caty3}, \[Theta]Vertical = ArcCos@First[catalog][[6]]; catT = Transpose@Select[catalog, #[[7]] == 4 &]; catx4 = Transpose[{catT[[3]]/Sin[\[Theta]Vertical], catT[[1]]}]; catT = Transpose@Select[catalog, #[[7]] == 2 &]; catx2 = Transpose[{-catT[[3]]/Sin[\[Theta]Vertical], catT[[1]]}]; catT = Transpose@Select[catalog, #[[7]] == 1 &]; caty1 = Transpose[{-catT[[3]]/Sin[\[Theta]Vertical], catT[[1]]}]; catT = Transpose@Select[catalog, #[[7]] == 3 &]; caty3 = Transpose[{catT[[3]]/Sin[\[Theta]Vertical], catT[[1]]}]; {{catx2, catx4}, {caty1, caty3}}]
When executed the result is:
In[28]:=
{vxApparent, vyApparent} = ApparentVelocityBeamGrouped[catalog]; gx = ListPlot[vxApparent, Sequence[ PlotRange -> {{-20, 20}, All}, Axes -> False, Frame -> True, AspectRatio -> 1.2, FrameLabel -> {"wind velocity (m/s)", "altitude (m)"}, PlotLabel -> Style[ "\!\(\*SubscriptBox[\(v\), \(x\)]\): East", FontSize -> 10]]]; gy = ListPlot[vyApparent, Sequence[ PlotRange -> {{-20, 20}, All}, Axes -> False, Frame -> True, AspectRatio -> 1.2, FrameLabel -> {"wind velocity (m/s)", "altitude (m)"}, PlotLabel -> Style[ "\!\(\*SubscriptBox[\(v\), \(y\)]\): North", FontSize -> 10]]]; GraphicsRow[{gx, gy}, ImageSize -> 600]
Out[31]=

There is an obvious offset between the beam 4 data and the reversed beam 2 data, see left figure. There is also an offset between the beam 3 data and the reversed beam 1 data, see right figure. These offsets are caused by non-zero vertical winds.

The relationship between wind velocity and the observed Doppler velocity V is

where θ is the beam elevation angle and ϕ is the beam azimuth angle. If vz≠0 then the Doppler velocity will contain a contribution -vzsinθ. This complicates the velocity inversion process.

All the information necessary for a more refined estimate of the wind velocity profile (vx(z),vy(x),vz(z)) is in the data catalog. Suppose that at an altitude of interest the observed Doppler velocities are arranged into the vector . The forward relationship between and the wind velocity is where A is the matrix of beam pointing directions corresponding to the entries in . Specifically, they are the data from columns 4-6 of the data catalog. Provided that A contains data from at least 3 unique beam pointing directions, then the least squares solution for the wind velocity is . In general it is best to require data from all four beams at an altitude before estimating the wind velocity vector. This reduces biasing effects. If θ is the beam elevation angle (80 deg for WiPPR) and ϕ the beam azimuthal angle, then the rows of the matrix A are of the form.

In the case where there are exactly four Doppler velocity measurements (V1,V2,V3,V4) at an altitude on beams (1,2,3,4) respectively then the least squares solution for the wind velocity (vx,vy,vz) takes a very simple form

where θ is the radar beam elevation angle.

The least squares wind vector estimation algorithm is:
In[32]:=
LeastSquaresWindEstimates[catalog_, minUniqueBeams_ : 4] := Module[{zunique, tmp, cat, A, Aunique, Nbeams, V, vxL, vyL, vzL, z, vx, vy, vz}, zunique = Sort@DeleteDuplicates@Transpose[catalog][[1]]; tmp = Reap[ Do[cat = zunique[[nz]] /. GroupBy[catalog, First]; A = Map[Take[#, {4, 6}] &, cat]; Aunique = DeleteDuplicates@A; Nbeams = Length[Aunique]; If[Nbeams == minUniqueBeams, V = Transpose[cat][[3]]; {vxL, vyL, vzL} = -Inverse[Transpose[A] . A] . Transpose[A] . V; Sow[{zunique[[nz]], vxL, vyL, vzL}]], {nz, 1, Length[zunique]}]]; {z, vx, vy, vz} = Transpose@tmp[[2]][[1]]; {Transpose[{vx, z}], Transpose[{vy, z}], Transpose[{vz, z}]}]
The application of the algorithm to the data in the contact catalog lead to the following figure. The least squares solution is shown in black:
In[33]:=
{vxLS4, vyLS4, vzLS4} = LeastSquaresWindEstimates[catalog, 4]; gx = ListPlot[vxApparent, Sequence[ PlotRange -> {{-20, 20}, All}, Axes -> False, Frame -> True, AspectRatio -> 1.2, FrameLabel -> {"wind velocity (m/s)", "altitude (m)"}, PlotLabel -> Style[ "\!\(\*SubscriptBox[\(v\), \(x\)]\): East", FontSize -> 10]], Epilog -> {Point[vxLS4]}]; gy = ListPlot[vyApparent, Sequence[ PlotRange -> {{-20, 20}, All}, Axes -> False, Frame -> True, AspectRatio -> 1.2, FrameLabel -> {"wind velocity (m/s)", "altitude (m)"}, PlotLabel -> Style[ "\!\(\*SubscriptBox[\(v\), \(y\)]\): North", FontSize -> 10]], Epilog -> {Point[vyLS4]}]; GraphicsRow[{gx, gy}, ImageSize -> 600]
Out[36]=
The vertical component of the wind velocity vector is shown in the following:
In[37]:=
gz = ListPlot[vzLS4, Sequence[ PlotRange -> {{-5, 5}, All}, Axes -> False, Frame -> True, AspectRatio -> 1.2, FrameLabel -> {"wind velocity (m/s)", "altitude (m)"}, PlotLabel -> Style[ "\!\(\*SubscriptBox[\(v\), \(z\)]\): up", FontSize -> 10], PlotStyle -> Black], ImageSize -> 250]
Out[37]=

Marshall Bradley, "Nocturnal Jet Wind Profiles Measured with WiPPR" from the Wolfram Data Repository (2022)  

Data Resource History

  • Date Created:

Source Metadata

  • Citation:
    • Wind Profiling Portable Radar (WiPPR) by Marshall Bradley, Mark Henderson, Michael Vollmuth and Quinn Smith, LogLinear Group, January 2022.

Publisher Information

  • Prepared for the Wolfram Data Repository By: Marshall Bradley
  • Publisher of Record: Marshall Bradley