Method for measuring microphysical characteristics of natural precipitation using particle image velocimetry
A method and video sensor for precipitation microphysical features measurement based on particle image velocimetry. The CCD camera is placed facing towards the light source, which forms a three-dimensional sampling space. As the precipitation particles fall through the sampling space, double-exposure images of precipitation particles illuminated by pulse light source are recorded by CCD camera. Combined with the telecentric imaging system, the time between the two exposures are adaptive and can be adjusted according to the velocity of precipitation particles. The size and shape can be obtained by the images of particles; the fall velocity can be calculated by particle displacement in the double-exposure image and interval time; the drop size distribution and velocity distribution, precipitation intensity, and accumulated precipitation amount can be calculated by time integration. This invention provides a method for measuring the shape, size, velocity, and other microphysical characteristics of various precipitation particles.
1. A method for measuring microphysical characteristics of precipitation particles using particle image velocimetry comprising following steps:
a. a CCD camera is placed facing towards the light source, which forms a three-dimensional sampling space;
b. auto-correlation images of precipitation particles illuminated by pulse light source are recorded by CCD camera;
c. horizontal and vertical dimensions of particles are inferred from the images based on the scaling relation between the images and the size;
d. shape, axis ratio, as well as orientation of particles are calculated by the image processing algorithm;
e. vertical and horizontal velocities of particles are calculated according to the particle displacement in double-exposure image and interval time;
f. phases of precipitation particles are classified according to the shape of particles, and the relationship between the velocity and size of particles;
wherein the defocused and noisy images are restored using the median filtering and point spread function, the shape, size, axis ratio, orientation of particles are calculated based on grayscale and contour of images, and the fall velocity of particles are calculated from auto-correlation images using the particle grayscale tracking method.
2. The method according to claim 1 , wherein the CCD camera takes images of precipitation particles with a frame frequency of 50±25 frame per second.
3. The method according to claim 1 , wherein the CCD camera and pulse light source are synchronized by a microprocessor and take a double-exposure of each precipitation particle on one frame, where each exposure time is no more than 10 μs.
4. The method according to claim 3 , wherein the time between two exposures are adaptive and can be adjusted according to the velocities of different particles with scales of 1 ms˜3 ms, so that the two exposures of the same particle can be recorded in one image.
5. The method according to claim 1 , wherein the spectral distribution of shape, velocity, size, axis ratio, orientation, and phase of precipitation particles are obtained by time integration.
6. The method according to claim 1 , wherein the background of imaging system is detected real-time using the three-dimensional sampling space, and the background image is subtracted from the raw image to obtain a difference image with flat background and clear particles.
7. The method according to claim 1 , wherein non-precipitation images are identified and eliminated by time series analysis of sequential frames using the following criteria: when images of the same size are detected in the same position in several sequential frame images, they are eliminated as noises, water spots on lens or other non-precipitation images.
8. The method according to claim 1 , wherein the raw image of particles are processed in the following steps:
a. Filter denoising: the noises in raw images are eliminated by median filtering;
b. Image restoration: the raw images with fuzzy edge are restored using the point spread function and the image degradation model is:
g ( x,y )= f ( x,y )* h ( x,y )+ n ( x,y ) (1)
where g(x, y) is the blur image, f(x, y) is the ideal image, n(x, y) is the additive noise, h(i, j) is the point spread function, which is defined as:
h
(
i
,
j
)
=
{
1
/
(
π
R
2
)
i
2
+
j
2
≤
R
0
i
2
+
j
2
>
R
(
2
)
where the R is the defocusing radius, x, y, i, j are the coordinates of pixels in images;
c. Image binarization and edge detection: the binarization images of particles are obtained by a global threshold value based on the grayscale histogram, and the edges of particles are detected using the binarization threshold.
9. The method according to claim 1 , wherein the microphysical characteristics are extracted using the following algorithms:
a. Area: the summation of pixels with grayscale values that are larger than the threshold:
A
R
=
∑
(
x
,
y
)
∈
R
1
(
3
)
b. Perimeter: the summation of pixel values on the contour of particle's image:
P
R
=
∑
(
x
,
y
)
∈
C
1
(
4
)
c. Long axis: the long axis of minimum enclosing rectangle of particle's image:
E
A
=
max
p
i
,
p
j
∈
C
⌊
D
(
p
i
,
p
j
)
⌋
(
5
)
d. Short axis: the short axis of minimum enclosing rectangle of particle's image:
E
B
=
max
p
i
,
p
j
∈
C
,
p
i
,
p
j
⊥
E
A
⌊
D
(
p
i
,
p
j
)
⌋
(
6
)
e. Equivolumetric diameter: the diameter of spherical particle with the same volume of non-spherical particle:
D
eq
=
(
6
V
π
)
-
1
/
3
=
(
E
A
2
E
B
)
-
1
/
3
(
7
)
f. Axis ratio: the ratio of short axis and the long axis of particle:
α
=
E
B
E
A
(
8
)
g. Canting angle: the angle between the symmetry axis and vertical direction:
θ=∠{right arrow over ( E B )}/{right arrow over (Vertical)} (9)
h. Velocities: the moving velocity of particles in horizontal direction and vertical direction:
{
v
x
=
Δ
x
Δ
t
=
x
2
-
x
1
Δ
t
v
y
=
Δ
y
Δ
t
=
y
2
-
y
1
Δ
t
(
10
)
i. The probability of a certain particle with a certain size and fall velocity being captured by the CCD camera:
Prob
(
D
)
=
H
-
D
V
-
V
(
D
)
·
T
interval
V
(
D
)
·
T
frame
(
11
)
where H is the height of sampling size, D V is the vertical size of drops, V (D) is the vertical velocity of drops, T interval is the exposure interval in each frame, T frame is the exposure time of each frame;
j. The raindrop size distribution N(D):
N
(
D
)
=
Num
(
D
)
Vol
·
dD
·
Prob
(
D
)
(
12
)
where N(D) is the sampled number concentration of particles per volume per size (m −3 mm −1 ); Num(D) is the number of particles of each bin, which can be obtained by counting particles and summing volume increments during the time interval; Vol is the effective sampling volume; dD is the interval size of particle.
10. The method according to claim 6 , wherein non-precipitation images are identified and eliminated by time series analysis of sequential frames using the following criteria: when images of the same size are detected in the same position in several sequential frame images, they are eliminated as noises, water spots on lens or other non-precipitation images.
11. The method according to claim 6 , wherein the raw image of particles are processed in the following steps:
a. Filter denoising: the noises in raw images are eliminated by median filtering;
b. Image restoration: the raw images with fuzzy edge are restored using the point spread function and the image degradation model is:
g ( x,y )= f ( x,y )* h ( x,y )+ n ( x,y ) (1)
where g(x, y) is the blur image, f(x, y) is the ideal image, n(x, y) is the additive noise, h(i, j) is the point spread function, which is defined as:
h
(
i
,
j
)
=
{
1
/
(
π
R
2
)
i
2
+
j
2
≤
R
0
i
2
+
j
2
>
R
(
2
)
where the R is the defocusing radius, x, y, i, j are the coordinates of pixels in images;
c. Image binarization and edge detection: the binarization images of particles are obtained by a global threshold value based on the grayscale histogram, and the edges of particles are detected using the binarization threshold.
12. The method according to claim 6 , wherein the microphysical characteristics are extracted using the following algorithms:
a. Area: the summation of pixels with grayscale values that are larger than the threshold:
A
R
=
∑
(
x
,
y
)
∈
R
1
(
3
)
b. Perimeter: the summation of pixel values on the contour of particle's image:
P
R
=
∑
(
x
,
y
)
∈
C
1
(
4
)
c. Long axis: the long axis of minimum enclosing rectangle of particle's image:
E
A
=
max
p
i
,
p
j
∈
C
⌊
D
(
p
i
,
p
j
)
⌋
(
5
)
d. Short axis: the short axis of minimum enclosing rectangle of particle's image:
E
B
=
max
p
i
,
p
j
∈
C
,
p
i
p
j
⊥
E
A
⌊
D
(
p
i
,
p
j
)
⌋
(
6
)
e. Equivolumetric diameter: the diameter of spherical particle with the same volume of non-spherical particle:
D
eq
=
(
6
V
π
)
-
1
/
3
=
(
E
A
2
E
B
)
-
1
/
3
(
7
)
f. Axis ratio: the ratio of short axis and the long axis of particle:
α
=
E
B
E
A
(
8
)
g. Canting angle: the angle between the symmetry axis and vertical direction:
θ=∠{right arrow over ( E B )}/{right arrow over (Vertical)} (9)
h. Velocities: the moving velocity of particles in horizontal direction and vertical direction:
{
v
x
=
Δ
x
Δ
t
=
x
2
-
x
1
Δ
t
v
y
=
Δ
y
Δ
t
=
y
2
-
y
1
Δ
t
(
10
)
i. The probability of a certain particle with a certain size and fall velocity being captured by the CCD camera:
Prob
(
D
)
=
H
-
D
V
-
V
(
D
)
·
T
interval
V
(
D
)
·
T
frame
(
11
)
where H is the height of sampling size, D V is the vertical size of drops, V(D) is the vertical velocity of drops, T interval is the exposure interval in each frame, T frame is the exposure time of each frame;
j. The raindrop size distribution N(D):
N
(
D
)
=
Num
(
D
)
Vol
·
dD
·
Prob
(
D
)
(
12
)
where N(D) is the sampled number concentration of particles per volume per size (m −3 mm −1 ); Num(D) is the number of particles of each bin, which can be obtained by counting particles and summing volume increments during the time interval; Vol is the effective sampling volume; dD is the interval size of particle.