A rain-type adaptive optical flow method and its application in tropical cyclone rainfall nowcasting

Abstract

A rain-type adaptive pyramid Kanade–Lucas–Tomasi (A-PKLT) optical flow method for radar echo extrapolation is proposed. This method introduces a rain-type classification algorithm that can classify radar echoes into six types: convective, stratiform, surrounding convective, isolated convective core, isolated convective fringe, and weak echoes. Then, new schemes are designed to optimize specific parameters of the PKLT optical flow based on the rain type of the echo. At the same time, the gradients of radar reflectivity in the fringe positions corresponding to all types of rain echoes are increased. As a result, corner points that are characteristic points used for PKLT optical flow tracking in the surrounding area will be increased. Therefore, more motion vectors are purposefully obtained in the whole radar echo area. This helps to describe the motion characteristics of the precipitation more precisely. Then, the motion vectors corresponding to each type of rain echo are merged, and a denser motion vector field is generated by an interpolation algorithm on the basis of merged motion vectors. Finally, the dense motion vectors are used to extrapolate rain echoes into 0–60-min nowcasts by a semi-Lagrangian scheme. Compared with other nowcasting methods for four landfalling typhoons in or near Shanghai, the new optical flow method is found to be more accurate than the traditional cross-correlation and optical flow methods, particularly showing a clear improvement in the nowcasting of convective echoes on the spiral rainbands of typhoons.

Keywords

optical flow method / radar echo classification / adaptive / typhoon / nowcasting

1 Introduction

A typhoon is a tropical cyclone (TC) that occurs in the western Pacific Ocean or northern Indian Ocean. Extreme disaster events are frequently caused by the heavy rainfall associated with typhoon landfalls (Chen et al., 2010). Landfalling typhoons can cause strong wind and heavy rainfall/flood damage and lead to human and economic losses (Zhang et al., 2010; Chen et al., 2019; Yu and Chen, 2019). The horizontal structure of a tropical cyclone usually includes an eye and eyewall, which is a band of convective precipitation near the strongest winds that contains several spiral rainbands (Houze, 2010). Rainfall regions associated with a landfalling typhoon can be categorized as follows: core region rain, spiral rainband, rain produced by meso- or microscale systems, unstable rain in the region south of the core where upper-level cold air and lower-layer warm air superimpose each other, rain associated with a squall line, and remote rain (Chen and Li, 2004).

Although tracking and intensity forecasts for TCs have made significant progress in recent years (Emanuel and Center, 2018), the heavy rainfall nowcasting of landfalling TCs remains a challenging issue. Since the 1940s, hurricanes have been initially observed by coastal radars (Wexler, 1945, 1947), and weather radars have become a standard and operational instrument in observing TCs. Using extrapolation methods, radars play an important role in rainfall nowcasting for landfalling TCs. The movement of spiral rainbands in TCs is affected by multiple factors, such as the movement of TCs and convective precipitation evolution, which increase the difficulties of rainfall forecasting. Therefore, it is necessary to research how to obtain accurate motion vectors for different types of precipitation in TCs to improve rainfall nowcasting.

Nowcasting refers to techniques dedicated to making forecasts over relatively short periods, generally with a lead time within 2 h (or 6 h) in China (Yu et al., 2012; Zheng et al., 2015; Yu and Zheng, 2020). At present, radar echo extrapolation remains one of the most important methods for nowcasting. Commonly used extrapolation approaches include centroid tracking, cross-correlation, and optical flow methods. The centroid tracking method (Crane, 1979; Dixon and Wiener, 1993; Johnson et al., 1998) regards storms as three-dimensional cells for recognition, analysis and tracking and implements nowcasting by fitting and extrapolating radar echoes. It is mainly used for storm cell tracking and prediction and includes Storm Cell Identification and Tracking (SCIT) (Johnson et al., 1998) and Thunderstorm Identification, Tracking, Analysis and Nowcasting (TITAN), which track and extrapolate storms linearly by identifying the centroids of a weather system (Dixon and Wiener, 1993). The Tracking Radar Echoes by Correlations (TREC) (Rinehart and Garvey, 1978;  Tuttle and Foote, 1990) approach obtains the motion vectors of radar echoes by calculating the optimal spatial correlation coefficient of different areas on time continuous radar echoes and predicting the position of radar echoes in the future with these motion vectors (Rinehart and Garvey, 1978; Chen et al., 2009). Centroid tracking methods are often used in convective precipitation, while cross-correlation methods can be used for both convective and stratiform precipitation. For locally generated precipitation, where its intensity and shape change rapidly with time, the quality of motion vectors given by cross-correlation methods will decrease, and the tracking failure will increase significantly.

The optical flow method has been used to improve image nowcasting (Gibson, 1950; Lucas and Kanade, 1981; Shi and Tomasi, 1994; Bowler et al., 2004; Cao et al., 2015; Chen et al., 2017; Woo and Wong, 2017). This approach focuses on image changes rather than selecting invariant features to track the motion of an image. When using the optical flow method to calculate motion vectors, the temporal and spatial variation of the image (e.g., radar echo) should be considered (Cao et al., 2015). The optical flow method uses strict constraints in the calculation process. Therefore, an accurate global movement trend of precipitation can still be obtained with the optical flow method, even if the movement and shape of the precipitation system image change rapidly. In this respect, the optical flow method is superior to the cross-correlation method (Han et al., 2008; Cao et al., 2015; Chen et al., 2017). The pyramid Kanade–Lucas–Tomasi (PKLT) optical flow method is better than the cross-correlation methods when tracking the movement of thunderstorms or mesoscale rainbands (squall lines, etc.) whose movement speed changes rapidly with time (Cao et al., 2015; Chen et al., 2017; Woo and Wong, 2017). However, studies have implied that there are disadvantages in the PKLT optical flow method (Cao et al., 2015; Chen et al., 2017). When multiple types of precipitation are mixed together, it is difficult to obtain high-quality motion vectors in various precipitation areas with the PKLT optical flow method. For some tropical precipitation cases in South China, the cross-correlation method performs slightly better than the optical flow method (Cao et al., 2015).

To improve the accuracy and efficiency of precipitation echo nowcasting, based on the traditional PKLT optical flow method and a rain-type classification algorithm for radar echoes, a new optical flow method scheme, a rain-type adaptive PKLT (A-PKLT), is proposed. The organization of this study is as follows. The methodology of the A-PKLT optical flow method is illustrated in Section 2. The data source and the selected typhoon cases are described in Section 3. The rain-type classification algorithm and optical flow parameter optimization are introduced in Section 4. The detailed schemes of optical vector generation using the A-PKLT method are described in Section 5. Comparisons of the A-PKLT, PKLT, and TREC nowcasts for 4 typhoons and a statistical analysis of their performance of nowcast skill are illustrated and discussed in Section 6. Finally, the main findings are summarized and discussed in Section 7.

2 Method

2.1 Limitations of traditional optical flow method

Optical flow is the instantaneous motion of a point in a moving object on an observation image plane. The optical flow method (Gibson, 1950) is designed to find the corresponding relationship between the previous frame and the current frame in the image sequence by finding the change and correlation of points in the time domain. To calculate the motion vectors of objects between adjacent frames, the PKLT optical flow method needs to abide by three assumptions: 1) the observed brightness of any object point is constant over time; 2) image points do not move very fast relative to the frame rate; and 3) the tracked object point keeps the same motion vector as all those of its neighborhood points.

In the calculation, the PKLT optical flow method uses the spatiotemporal gradient function of a time-varying grayscale to calculate the motion vectors of points. First, the grayscale values obtained after preprocessing are used as the input, and the eigenvalue and the corresponding corner point threshold are calculated by the Shi-Tomasi algorithm. The corner point here refers to the intersection of multiple outlines. The gradient of the corner point changes greatly in direction. The position of the corner point can be determined according to the eigenvalue of the point and the threshold value of the corner points. After that, specific basic optical flow equations are established in the neighborhood around the corner point according to the three assumptions of PKLT. The least square method is used to solve the equations to obtain the motion vectors of the corner points.

According to the Shi-Tomasi algorithm, the corner point threshold is determined by the product of the largest point eigenvalue in the image (radar echo) and the quality coefficient, while the point eigenvalue is determined by the grayscale gradient of the point in the west–east and south–north directions. The gradient is proportional to the eigenvalue. When the point eigenvalue is greater than or equal to the corner point threshold, the point is considered the corner point. However, when the gradients of the points (such as radar reflectivity grid points) are not well distributed, the eigenvalues of some points will be relatively large, which will increase the threshold value of corner points. At this time, without reducing the quality coefficient, it is difficult to obtain corner points in the area with small grayscale gradients. For example, in heavy rainfall cases, there are often large gradients in the fringe positions of radar reflectivity and small gradients in the interior of strong echoes. In this case, if the corner points are increased by reducing the quality coefficient (corner point threshold), the accuracy of the optical flow calculation may be decreased. However, the lack of sufficient representative corner points will cause the PKLT optical flow method to fail to calculate the motion vectors in some critical areas (such as convective echoes related to heavy rain). Thus, the effective motion estimation of some parts of the echo area is lost, which will reduce the accuracy and precision of the motion vectors in the entire radar echo. Therefore, it is necessary to increase the number of corner points in the key precipitation areas to obtain more detailed and accurate motions of precipitation echoes.

2.2 Rain-type adaptive optical flow method

It is difficult for traditional PKLT optical flow method to obtain corner points in the areas of radar echoes which have small gradients. The lack of sufficient representative corner points will cause the traditional PKLT optical flow method fails to calculate the motion vectors in some critical areas (such as convective echoes related to heavy rain). A-PKLT optical flow method is able to add corner points in different types of rainfall areas by calculating the corner points of different types of precipitation echoes respectively. This can make up for the deficiency of traditional PKLT optical flow method.

To solve the above problems existing in the traditional PKLT optical flow method, the following improvements are made in the A-PKLT optical flow method. The fundamental procedure is as follows. First, a rain-type classification algorithm is used to classify the types of precipitation areas in radar reflectivity. Six types of rainfall echoes can be identified with the algorithm: stratiform, convective, surrounding convective, isolated convective core, isolated convective fringe, and weak echoes. Then, after extracting different types of precipitation echoes, the gradients of radar reflectivity in the fringe positions corresponding to each type of precipitation are increased. As a result, the number of corner points (characteristic points used for PKLT optical flow tracking) in the surrounding area will be increased. Therefore, more motion vectors can be obtained purposefully for each type of precipitation. It is helpful to describe the motion characteristics of precipitation more precisely.

The design flow chart of the radar echo extrapolation method based on the A-PKLT optical flow method is shown in Fig. 1, which includes “Data preprocessing”, “Rain-type classification”, “Optical flow parameter optimization”, “PKLT motion vector calculation (for all rain types)”, “Motion vector merging”, “Dense motion vector calculation”, and “Radar echo extrapolation”. The procedure of traditional PKLT optical flow method is shown in the dashed rectangle. The difference between the A-PKLT optical flow method and the traditional PKLT optical flow method is shown in Fig. 1.

3 Data

3.1 Radar data and preprocessing

The high-resolution reflectivity and correlation coefficient of the Nanhui Weather Surveillance Radar-1988 Doppler (WSR-88D) dual polarization Doppler weather radar and the reflectivity data of the Qingpu WSR-98D Doppler radar data are used in the algorithm. The level II data include reflectivity, mean radial velocity, and spectrum width for both radars and the differential reflectivity, correlation coefficient, and differential phase from the Nanhui WSR-88D. The high-resolution reflectivity reflects the scale and density distribution of precipitation particles inside the meteorological target, which is used to represent the intensity of the precipitation (Crum et al., 2013). The correlation coefficient (CC) represents the correlation between the detection pulse changes of the horizontal channel and vertical channel on the dual polarization Doppler weather radar, and the value ranges from 0 to 1 (Yang et al., 2017).

Radar data preprocessing process includes radar data acquisition, data quality control, coordinate conversion and interpolation, grayscale conversion, and denoising processing. Radar data acquisition incorporates the reflectivity data of the 0.5° elevation tilt at the current time (T_{\mathrm{0}}), and the previous scans (configurable, generally 1–3 scans before T_{\mathrm{0}}) are processed.

A data quality control method based on the correlation coefficient is used to eliminate the nonprecipitation echo, but only for WSR-88D. The correlation coefficient can provide a measure of the consistency of the shapes and sizes of targets within the radar beam. A higher value shows a higher consistency in the size and shape of radar targets, while a lower value indicates greater variability in shapes and sizes. According to research results indicating that in precipitation including hail, the correlation coefficient seldom falls below 0.9 (Dusan et al., 2006), nonprecipitation radar echoes are eliminated if the CC value is lower than 0.9, and the corresponding reflectivity data are set as an invalid value of 0. For the Qingpu WSR-98D Doppler radar, the scattered and isolated echoes are filtered, and then the ground clutter echoes with low height are removed through a comparison of the reflectivity of the lowest three elevation angles and their near-zero radial velocity.

The quality-controlled reflectivity data are converted from polar coordinates to rectangular coordinates. The four closest reflectivity data in both the radial and azimuth directions are used for interpolation with an inverse distance weighting method. As for grayscale conversion, all radar reflectivity data after quality control are normalized and converted to grayscale values (value range is 0–255) using Eq. (1):

(1)

where ref is the reflectivity value that should be converted and ref\text{'} is the output reflectivity grayscale value. Symbol min means to take the minimum value. When the value of ref is less than 1, ref\text{'} is set to zero.

As for denoising, using the same structural element, the image is processed with an erosion operation and again with a dilation operation: this is called an open operation. The open operation is typically used to eliminate small noise points. Using an open operation, the isolated echo points are removed from the radar reflectivity data. All radar reflectivity grayscale data are processed with erosion operations and then processed with dilation operations to remove isolated echo points with a denoising Eq. (2):

(2)A\circ \mathrm{B}=\left(A\Theta B\right)\oplus B,

where A is the radar reflectivity grayscale data and B is the structural element. \Theta is the symbol for the erosion operation. \oplus is the symbol for the dilation operation. The structural element can be generally defined as a rectangular structure, elliptical structure or cross structure. The size of the structural element can also be customized. Different custom parameters will affect the denoising effect. The elliptical structure is used as the structural element, and the size is 3×3.

3.2 Typhoon cases

Four typhoon cases are used to verify the echo extrapolation nowcast based on the new optical flow method: Ampil (No. 1810), Jongdari (No. 1812), and Rumbia (No. 1818), which made landfall in Shanghai in 2018, and Lekima (1909), which made landfall in Zhejiang Province in 2019. In Table 1, the date and time, location, center pressure, maximum wind speed near the center at landfall, maximum rainfall and maximum hourly rainfall in Shanghai, and data period of all typhoons in this study are given.

Ampil (No. 1810) made landfall in eastern Chongming, Shanghai, at 0430 UTC 22 July, 2018 with a maximum total rainfall of 118.1 mm and a maximum hourly rainfall of 46.5 mm on the east coast of Shanghai. A total of 317 scans of the Shanghai Nanhui WSR-88D radar are selected for the period of 1107 UTC 21 July – 1759 UTC 22 July, 2018 when Ampil was close to Shanghai.

After making an initial landfall in Japan, Jongdari (No. 1812) turned westward into the East China Sea and made landfall on the coast of Jinshan District, Shanghai, at 0230 UTC 03 August, 2018. Jongdari brought heavy rain in central and southern Shanghai, with a maximum rainfall of 104.4 mm in Songjiang District within 20 h and a maximum wind gust of 26.4 m/s in southern Pudong District. A total of 356 scans of the Shanghai Nanhui WSR-88D radar are selected for the period of 0204 UTC 2 August, – 1055 UTC 3 August, 2018.

Rumbia (1818) made landfall on the southern coast of Shanghai Pudong New District in the early morning of 17 Aug. From 0000 UTC 16 August, to 0700 UTC 17 Aug; Shanghai received heavy rain, a maximum rainfall of 154.4 mm in Minhang and a maximum hourly rainfall of 46.6 mm in Xujing, Qingpu (0500–0600 UTC 17 Aug). A total of 407 scans of the Shanghai Qingpu WSR-98D radar are used for the period of 00030 UTC 16 August – 1706 UTC 17 August, 2018.

Lekima (No. 1909) was the strongest landfalling typhoon in China in 2019, causing severe wind and rain impacts in East China (Zhou and Gao, 2019). Lekima made landfall on the coast of Wenling City, Zhejiang Province at approximately 1745 UTC on 9 Aug. The maximum wind speed near the center at landfall was 52 m/s (super typhoon), and the lowest center pressure was 930 hPa. As a supertyphoon, Lekima ranked fifth among the landfalling typhoons in Chinese mainland since 1949. Lekima moved slowly after landfall and lasted for 44 hours on land, including Zhejiang Province for 20 hours, making it the longest supertyphoon in Zhejiang Province. Lekima received a maximum total rainfall of 833 mm in Kuocang Mountain, which was the second heaviest daily rainfall record for landfalling typhoons in Zhejiang. Lekima caused disasters such as floods in urban and rural areas, floods in small and medium-sized rivers, flash floods, and landslides in several provinces.

The overall rainfall in Shanghai exceeded 150 mm, while Minhang Pujiang recorded a maximum rainfall of 276.1 mm from 0000 UTC 9 August, to 0100 UTC 11 August. It is not unusual for spiral rainbands in tropical cyclones to become more convective, with echoes reaching intensities comparable to those of thunderstorms (Fabry, 2015). One of the spiral rainbands of Lekima brought a maximum hourly rainfall of 112.2 mm from 0600 to 0700 UTC on 10 August, at Minhang Pujiang (Fig. 2(a)).

A total of 481 scans of the Shanghai Nanhui WSR-88D radar are selected for the period of 0000 UTC 9 August, – 0000 UTC 11 August, 2019 when Lekima was making landfall and moving from Zhejiang Province to Jiangsu Province (Fig. 3).

4 Echo classification and optical flow parameters

4.1 Rain echo classification

Steiner et al. (1995) (hereafter SHY95) developed algorithms to separate radar echoes into convective and stratiform regions, statistically summarize the vertical structure of radar echoes, and determine precipitation rates and amounts at high spatial resolution. Then, the SHY95 algorithm was improved to classify the type of precipitation areas in quality-controlled radar reflectivity data by Powell et al. (2016) (hereafter PHB16). Six types of rainfall echoes could be identified with the algorithm: stratiform, convective, surrounding convective, isolated convective core, isolated convective fringe, and weak echoes. Among them, surrounding convective echoes are defined as “uncertain” echoes surrounding convective cores, and the isolated convective core is the strongest echo in small echo objects. It often represents the cores of developing shallow and isolated convection, and the isolated convective fringe echoes are related to the weaker echoes in small echo objects. This can include weak, decaying convection and echoes surrounding isolated convective cores. Based on the rain-type classification algorithm (decision tree) shown in Fig. 5 of PHB16, the precipitation echo regions are classified into 6 types of rain echoes: convective, stratiform, surrounding convective, isolated convective core, isolated convective fringe, and weak echoes.

4.2 Optical flow parameter optimization

Since the optical flow equation only retains the first-order approximation term after expansion with Taylor’s formula, when the object to be tracked moves faster, the sum value of higher-order terms is no longer approximate to zero. As a result, the accuracy of motion vector calculation will be decreased. If the corner point to be tracked between two frames moves so fast that it is out of the optical flow search window, the accuracy of motion vector calculation will be decreased. With a pyramid layered method, the target image will be decomposed layer by layer. The image with the lowest resolution is placed on the top layer, and the original image is placed on the bottom layer. The motion vectors are modified from the top layer to the bottom layer. In this way, the motion vectors can be estimated more accurately, even if the target moves fast. Generally, two or three pyramid layers are used in the PKLT method. For typical image sizes, it makes no sense to go above four pyramid layers (Bouguet, 2000).

Considering the motion characteristics of the radar echo, a reasonable motion speed range should be set to filter out unreliable motion that is too fast or too slow. Stratiform precipitation is often affected by advection, and its moving speed is relatively stable, while convective precipitation is affected not only by advection but also by rapid changes in its intensity (Dai, 2013; Sun et al., 2015). According to an analysis of the historical radar echo data of the WSR-88D radar in Shanghai Nanhui, it is generally considered that the maximum distance a precipitation echo moves between two radar scans (approximately 6 min) is 10–15 km.

1) Number of pyramid layers

When the echo moves fast, the accuracy of the PKLT optical flow method can be improved by using a higher number of pyramid layers. The following empirical parameters are selected here: three-layer pyramids are used for convective and surrounding convective echoes, and two-layer pyramids are used for stratiform, isolated convective core, and isolated convective fringe echoes.

2) Speed threshold of motion vectors

According to the differences in stratiform, convective, surrounding convective, isolated convective core, isolated convective fringe, and weak echoes, more reliable motion vectors are obtained by setting different velocity ranges for motion vectors. Eq. (3) is used to calculate the velocity (distance) Vcorresponding to each corner point separately.

(3)V=\sqrt{u^{\mathrm{2}}+v^{\mathrm{2}}},

where u,v are the motion vectors of the corner point in the west–east and south–north directions.

Equation (4) is used to determine the maximum and minimum velocity thresholds of the motion vectors corresponding to the corner points in various precipitation areas, which are used to filter out the corner points with motion speeds that are too fast and too slow.

(4)\begin{array}{c}\hfill V_{\max }=\min (T_{\max },\hspace{0.17em}q\mathrm{2}+\mathrm{3}\times (q\mathrm{3}-q\mathrm{1}\left)\right)\hfill \\ \hfill V_{\min }=\max (\mathrm{0},\hspace{0.17em}q\mathrm{2}-\mathrm{2}\times (q\mathrm{3}-q\mathrm{1}\left)\right),\hfill \end{array}

where V_{\max }, V_{\min } and T_{\max } represent the maximum velocity threshold, the minimum velocity threshold, and the maximum velocity limit constant, and q\mathrm{1}, q\mathrm{2}, and q\mathrm{3} are the first, second, and third quartiles of velocity for all corner points. When the velocity V of a corner point is greater than V_{\min } and less than V_{\max }, the corresponding motion vectors u,v and corner point coordinates (x,y) are reserved.

3) Parameter comparisons

The number of pyramid layers and the T_{\max }value should be within a reasonable range. Considering the image size and reasonable moving speed of a radar echo, the number of pyramid layers should be 2 or 3, and the T_{\max }value should be between 10 and 15. According to the test, if the parameters are selected within this range, the calculation results of the optical flow method will be relatively reliable.

Two parameter schemes are tested and evaluated in different types of precipitation systems (e.g., spiral rainbands, mesoscale convective systems, and afternoon local convective storms). These parameters are used to optimize the motion vector calculation results of A-PKLT optical flow method.

Parameter scheme 1: For convective regions and surrounding convective regions, a 3-layer pyramid and T_{\max }=\mathrm{15} are selected. For stratiform regions, a 2-layer pyramid and T_{\max }=\mathrm{10} are selected. For isolated convective cores and isolated convective fringe regions, a 2-layer pyramid and T_{\max }=\mathrm{12} are used.

Parameter scheme 2: A 2-layer pyramid and T_{\max }=\mathrm{10} are selected for all types of precipitation echoes.

The results show that the extrapolation accuracy of parameter scheme 1 is slightly higher than that of parameter scheme 2. Therefore, parameter scheme 1 is utilized in the A-PKLT optical flow method.

5 Motion vectors

First, based on the rain-type classification results, the reflectivity grayscale data are extracted iteratively by rain type, and five types of rain echoes (except for weak echoes) are extracted. Then, different optical flow calculation parameters and motion vector filtering parameters are used for different types of precipitation echo areas to calculate the motion vectors iteratively. When all the iterations are completed, the motion vectors corresponding to all five types of rain echoes are merged. After declustering and interpolation, the A-PKLT dense motion vectors are obtained. The motion vector procedure includes the following steps.

5.1 Corner points

Five types of classified rain echoes, including stratiform, convective, surrounding convective, isolated convective core, and isolated convective fringe echoes, are used for corner point calculations. The radar reflectivity grayscale data of the first N frames and the first N+1 frames before the current time (T_{\mathrm{0}}) are extracted according to each type of rain echo. Then, the Shi-Tomasi algorithm is used to recognize and detect the corner point coordinates in the radar reflectivity grayscale data of the first N frames, and the corner point coordinates in the reflectivity grayscale data of the first N frames before the current time (T_{\mathrm{0}}) are iteratively calculated for all types of rain echoes.

The corner point detection process is as follows: 1) input the radar reflectivity grayscale data of the first N frames, 2) set up a window to slide on the image plane, 3) calculate the grayscale change of points in the sliding window (15 px×15 px) using Eq. (5), and 4) calculate the positions of corner points according to the grayscale change of points in the sliding window.

(5)E(u,v)={\displaystyle \sum _{x,y}w(x,y)}{\left[I\right(x+u,y+v)-I(x,y\left)\right]}^{\mathrm{2}},

where (u,v) represents the offset, (x,y) is the point coordinate within the window, w(x,y)is the window function, and E(u,v) represents the grayscale change.

After using Taylor’s formula to expand Eq. (5):

(6)E(u,v)\cong [u,v]M\left[\begin{array}{c}\hfill u\hfill \\ \hfill v\hfill \end{array}\right],

(7)M={\displaystyle \sum _{x,y}w(x,y)\left[\begin{array}{cc}\hfill I_x{I}_x\hfill & \hfill I_x{I}_y\hfill \\ \hfill I_x{I}_y\hfill & \hfill I_y{I}_y\hfill \end{array}\right]},

where M is the covariance matrix, and I_x and I_y are the gradients of the point in the west–east and south–north directions.

The corner detector response function is defined as:

(8)R=\min ({\lambda }_{\mathrm{1}},{\lambda }_{\mathrm{2}}),

where {\lambda }_{\mathrm{1}} and {\lambda }_{\mathrm{2}} are the eigenvalues of the covariance matrix corresponding to the point. This is defined as the change component of the two orthogonal directions after the matrix is diagonalized. R is the minimum value of {\lambda }_{\mathrm{1}} and {\lambda }_{\mathrm{2}}, and R is the eigenvalue of the point.

The quality coefficient (generally 0.01–0.1) is used to determine the minimum quality level of the acceptable corner points. The quality coefficient is multiplied by the maximum R value of all points as the corner point threshold, and points with R values greater than or equal to the corner point threshold are considered corner points.

5.2 Motion vector estimation

As mentioned in Section 2.1, the PKLT optical flow method needs to abide by three assumptions. From assumption 1, the following equation can be established:

(9)I(x,y,t)=I(x+\mathrm{d}x,\hspace{0.17em}y+\mathrm{d}y,\hspace{0.17em}t+\mathrm{d}t).

When the right part of Eq. (9) is expanded by a Taylor formula, and the high-order terms are ignored:

(10){\displaystyle \frac{\partial I}{\partial x}}u+{\displaystyle \frac{\partial I}{\partial y}}v+{\displaystyle \frac{\partial I}{\partial t}}=\mathrm{0.}

From assumption 3, assuming that the tracked object point and all points in its neighborhood maintain the same motion vector, the equation is established as:

(11)\begin{array}{c}\hfill I_x\left(q_{\mathrm{1}}\right)u+I_y\left(q_{\mathrm{1}}\right)v+I_t\left(q_{\mathrm{1}}\right)=\mathrm{0}\hfill \\ \hfill ⋮\hfill \\ \hfill I_x\left(q_n\right)u+I_y\left(q_n\right)v+I_t\left(q_n\right)=\mathrm{0},\hfill \end{array}

where I_x\left(q_n\right), I_y\left(q_n\right), and I_t\left(q_n\right) are the partial derivatives of the point, which are solved by the least squares method:

\left[\begin{array}{c}\hfill u\hfill \\ \hfill v\hfill \end{array}\right]={\left[\begin{array}{cc}\hfill {\displaystyle {\sum }_n{I}_x{\left(q_n\right)}^{\mathrm{2}}}\hfill & \hfill {\displaystyle {\sum }_n{I}_x\left(q_n\right)I_y\left(q_n\right)}\hfill \\ \hfill \hfill \\ \hfill {\displaystyle {\sum }_n{I}_y\left(q_n\right)I_x\left(q_n\right)}\hfill & \hfill {\displaystyle {\sum }_n{I}_y{\left(q_n\right)}^{\mathrm{2}}}\hfill \end{array}\right]}^{-\mathrm{1}}

(12)\left[\begin{array}{c}-{\displaystyle \sum {}_n{I}_x\left(q_n\right)I_t\left(q_n\right)}\\ \hfill -{\displaystyle \sum {}_n{I}_y\left(q_n\right)I_t\left(q_n\right)}\hfill \end{array}\right].

Finally, the motion vectors of the corner points are calculated.

5.3 Motion vector merging

The motion vectors u,\hspace{0.17em}v and corner point coordinates (x,y) of five types of precipitation echoes are merged iteratively into four arrays in which the values of u,\hspace{0.17em}v,\hspace{0.17em}x,\hspace{0.17em}y are stored separately.

5.4 Dense motion vectors

1) Declustering

After the motion vectors of all types of rain echoes are merged, a box for declustering calculation (a default box of 20 px×20 px) is established. Only the median of motion vector (u,v) and its coordinates (x,y) in one grid will be stored (Fig. 4).

Through the methods above, more reliable corner points and corresponding motion vectors are obtained. Precipitation echoes often show large reflectivity gradients in the fringe regions and small reflectivity gradients in the interior regions of the precipitation system. Due to the limitations of the traditional PKLT corner detection algorithm, it is difficult to obtain the corner points in the interior regions with small reflectivity gradients, resulting in fewer motion vectors in the interior areas that are always associated with heavier rainfall. Figure 5 shows the sparse motion vectors calculated by the A-PKLT optical flow method (b) and the traditional PKLT optical flow method (a) for Typhoon Lekima. More motion vectors in the interior areas of two spiral rainbands that have strong echoes are obtained when using the A-PKLT optical flow method than when using the traditional PKLT optical flow method. A more precise description of radar echo motion characteristics will be obtained in this region when calculating motion vectors with the new scheme.

2) Interpolation

As the result of the PKLT optical flow method is sparse optical flow, the reliable corner points and motion vectors obtained in the previous steps are used as the original data for the interpolation of motion vectors to obtain the dense motion vectors. The inverse multiquadric radial basis function is used to calculate the interpolation weights, and the K-nearest neighbor algorithm (KNN) is used to speed up the interpolation (only the nearest K motion vectors are taken for calculation). The weight calculation function is as follows:

(13)w={\displaystyle \frac{\mathrm{1}}{\sqrt{\mathrm{1}+{\left(ϵr\right)}^{\mathrm{2}}}}},

where ϵ is the adjustment factor and r is the Euclidean distance from the interpolated point to the corresponding corner point.

The advantage of inverse multiquadric radial basis function weight interpolation is that it can obtain smooth motion vectors when the original motion vectors do not change suddenly in a short distance. After being calculated by Eq. (4) in Section 4.2, only the original motion vectors with reasonable speed are retained. Under this premise, the interpolated motion vector field can be considered as smooth after using the KNN algorithm and inverse multiquadric radial basis function weight interpolation to interpolate.

5.5 Echo extrapolation

Finally, with dense motion vectors, the radar echo is extrapolated by a semi-Lagrangian scheme (Germann and Zawadzki, 2002; Zhang et al., 2014), which can solve the problem whereby the echo cannot be rotated when extrapolating linearly. The predicted radar echo \Psi at time t_{\mathrm{0}}+\tau (initial+forecast time) is given by a semi-Lagrangian advection formula:

(14)\psi (t_{\mathrm{0}}+\tau ,x)=\Psi (t_{\mathrm{0}},x-\sigma ),

where \Psi is the observed radar echo and \sigma is the motion vector within the prediction time.

The extrapolation process within time \tau is divided into K time steps {\mathrm{\Delta }}_t. For each time step, the motion vector can be calculated by the following formula:

(15){\alpha }^{(K+\mathrm{1})}={\mathrm{\Delta }}_{t}u(t_{\mathrm{0}},x-{\displaystyle \frac{{\alpha }^{\left(K\right)}}{\mathrm{2}}}),

where u(t_{\mathrm{0}},x-{\displaystyle \frac{\alpha }{\mathrm{2}}}) is the velocity of the echo at x-{\displaystyle \frac{\alpha }{\mathrm{2}}}, and K is the number of iterations. The initial value \alpha is zero, and the final motion vectors are the sum of multiple vectors in every single time step.

6 Results

Four typhoon cases are used to verify the echo extrapolation nowcast based on the new optical flow method using the critical success index (CSI) score. The CSI is a biased score that is dependent upon the frequency of the event (Schaefer, 1990).

(16)CSI={\displaystyle \frac{Ns}{Ns+Nf+Nm}},

where Ns, Nf, and Nm represent a hit, false, and miss, respectively, defined by events being predicted and verified by the observation, events being predicted but not verified by the observation, and events being not predicted but verified by the observations. Here, the CSI is used as an evaluation score for radar echo extrapolation based on the A-PKLT, PKLT, and TREC methods within a lead time ranging from 6 to 60 min. The higher the CSI is, the better the extrapolation method is.

Figure 6 shows the average CSI of the reflectivity nowcasts at each reflectivity level of 20–55 dBZ within the 6–60-min lead time of the A-PKLT, PKLT and TREC methods during the influence of the four typhoons: a) Ampil, b) Jongdari, c) Rumbia, and d) Lekima. In general, 1) as the nowcasted reflectivity increases, all methods’ CSI scores for the four typhoons continue to decrease. For example, in Typhoon Ampil, the average CSI of the 6– 60-min reflectivity nowcast of the three methods drops from 0.42–0.53 at 20 dBZ down to 0.23–0.32 at 35 dBZ and approximately 0.02 at 50 dBZ. 2) Both PKLT and A-PKLT perform significantly better than TREC, which is similar to some previous studies in which the optical flow method is better than the cross-correlation methods when tracking the movement of mesoscale rainbands whose movement speed changes rapidly with time (Cao et al., 2015; Chen et al., 2017; Woo and Wong, 2017). At higher reflectivity levels (above 30 dBZ), A-PKLT is slightly better than PKLT. 3) In different typhoon cases, the three methods perform differently, but all show their best performance in Typhoon Lekima. Among them, the A-PKLT performs the best. The mean CSI of the A-PKLT’s 6–60-min 20 dBZ reflectivity nowcast for Typhoon Lekima reaches approximately 0.7, 0.44 for 35 dBZ, and 0.28 for 40 dBZ, which are significantly higher than those for the other three typhoon cases. The largest gap between the PKLT-based methods (PKLT and A-PKLT) and TREC, however, is found in Typhoon Ampil, while the smallest gap is found in Typhoon Lekima.

Figure 7 shows the 0.5° elevation angle, (a) the reflectivity observed by the Nanhui WSR-88D radar at 0123 UTC on 22 July, 2018, and the 60-min reflectivity extrapolation at 0027 UTC on 22 July, 2018 using three extrapolation methods: (b) TREC, (c) PKLT, and (d) A-PKLT for Typhoon Ampil. A-PKLT's reflectivity nowcast for Typhoon Ampil at this moment is more accurate (Fig. 7(d)) than TREC and PKLT, showing a strong echo region of 40–50 dBZ over the radar site, which is closer to the observation (Fig. 7(a)), while the strong echo regions of PKLT and TREC’s nowcasts are still offshore.

Compared with Typhoon Ampil with only one major spiral rainband (Fig. 7(a)), Jongdari’s precipitation echoes (Fig. 8(a)) are relatively scattered, and no major spiral rainband is found, while Rumbia (Fig. 8(b)) has two major spiral rainbands and larger rain echo coverage than the other two typhoons. On the other hand, according to Fig. 6, the PKLT and A-PKLT nowcasts for Typhoon Rumbia are better than those for Ampil and Jongdari, while the TREC nowcasts show similar performances in both Jongdari and Rumbia. This may be related to the larger and stronger precipitation echoes of Rumbia (Fig. 8(b)) than those of the other two typhoons.

Table 2 shows the average CSI scores of the 6–60-min (lead time) A-PKLT 20–55 dBZ reflectivity nowcast for Typhoon Lekima. In general, the CSI scores decrease with time, similar to the other three typhoons in 2018. The CSI of 30 dBZ continues to decrease from 0.75 at 6 min to 0.43 at 60 min, and the CSI of 45 dBZ decreases from 0.36 at 6 min to 0.08 at 60 min. As the forecasted reflectivity increases, the CSI score also decreases. Taking the 30-min lead time forecast as an example, the CSI drops from 0.71 at 20 dBZ and 0.44 at 35 dBZ to 0.05 at 50 dBZ. Compared with the three other typhoons (Ampil, Jongdari, and Rumbia) in 2018, A-PKLT’s nowcasts for Typhoon Lekima (No. 1909) are even better, which may be related to the larger and stronger precipitation area of Lekima (Fig. 9(a)) than the other three typhoons. For all 20– 40 dBZ and lead times over 36 min of 45–50 dBZ nowcasts, the CSI scores were higher than those of typhoons in 2018. This may be associated with the larger precipitation echo region, the stronger intensity and longer lifetime related to the convective echo in Lekima.

Table 3 gives the CSI difference of A-PKLT relative to PKLT in the case of Typhoon Lekima. Compared with PKLT, most of the CSI scores of A-PKLT have an average increase of above 0.02. In the 6–30-min lead time nowcast for strong reflectivity (40–50 dBZ), the CSI has an increase of 0.04–0.06, which is better than that for the other three typhoons in 2018.

Based on the comparisons among methods for four typhoons, the optical flow-based methods (PKLT and A-PKLT) perform significantly better than the traditional cross-correlation method (TREC), especially for those typhoons whose spiral rainband structure is clear and whose other types of precipitation echoes are smaller (e.g., Ampil and Rumbia). When there is a larger range of other types of precipitation echoes rather than nonspiral rainbands, the advantages of optical flow methods against cross-correlation methods are reduced (e.g., Jongdari and Lekima). On the other hand, when making nowcasts for typhoons with two or more long and strong spiral rainbands (e.g., Rumbia and Lekima), the advantages of the A-PKLT algorithm over PKLT gradually become prominent, especially for the convective region with reflectivity over 35 dBZ. Since the main factor affecting extrapolation-based nowcasting is the change in precipitation system intensity, when the spiral rainband is increasingly longer, the rainband’s major characteristics may be slightly longer than those of the smaller and weaker spiral rainbands. Because the A-PKLT algorithm is designed to optimize the recognition of convection regions based on optical flow technology, the A-PKLT rain-type adaptive extrapolation method can accurately identify and extrapolate longer and stronger spiral rainbands in typhoons.

To further explore the advantages of the A-PKLT algorithm, its performance at different stages of Typhoon Lekima is analyzed. When calculating the statistics of CSI performance during the major spiral rainband developed and arriving at Shanghai from 0423 UTC to 0807 UTC on 10 August, the average CSI scores of A-PKLT (Fig. 10(a)) for 45–55 dBZ are significantly higher than the 2-day average, with CSI scores of 0.31 to 0.24 for 45–50 dBZ nowcasts, which are much higher than those for the entire study period of Lekima. Based on Lekima’s radar images (Fig. 3(e)), this period is when the major spiral rainband in Shanghai developed in its mature stage. With the same period as Fig. 10(a), the lead time-based CSI scores of the 45 dBZ and 50 dBZ 60-min forecasts reach 0.20 and 0.13, respectively, much higher than the average of 0.08 and 0.02 over the two days in Lekima and those for the three typhoons in 2018. At same time, the CSI gains of A-PKLT over PKLT (Table 4) are also much higher than the 2-day average (in Table 3), mainly at higher reflectivity levels (40–55 dBZ). The maximum CSI gain (0.09) is found within the 24–36 min 50 dBZ nowcasting, revealing that A-PKLT performs better for those strong and long-lived rainbands.

Figure 10(b) shows hourly average CSI gains of A-PKLT relative to PKLT of 40, 45 and 50 dBZ reflectivity nowcasts during 0000 UTC 9 August – 2300 UTC 10 August. According to the reflectivity patterns shown in Fig. 3, the A-PKLT method has a greater improvement over the traditional PKLT optical flow method: 1) when Typhoon Lekima shows a spiral rainband structure, especially when there is a strong convective spiral rainband (reflectivity reaching 40–50 dBZ); 2) when the spiral rainband strong echo zone shows a significant different movement speed or/and moving direction with typhoon overall rain zone, e.g., 0600–1800 UTC of 10 August, 2019.

7 Discussion and Conclusions

In this paper, a rain-type adaptive pyramid Kanade–Lucas–Tomasi (A-PKLT) optical flow method for radar echo extrapolation is introduced. After using the A-PKLT rain echo classification algorithms, six types of rain echoes: stratiform, convective, surrounding convective, isolated convective core, isolated convective fringe, and weak echo can be identified. New schemes are designed to adaptively optimize some parameters of the PKLT optical flow based on rain-type radar echoes. By giving the specific motion speed threshold and pyramid layers to each rain type and increasing the gradients of radar reflectivity in the fringe positions, more corner points and more motion vectors are captured. It is helpful to describe the motion characteristics of precipitation (such as convective precipitation regions) more precisely.

By comparing the 6–60-min precipitation echo nowcast for a series of typhoon cases that made landfall in or passed by Shanghai in recent years, the A-PKLT optical flow method is verified to be more accurate than the traditional cross-correlation and traditional optical flow methods. The comparisons reveal that the optical flow methods (PKLT and A-PKLT) are significantly better than the traditional cross-correlation method (TREC) for typhoons whose spiral rainband structure is clear and whose other types of precipitation echoes are smaller (e.g., Ampil and Rumbia). When there is a larger range of other types of precipitation echoes rather than nonspiral rainbands, the advantages of optical flow methods are reduced (e.g., Jongdari and Lekima). On the other hand, when providing nowcasts for typhoons with two or more long and strong spiral rainbands (e.g., Rumbia and Lekima), the advantages of the A-PKLT algorithm over PKLT gradually become prominent, especially for reflectivity above 35 dBZ. Because the A-PKLT algorithm is designed to optimize the recognition of convection regions based on optical flow technology, this rain-type adaptive optical flow extrapolation method can effectively improve nowcasting skills for heavy rainfall related to tropical cyclones, especially for strong and long-lived convective spiral rainbands.

The A-PKLT optical flow method has shown superior performance compared with other extrapolation-based methods. However, there are many factors that lead to heavy rainfall, such as extratropical cyclones, Meiyu fronts, mesoscale convective systems, or local convective weather in the early afternoon. Optical flow parameter optimization schemes require more precipitation systems to verify and improve results. On the other hand, as an extrapolation-based method, its performance decreases rapidly with lead time (after 30 min) because it does not account for the initiation, growth, and dissipation of precipitation patterns and the uncertainty of the displacement (Browning and Collier, 1989; Wilson, 2004; Dai, 2013). Indeed, it is impossible for this optical flow method to predict the initiation, growth, and decay of precipitation due to a lack of physical mechanisms, which is also a common defect for all extrapolation algorithms. To overcome this issue, in the future, some new approaches, such as extrapolation nowcast and downscaled numerical weather prediction model forecast blending techniques and machine learning and/or artificial intelligence methods, can be used to improve precipitation evolution nowcasting systems (Yu et al., 2012; Zheng et al., 2015; Yu and Zheng, 2020). The optical flow technique has potential applications in both precipitation and wind nowcasting (Bechini and Chandrasekar, 2017).

References

[1] Bechini R, Chandrasekar V (2017). An enhanced optical flow technique for radar nowcasting of precipitation and winds. J Atmos Ocean Technol, 34(12): 2637–2658

[2] Bouguet J Y (2000). Pyramidal Implementation of the Lucas Kanade Feature Tracker description of the algorithm. Intel Corporation: 1–9

[3] Bowler N E H, Pierce C E, Seed A (2004). Development of a precipitation nowcasting algorithm based on optical flow techniques. J Hydrol (Amst), 288(1–2): 74–91

[4] Browning K A, Collier C G (1989). Nowcasting of precipitating systems. Rev Geophys, 27(3): 345–370

[5] Cao C Y, Chen Y Z, Liu D H, Li C, Li H, He J J (2015). The optical flow method and its application to nowcasting. Acta Meteorol Sin, 73(3): 471–480 (in Chinese)

[6] Chen L, Dai J H, Tao L (2009). Application of an improved TREC algorithm (CoTREC) for precipitation nowcast. J Trop Meteorol, 25(1): 117–122 (in Chinese)

[7] Chen L S, Li Y (2004). An overview on the study of the tropical cyclone rainfall. In: Proc. Inter. Conf. on Storms, Brisbane, Australian Meteorological and Oceanographic Society: 112–113

[8] Chen L S, Li Y, Cheng Z Q (2010). An overview of research and forecasting on rainfall associated with landfalling tropical cyclones. Adv Atmos Sci, 27(5): 967–976

[9] Chen P Y, Yu H, Xu M, Lei X T, Zeng F (2019). A simplified index to assess the combined impact of tropical cyclone precipitation and wind on China. Front Earth Sci, 13(4): 672–681

[10] Chen Y Z, Lan H P, Chen X L, Zhang W H (2017). A nowcasting technique based on application of the particle filter blending algorithm. J Meteorol Res, 31(5): 931–945

[11] Crane R K (1979). Automatic cell detection and tracking. IEEE Trans Geosci Electron, 17(4): 250–262

[12] Crum T, Ciardi E J, Boettcher J B, Istok M, Stern A (2013). How the wsr-88d and its new dual polarization capability can benefit the wind energy industry. In: 93rd American Meteorological Society Annual Meeting

[13] Dai J H (2013). Analysis of thunderstorm development and evolution features and mechanism study for the Yangtze River Delta region. Dissertation for the Doctoral Degree. Nanjing: Nanjing University (in Chinese)

[14] Dixon M, Wiener G (1993). TITAN: Thunderstorm identification, tracking, analysis, and nowcasting-A radar-based methodoloy. J Atmos Ocean Technol, 10(6): 785–797

[15] Dusan S Z, Valery M M, Alexander V R (2006). Correlation coefficients between horizontally and vertically polarized returns from ground clutter. J Atmos Ocean Technol, 23(3): 381–394

[16] Emanuel K, Center L (2018). 100 Years of progress in tropical cyclone research. Meteorol Monographs, 59: 15.1–15.68

[17] Fabry F (2015). Radar Meteorology: Principles and Practice. Cambridge: Cambridge University Press

[18] Germann U, Zawadzki I (2002). Scale-dependence of the predictability of precipitation from continental radar images. Part I: description of the methodology. Mon Weather Rev, 130(12): 2859–2873

[19] Gibson J J (1950). The Ecological Approach to Visual Perception. Boston: Houghton Mifflin

[20] Han L, Wang H Q, Lin Y J (2008). Application of optical flow method to nowcasting convective weather. Acta Scientiarum Naturalium Universitatis Pekinensis, 44(5): 751–755

[21] Houze R Jr (2010). Clouds in tropical cyclones. Mon Weather Rev, 138(2): 293–344

[22] Johnson J T, MacKeen P L, Witt A, Mitchell E D, Stumpf G J, Eilts M D, Thomas K W (1998). The storm cell identification and tracking (SCIT) algorithm: an enhanced WSR-88D algorithm. Weather Forecast, 13(2): 263–276

[23] Lucas B D, Kanade T (1981). An iterative image registration technique with an application to stereo vision. Proceedings of Imaging Understanding Workshop, 121–130

[24] Powell S W, Houze R A, Brodzik S R (2016). Rainfall-type categorization of radar echoes using polar coordinate reflectivity data. J Atmos Ocean Technol, 33(3): 523–538

[25] Rinehart R E, Garvey E T (1978). Three-dimensional storm motion detection by conventional weather radar. Nature, 273(5660): 287–289

[26] Schaefer J T (1990). The critical success index as an indicator of warning skill. Weather Forecast, 5(4): 570–575

[27] Shi J B, Tomasi C (1994). Good features to track. In: Proceedings IEEE Computer Society Conference on Computer Vision and Pattern Recognition (CVPR) Seattle: 593–600

[28] Steiner M, Houze R Jr, Yuter S E (1995). Climatological characterization of three-dimensional storm structure from operational radar and rain gauge data. J Appl Meteorol, 34(9): 1978–2007

[29] Sun M, Dai J H, Yuan Z H, Tao L (2015). An analysis of a back-propagating thunderstorm using the three-dimensional wind fields retrieved by the dual-Doppler radar data. Acta Meteorol Sin, 73(2): 247–262 (in Chinese)

[30] Tuttle J D, Foote G B (1990). Determination of the boundary layer airflow from a single Doppler radar. J Atmos Ocean Technol, 7(2): 218–232

[31] Wexler H (1945). The structure of the September, 1944, Hurricane when off Cape Henry, Virginia. Bull Am Meteorol Soc, 26(5): 156–159

[32] Wexler H (1947). Structure of hurricanes as determined by radar. Ann N Y Acad Sci, 48(8): 821–844

[33] Wilson J (2004). Precipitation nowcasting: past, present and future. In: Sixth Int. Symp. on Hydrological Applications of Weather Radar, Melbourne, VIC, Australia, Centre for Australian Weather and Climate Research, 8

[34] Woo W C, Wong W K (2017). Operational application of optical flow techniques to radar-based rainfall nowcasting. Atmosphere, 8(12): 48

[35] Yang L H, Liu X Y, Chen Y H (2017). Product analysis of correlation coefficient of dual polarization radar. Guangdong Meteorol, 39(3): 69–72 (in Chinese)

[36] Yu H, Chen L S (2019). Impact assessment of landfalling tropical cyclones: introduction to the special issue. Front Earth Sci, 13(4): 669–671

[37] Yu X D, Zheng Y G (2020). Advances in severe convective weather research and operational service in China. J Meteorol Res, 34(2): 189–217

[38] Yu X D, Zhou X G, Wang X M (2012). The advances in the nowcasting techniques on thunderstorms and severe convection. Acta Meteorol Sin, 70: 311–337 (in Chinese)

[39] Zhang L, Wei M, Li N, Zhou S H (2014). Improved optical flow method application to extrapolate radar echo. Sci Technol Engineering, 14(32): 133–137 (in Chinese)

[40] Zhang Q H, Wei Q, Chen L S (2010). Impact of landfalling tropical cyclones in Chinese mainland. Sci China Earth Sci, 53(10): 1559–1564

[41] Zheng Y G, Zhou K H, Sheng J, Lin Y J, Tian F Y, Tang W Y, Lan Y, Zhu W J (2015). Advances in techniques of monitoring, forecasting and warning of severe convective weather. J Applied Meteorol Sci, 26(6): 641–657 (in Chinese)

[42] Zhou G B, Gao S Z (2019). Analysis of the August 2019 atmospheric circulation and weather. Meteorol Mont, 45(11): 1621–1628 (in Chinese)

RIGHTS & PERMISSIONS

Higher Education Press