Precipitation Nowcasting with Three-Dimensional Space–Time Extrapolation of Dense and Frequent Phased-Array Weather Radar Observations (original) (raw)

1. Introduction

Conventional weather radars with parabolic antennas have been used to observe precipitation for several decades. They scan the whole sky every 5–10 min at about 15 elevation angles. These data have been used for real-time monitoring and short-range forecasting (nowcasting) of precipitation (e.g., Pierce et al. 2012). Conventional nowcasting techniques of radar precipitation (e.g., Li et al. 1995; Bowler et al. 2006) predict the horizontal translation of precipitating areas at scales of O(100) km up to 1–6 h. This was successful at these spatiotemporal scales because of the quasilinear nature of the horizontal motion. However, severe local storms may develop on short time scales of ~10 min (e.g., Kato and Maki 2009; Kim et al. 2012), and more rapid monitoring is essential for these local storms. For example, Kramar et al. (2005) reported that the Tracking Radar Echoes by Correlation (TREC; Rinehart and Garvey 1978) algorithm failed to compute two-dimensional (2D) motion vectors with intense vertical motions inside supercells.

The phased-array weather radar (PAWR) is a new-generation radar that enables a rapid three-dimensional (3D) volume scan (Zrnic et al. 2007; Bluestein et al. 2010). Recently, a PAWR that can scan the whole sky every 10–30 s at about 100 elevation angles has been developed in Japan (Yoshikawa et al. 2013; Ushio et al. 2015). Although Rinehart and Garvey (1978) have already presented the idea of 3D nowcasting by extrapolation, sufficiently dense and rapid radar observations were not available for many years. To accurately capture a 3D evolution of convective cells at 1–10-km scales in 1–10 min, spatiotemporal resolutions must be <1 km and <1 min. PAWRs satisfy these requirements and would be essential in nowcasting fast-developing local severe storms.

This study aims to demonstrate the advantage of PAWR in terms of 3D motion tracking for precipitation nowcasting at spatial scales of 1–10 km up to several minutes. These spatial and temporal scales are orders of magnitude smaller than those treated by the current operational nowcasting systems based on the 2D motion of conventional parabolic radars. For that purpose, 2D and 3D extrapolation systems are tested (Fig. 1) with both simulated and real PAWR observations. The simulated radar images are free from nonmeteorological echoes and provide the baseline for the proof of concept and a verification of the newly developed system. Experiments with the real radar observations will clarify possible benefits and limitations of the 3D extrapolation system under realistic conditions.

The paper is organized as follows. Section 2 describes the 3D nowcasting system, and section 3 explains the details of the simulated radar images and real PAWR observations. Section 4 describes the experimental settings, and section 5 presents the results. Section 6 provides a discussion, and section 7 concludes the paper.

2. Method

a. Motion vector

A 3D extension of a conventional motion vector algorithm is developed with several modifications for PAWR. In particular, fractional motion vectors are computed to obtain smooth motions because motions computed from frequently updated radar images tend to be discontinuous, as described below.

The new algorithm is based on the TREC algorithm and computes 3D motion vectors

at each

, where

denotes the 3D grid index and t denotes the time. Among various motion vector schemes, TREC is selected because its simple implementation would be suitable for demonstrating the advantage of 3D nowcasting. Instead of the correlation coefficients (CCs) in TREC, the mean absolute difference (MAD) of radar reflectivity (dB_Z_) is used to reduce the computational cost, although either MAD or CC results in similar nowcasting accuracy. The MAD calculates the similarity between two subimages at t = t_0 and t_0 + Δ_t (where Δ_t is the time interval) around a reference position

assuming persistency of the image and translation by

, where δ_i_, δ_j_, and δ_k_ are integer values:

e1

e1

Here

denotes the array of radar reflectivity values of size

(matching box);

denotes the effective number of elements in

after removing missing data (see below); and

,

, and

represent the maximum search distances in the i, j, and k directions, respectively. In this study, motion vectors are calculated at each grid point with overlapped reference boxes.

Usually, a motion vector is calculated simply by

, where the offset

gives the best match, and

denotes the grid spacing. However, this leads to spatial discontinuities because

is discrete and

becomes large for PAWR [e.g., (100 m)/(30 s) = 3.3 m s−1]. Therefore, a weighted mean of possible offsets

is introduced instead of

so that

is fractional whereas

takes only integer values. Taking a weighted mean also filters erroneous local minima of MADs. Here, offsets

with small MADs are used for calculating

:

e2

e2

e3

e3

e4

e4

e5

e5

e6

e6

where the tuning parameter

is chosen to focus on the bottom 10% of the variable range of MAD at

. Using fractional

instead of discrete

produces spatially smoother motion vectors.

The real PAWR data contain missing data because of quality reasons such as ground clutter, shielding, sidelobe contamination, and rain attenuation (Ruiz et al. 2015). When calculating MAD, the missing data are simply ignored. If more than 95% of the matched pairs of pixels are missing at given and , is considered to be missing. If more than 1% of the MADs are missing in the calculation of , is not calculated and is filled with the domain-averaged motion vector at that altitude.

To spatially smooth out random observation noise, Lanczos filtering with a critical scale of three grids is applied in the horizontal and vertical to the motion vectors in the real cases. Finally, the continuity of TREC vectors (COTREC) technique (Li et al. 1995), which constrains the motion vectors by the 2D Boussinesq mass continuity equation, adjusts the horizontal motion vectors to be nondivergent before advecting. Although the flow is not necessarily nondivergent around convective cells, COTREC improves the forecast performance because COTREC avoids erroneous motion vectors as well as small-scale divergent features that are not predictable on the time scales of interest. The vertical motion is not adjusted by COTREC because the nondivergent assumption does not hold.

b. Advection

Reflectivity Z is converted to rain rate after the computation of the motion vectors and before the advection. Although an optimal Z_–_R relationship depends on various factors such as the phase of the hydrometeor, a standard Z_–_R relationship is assumed in this study; optimizing the Z_–_R relationship is beyond the scope of this work. The rain rate is advected by the fifth-order weighted essentially nonoscillatory scheme (WENO; Liu et al. 1994). The motion vector field is also advected by WENO in the idealized cases, and by the first-order upwind scheme in the real cases. This is because WENO outperforms the first-order upwind scheme in the idealized case, and the first-order upwind scheme outperforms WENO in the real cases. The diffusive nature of the first-order upwind scheme smooths out noises and would be beneficial in the real cases. Time integration is performed by the second-order Adams–Bashforth scheme. For rain, the lateral, bottom, and top boundary conditions are fixed in time for the inflow conditions (no rain), whereas for the outflow boundary conditions, zero gradient conditions are imposed. For the motion vector field, the lateral and top boundary conditions are fixed to the mean motion vector at each altitude, and for the bottom boundary condition, a zero gradient condition is imposed to allow downward motions of precipitation at the surface.

3. Data

We selected two heavy rainfall events that occurred on 13 July 2013 in Kyoto (Ruiz et al. 2015) and 11 September 2014 in Kobe, Japan (hereafter, Case2013 and Case2014, respectively). The events were typical localized heavy rainfalls during the summer in Japan and main PAWR targets. However, the two cases differ in several aspects, which will be described in section 3b. We use two different datasets: simulated radar images for Case2013 and the real PAWR observations for the two cases. Both are interpolated to Cartesian coordinates with a grid spacing of Δ_x_ = Δ_y_ = Δ_z_ = 100 m. For the simulated case, the model data are linearly interpolated in the vertical. The real PAWR data are converted by the bilinear interpolation; first, data are interpolated to constant z levels for each azimuth slice, and then interpolated to the regular x, y grid. The time interval is 30 s.

a. Simulated radar images

The Weather Research and Forecasting (WRF) Model, version 3.6 (Skamarock et al. 2008), is used to produce synthetic radar images for Case2013. The model is configured with five one-way-nested computational domains with 62.5-, 12.5-, 2.5-, 0.5-, and 0.1-km horizontal grid spacings. The finest domain is centered at 34.8855°N, 135°E and has 1000 × 600 × 175 grid points, covering the cities of Kobe, Osaka, and Kyoto in Japan (Fig. 2a). The model uses the WRF single-moment 6-class microphysics scheme, the Mellor–Yamada–Nakanishi–Niino level 3 planetary boundary layer scheme, the Rapid Radiative Transfer Model longwave radiation scheme, and Dudhia’s shortwave radiation scheme. The Kain–Fritsch cumulus parameterization scheme is used in the 62.5- and 12.5-km-resolution domains only. The WRF precipitation output is converted to reflectivity using the NCAR Command Language (UCAR/NCAR/CISL/VETS 2014), which assumes spherical particles of constant density with exponential size distributions.

Fig. 2.

Fig. 2.

Fig. 2.

(a) Topography in the 100-m-resolution WRF simulation and the domains of I_2D/3D/P, R1_2D/3D/P, and R2_2D/3D. (b) The WRF-simulated rain rate in the domain of I_2D/3D/P at 0400:00 UTC 13 Jul 2013. The vertical cross section is at y = 32 km and the horizontal plot at z = 2 km. (c) The PAWR rain rate in the domain of R1_2D/3D/P at 0500:10 UTC 13 Jul 2013. The vertical cross section is at y = 17 km and the horizontal plot at z = 2 km. (d) The PAWR rain rate in the domain of R2_2D/3D at 2327:02 UTC 10 Sep 2014. The vertical cross section is at y = 13 km and the horizontal plot at z = 2 km. The bottom-left corner of each nowcasting domain is defined as the origin for each experiment. The dashed lines in (b)–(d) denote the positions of the cross sections.

Citation: Weather and Forecasting 31, 1; 10.1175/WAF-D-15-0063.1

The initial and boundary conditions of the outermost domain are produced from the National Center for Environmental Prediction Final Analysis results. The coarsest domain starts at 0000:00 UTC 13 July 2013. For nowcasting, 61 snapshots of radar reflectivity in the finest domain every 30 s from 0400:00 UTC are used as the input. Two consecutive frames are used to detect motion vectors, so that 60 initial conditions are available for nowcasting.

Figure 2b shows the rain rate in the finest domain at 0400:00 UTC, the first snapshot. A strong convective rainband ranging from west to east and moving southeastward is simulated.

b. PAWR observations

This study uses the data from a 9.4-GHz X-band PAWR implemented at Osaka University’s Suita Campus. This PAWR was developed by Osaka University, Toshiba Corporation, and the National Institute of Information and Communications Technology (NICT), and has been operated since summer 2012 (Yoshikawa et al. 2013; Ushio et al. 2015). The PAWR observes reflectivity, Doppler velocity, and spectrum width (no polarimetric parameters) at 360 azimuthal angles, a 100-m radial resolution, and 110 elevation angles, and scans the 60-km range every 30 s. Ruiz et al. (2015) developed a quality control algorithm for this particular PAWR.

For Case2013, a 70 km × 40 km region northeast of the radar site is used, and the radar is located at 34.8233°N, 135.5229°E, or (x, y) = (10 km, 0 km) in the computational domain (Fig. 2a). Sixty-one images from 0500:10 UTC are used. To avoid erroneous motion vectors resulting from the observation noise, reflectivities >20 dB_Z_ are used in Eq. (1). Figure 2c shows the rain rate at 0500:10 UTC, where the cell that produces the heavy rainfall over Kyoto lies 20–30 km northeast of the radar site.

For Case2014, a 35 km × 25 km region southwest of the radar site is used, and the radar is located outside the computational domain [equivalent to (x, y) = (47.5 km, 27.5 km)]. Thirty-one images from 2317:02 UTC 10 September 2014 (0817:02 Japanese standard time 11 September 2014) are used. Because of different noise characteristics, reflectivities >25 dB_Z_ are used in Eq. (1). An isolated convective system developed rapidly over Kobe and affected people’s lives (Fig. 2d).

The two cases have different characteristics; Case2013 has a rainband associated with an island-scale convergence zone with a clear southwestward motion, whereas Case2014 is an isolated convective system with slower horizontal motions. In addition, convections in Case2013 are deeper and more intense than those in Case2014, leading to a clearer evolution of convective cores (Figs. 2c,d). Case2014 may reveal limitations of the proposed method because Case2014 seems more difficult than Case2013.

4. Experimental setup

Eight experiments, I_2D/3D/P, R1_2D/3D/P, and R2_2D/3D, are performed. I_2D/3D/P are the idealized experiments using the WRF-simulated radar images as the input into the nowcasting system, and R1_2D/3D/P and R2_2D/3D use the real PAWR data. I_2D/3D/P and R1_2D/3D/P are experiments for Case2013, whereas R2_2D/3D are for Case2014. I_P and R1_P represent persistence forecasts. In I_3D, the matching box size b x × b y × b z = 10 × 10 × 10 corresponds to 1 km × 1 km × 1 km, whereas in R1_3D and R2_3D, 40 × 40 × 10 corresponds to 4 km × 4 km × 1 km. In the 2D cases, b z = 1. Motion vectors are extrapolated near the lateral and upper/lower boundaries of the computational domains with a linear tapering because motion vectors are not available there. This does not have much of an effect on the results because we focus on the interior of the domains. In R1_2D/3D and R2_2D/3D, a larger matching box size than that in I_2D/3D is chosen to reduce the effects of small-scale noise. The maximum search distance corresponds to a feature velocity of 16.7 m s−1 in each direction. This value is manually optimized for the lower troposphere, but might be too small in the upper troposphere. The parameter sensitivity is beyond the scope of this paper.

5. Results

This study evaluates the nowcasting accuracy mainly at z = 2 km, following the operational radar-based precipitation products from the Japan Meteorological Agency.

a. Idealized experiment

First, the horizontal motion vectors at 0430:00 UTC computed by TREC/COTREC are compared with the WRF-simulated horizontal winds (Fig. 3). Although apparent motions of rain distributions do not necessarily agree with the actual winds, the horizontal motion vectors are generally consistent with the WRF winds; eastward motions are dominant within the validation domain. The southward winds around x = 54–55 km, y = 23–25 km in the WRF simulation are more consistent with the motions in I_3D than those in I_2D. In addition, I_2D produces northeastward motions around (x, y) = (58 km, 24 km), different from both the WRF-simulated winds and TREC/COTREC motions in I_3D. On the other hand, northeastward winds around y = 21 km in the WRF simulation are inconsistent with the motion vectors in I_2D/3D. Comparing the TREC and COTREC motions, the nondivergent constraint successfully ensures continuous motion vectors, although COTREC tends to increase the southward motions in this validation domain.

Fig. 3.

Fig. 3.

Fig. 3.

(a) The WRF-simulated horizontal winds at z = 2 km, the COTREC motion vectors at z = 2 km in (b) I_3D and (c) I_2D, and the TREC motion vectors at z = 2 km in (d) I_3D and (e) I_2D. The plots are valid at 0430:00 UTC. Color shading represents the rain rate at z = 2 km.

Citation: Weather and Forecasting 31, 1; 10.1175/WAF-D-15-0063.1

Figure 4 shows nowcasting initiated at 0430:00 UTC. In the WRF simulation, a local maximum of rain rate at (x, y) = (55.6 km, 22.9 km) at t = 0 s moves to (56.9 km, 23.3 km) at t = 120 s and intensifies from 94.6 to 126.6 mm h−1. The rain rate further intensifies at t = 240 s. I_3D captures the intensification; the local maximum lies at (56.9 km, 22.6 km) at t = 120 s, which is somewhat to the south of the true position, and the rain rate becomes 108.2 mm h−1. Because we simply advect the rain rate, only the vertical advection can explain the intensification at z = 2 km. In fact, I_2D does not show the intensification; the forecasted peak lies similarly at I_3D, but its intensity decreases in time as a result of numerical diffusion.

Fig. 4.

Fig. 4.

Fig. 4.

Horizontal maps of rain rate at z = 2 km for (top) the WRF simulation, (middle) I_3D, and (bottom) I_2D. Advection starts at 0430:00 UTC and lead times are (left) 0, (center) 120, and (right) 240 s. Vertical cross sections along the white dashed lines are shown in Fig. 5.

Citation: Weather and Forecasting 31, 1; 10.1175/WAF-D-15-0063.1

The vertical cross sections (Fig. 5) illustrate how the intensification occurs. In the WRF simulation, the top half of a convective core centered at (x, z) = (55.7 km, 3.5 km) at t = 0 s moves upward and the bottom half moves downward, and the area of high rain rate crosses z = 2 km (the horizontal dashed lines in Fig. 5). I_3D reasonably captures this motion, resulting in the intensification of rain rate at z = 2 km. This is consistent with the negative values of the vertical component of motion vectors below z = 3.5 km in the convective core. As shown in Fig. 5, the rain rate in the WRF simulation increases with time. However, the reflectivity change within 30 s is limited and does not affect the accuracy of the motion vectors. Incorporating growth/decay of the reflectivity (Li et al. 1995) may further improve the nowcasting performance at longer lead times, although we focus on the advantage of 3D motion tracking in this paper. The vertical winds in the WRF simulation show much stronger updrafts than the vertical motions of reflectivity. An apparent motion of a convective core is the sum of the actual vertical winds and the terminal velocity of hydrometeors. In an early stage of convection, the updraft prevails, whereas in the mature stage, the terminal velocity of hydrometeors starts to prevail. Therefore, advecting rain distributions by the actual 3D winds may fail to capture the correct vertical motions of convective cores.

Fig. 5.

Fig. 5.

Fig. 5.

Vertical cross sections of rain rate along the white dashed lines in Fig. 4 for (top) the WRF simulation, (middle) I_3D, and (bottom) vertical component of motion vectors in I_3D. Advection starts at 0430:00 UTC and lead times are (left) 0, (center) 120, and (right) 240 s. Horizontal dashed lines denote z = 2 km.

Citation: Weather and Forecasting 31, 1; 10.1175/WAF-D-15-0063.1

The advantage of 3D advection is further examined statistically by computing threat scores for the 60 consecutive initialization times. The results show that I_3D outperforms I_2D for two different threshold values of 1 and 10 mm h−1 (Fig. 6). The advantage of I_3D in the 10 mm h−1 threat score is maximized at the 10-min lead time, where I_3D leads I_2D by more than 3 min (30% of the lead time). Clearly, I_2D outperforms I_P, indicating that the 2D COTREC has forecasting skill.

Fig. 6.

Fig. 6.

Fig. 6.

Threat scores for the threshold values of 1 (black) and 10 mm h−1 (red) as a function of the forecast lead time for I_3D (solid), I_2D (dotted), and I_P (dashed–dotted). The analysis domain covers x = 10–100 km, y = 0–55 km, and z = 2 km. Here, Δ_x_ = 100 m. An average of 60 nowcasts are initiated every 30 s.

Citation: Weather and Forecasting 31, 1; 10.1175/WAF-D-15-0063.1

b. Real-case experiment 1: 13 July 2013, Kyoto

In the real PAWR observations, an intense rain area of more than 100 mm h−1 lies at (x, y) = (19.5 km, 18 km) at 0502:10 UTC (Fig. 7a). The 90-s forecasts indicate a clear advantage of 3D (Figs. 7b,c). R1_3D well captures the intense rain area, whereas in R1_2D, the rain rate is at most 60 mm h−1 in that area. The downward motion in R1_3D advects the intense rain area aloft down to z = 2 km, contributing to an accurate prediction of intensifying rainfall there. Clearly, R1_3D and R1_2D differ at around (x, y) = (21 km, 22 km) and (26 km, 20 km). In most cases, R1_3D is closer to PAWR than R1_2D.

Fig. 7.

Fig. 7.

Fig. 7.

Horizontal maps of rain rate at z = 2 km at 0502:10 UTC 13 Jul 2013 for (a) observations, (b) R1_3D, and (c) R1_2D. The lead time in (b) and (c) is 90 s. (d) As in Fig. 6, but for R1_3D (solid), R1_2D (dotted), and R1_P (dashed–dotted). The analysis domain covers x = 10–70 km, y = 0–35 km, and z = 2 km.

Citation: Weather and Forecasting 31, 1; 10.1175/WAF-D-15-0063.1

Figure 7d shows the threat scores for the 60 consecutive initialization times. Similarly to I_2D/3D, R1_3D outperforms R1_2D for the threshold values of both 1 and 10 mm h−1. The advantage of R1_3D in the 10 mm h−1 threat score maximizes at lead times of 3–4 min, where R1_3D leads R1_2D by more than 30 s (16.7% of the lead time). Not surprisingly, the threat scores for R1_2D/3D are generally lower than those for I_2D/3D. Although this can partly be caused by the differences in the computational domain, period, and distributions of the convective cells (cf. I_P and R1_P), R1_3D is closer to the persistence forecasts than I_3D (Figs. 6 and 7d), implying that R1_2D/3D suffer from various limitations in the radar data.

c. Real-case experiment 2: 11 September 2014, Kobe

Here, nowcasting experiments for Case2014 are performed. Figure 8 shows an example of nowcasting by R2_2D/3D starting at 2331:02 UTC 10 September 2014. At z = 2 km (Figs. 8a–d), the PAWR observations show that the rain rate decreases with time around (x, y) = (21 km, 12 km); R2_3D captures it whereas R2_2D misses it. However, R2_3D underestimates the rain rate near (x, y) = (18.5 km, 12.5 km) and (23.5 km, 14.5 km), and R2_2D is closer to the observation. Overall, R2_3D does not show a clear advantage over R2_2D at z = 2 km. By contrast, R2_3D clearly outperforms R2_2D at z = 1.5 km. In the observation (Figs. 8e,f), the rain rate increases around (x, y) = (18.5 km, 12.5 km) and (23.5 km, 14.5 km), and decreases around (x, y) = (21 km, 12 km). R2_3D captures the changes whereas R2_2D fails (Figs. 8g,h). The vertical cross sections in Figs. 8i–k clearly show that R2_3D captures the descending convective core at x = 18.5 km, which leads to the improvement at z = 1.5 km (Fig. 8g). Because the convective core initiates at z = 2 km, R2_3D and R2_2D do not differ much at z = 2 km (Figs. 8c,d).

Fig. 8.

Fig. 8.

Fig. 8.

Horizontal maps of rain rate at z = 2 km for (a) observation at 2331:02 UTC 10 Sep 2014 and (b) observation, (c) R2_3D, and (d) R2_2D at 2332:32 UTC 10 Sep 2014. The lead time in (c) and (d) is 90 s. (e)–(h) As in (a)–(d), but for z = 1.5 km. (i)–(k) Vertical cross sections that correspond to (a)–(c) and (e)–(g). The horizontal lines represent the positions of the cross sections.

Citation: Weather and Forecasting 31, 1; 10.1175/WAF-D-15-0063.1

The advantage of R2_3D over R2_2D is examined using the threat score differences for the 30 consecutive initialization times (Fig. 9). For the two threshold values of 1 and 10 mm h−1, the threat scores at z = 1.5 km are higher than those at z = 2 km except in the first minute for 10 mm h−1. The threat score difference is maximized at lead times of about 2 min at z = 2 km, whereas at about 5 min at z = 1.5 km. This time lag is consistent with the fact that the descending convective cores initiate around z = 2 km. R2_2D outperforms R2_3D at the end of the 10-min nowcasting partly because of the limitations of 3D nowcasting. However, both R2_2D and R2_3D perform poorly because of difficulties in estimating horizontal motions of slowly moving convections. These limitations will be discussed in the next section.

Fig. 9.

Fig. 9.

Fig. 9.

Threat score difference (R2_3D − R2_2D) for the threshold values of 1 (black) and 10 mm h−1 (red) at z = 1.5 (solid) and 2 km (dotted). The analysis domain is x = 0–35 km, y = 0–25 km. An average of 30 nowcasts are initiated every 30 s.

Citation: Weather and Forecasting 31, 1; 10.1175/WAF-D-15-0063.1

6. Discussion

Here, we discuss the possible benefits and limitations of the proposed system by comparing Case2013 and Case2014. From Figs. 4, 5, 7, and 8, it is clear that the 3D nowcasting has a clear advantage when predicting descending convective cores. Improvements in weak precipitation areas are also related to descending motions in many cases. Because the prediction of descending convective cores is essential for early warning of severe local storms, 3D nowcasting will be valuable.

However, the two cases illuminate some limitations of the proposed method. As described in section 3b, convections in Case2014 are shallower than those in Case2013; the echo-top height is around 6 km in Case2014, whereas that in Case2013 is above 10 km (Figs. 2c,d). In particular, strong convective cores initiate around z = 2 km in Case2014 (Fig. 8i). This affects the nowcasting performance above z = 2 km, or more generally, above the altitude of the convective cores; the 3D nowcasting has fewer chances of predicting precipitation growth/decay associated with descending convective cores.

The initial altitude of a convective core may also limit the predictable time scale. When a descending convective core approaches the ground, it becomes difficult to make accurate predictions because of the lower boundary conditions. Thus, the predictable time scale is limited by the initial altitude of the convective core divided by the descending speed. In fact, the threat score improvement by R2_3D does not last beyond 5 min at z = 1.5 km for Case2014 (Fig. 9), whereas I_3D and R1_3D for Case2013 show their advantage longer (Figs. 6 and 7d). Beyond that time scale, 2D extrapolation must be used, or merging the extrapolated field with a high-resolution NWP may be a straightforward choice.

Another difficulty involves the speed of horizontal motion. From Hovmöller analyses, the mean horizontal motion of convective systems (u, υ) in the x and y directions at z = 2 km in Case2014 is estimated as (4.4 m s−1, −1.1 m s−1), equivalent to (δi, δj) = (1.3, −0.33) for Δ_x_ = Δ_y_ = 100 m and Δ_t_ = 30 s. Although a fractional motion is computed to alleviate the discretization problem by Δ_x_/Δ_t_ [Eqs. (2)–(6)], this may affect the accuracy of the estimation; the mean motion vector by the proposed method in R2_3D becomes (u, υ) = (3.4 m s−1, −0.76 m s−1), which is slightly different from the observed motion. In contrast, the eastward motion in Case2013 exceeds 10 m s−1, and it seems more difficult to capture the horizontal motions accurately in Case2014 than in Case2013. To improve the horizontal motion vectors, taking a larger Δ_t_ may help. However, a larger Δ_t_ may degrade the vertical motions because convective cores change rapidly. To summarize, improving the 3D motion vectors for shallow and slowly moving convections will require further investigations. Data assimilation techniques may also help improve motion vectors and rain distributions using noisy observations.

The advantage of 3D nowcasting may also depend on the life stage of convections. In the mature and decaying stages, descending convective cores contribute to the improvements by 3D extrapolation, as shown in Figs. 4, 5, 7, and 8. During the early stages of convective activities, in contrast, it seems more difficult to improve predictions made by 3D nowcasting because descending convective cores do not exist, and predicting upward motions in the reflectivity field may not directly contribute to the prediction of surface precipitation.

Assuming a single Z_–_R relationship for the liquid, solid, and mixed phases may also limit the accuracy of the 3D nowcasting. Motion vectors might be affected by reflectivity patterns produced by the phase change. Estimation of the surface precipitation rate may also depend on the estimation of ice-phase precipitation aloft. Further investigations are needed for more accurate nowcasting.

7. Conclusions

This study has investigated the feasibility of precipitation nowcasting in updrafts and downdrafts inside cumulus convections for up to several minutes by taking advantage of frequent and dense PAWR observations. The 3D extension of a COTREC-based algorithm was developed with several noise reduction techniques for a real case, and comparisons between 2D and 3D motion vectors revealed that the 3D nowcasting outperformed the 2D nowcasting in both idealized and real-case experiments. In particular, the 3D motion vectors captured the downward motion of convective cores and reproduced an enhancement in the rain rate at 2-km altitude. Although cumulus convection is a nonlinear process, the radar reflectivity every 30 s shows almost continuous evolution and enables the quasilinear assumption for up to several minutes. This is promising for very short-range precipitation nowcasting, especially when merged with a superrapid-update numerical weather prediction. However, the present study examined only two rain events within a limited area. The 3D nowcasting technique shows its advantage for fast-moving deep convection patterns, whereas the performance improvement is limited for slowly moving shallow convections. Further statistical consideration is needed to validate the robustness of the results under various conditions.

Acknowledgments

The authors thank Chris Snyder and an anonymous reviewer for their helpful comments. This work was supported by CREST, Japan Science and Technology Agency. The PAWR data are open to public (http://pawr.nict.go.jp/). RK, YK, and YT developed a prototype of the 3D nowcasting system during the RIKEN AICS HPC Computational Science Internship Program in September 2014. The figures were produced by the GFD-DENNOU Library.

REFERENCES