The decommissioning of the Fukushima Daiichi (1F) nuclear power plant after the tsunami in 2011 is in progress  . Interior inspections of the primary containment vessels of units 1-3 have been conducted recently to assist in the removal of fuel debris. These inspections must determine the distribution of fuel debris and locate leakage points in the superstructure. Leakages must be located as water is being pumped into the vessels to cool the fuel debris, and it has been leaking from the vessels after having been contaminated by the radioactive debris. Optical inspections were conducted with a video camera in July 2017, and the video cameras were able to capture the first images of the fuel debris  . However, the spatial distribution of the fuel debris could not be accurately determined from the camera footage. Leakages were also not located, since visibility was limited in the murky water. Primary containment vessels will likely need to be inspected with a non-optical technique.
Ultrasonic measurements were considered for this purpose, given that ultrasound can be used in opaque liquids and ultrasonic transducers are generally robust against high radiation levels, as demonstrated during the decommissioning of the Three Mile Island (TMI-2) power plant  or during the treatment of radioactive wastewater  . Among the many ultrasonic measurement techniques available today, the present work focuses on echo particle-image velocimetry (echo-PIV)  . In this method, echo-intensity images are recorded with a sensor array at specific interval. Then, those images are processed with a PIV algorithm to construct velocity vector fields. Echo-PIV can be used to identify leakage points by visualising the flow of the fluid inside the containment vessel. In addition, the surface of stationary objects can be imaged by superposing the echo intensities in the obtained images in parallel to PIV processing. These images can be used to determine the distribution of debris within an underwater vessel.
The conventional echo-PIV method presented in  uses a convex sensor array and employs sector-scan imaging as shown in Figure 1(a). The scanning time for one image is determined by the range to the object of interest and the number of transmitters in the array:
where Ts is the scanning time, Tml is the time a transmitted pulse needs to reach are flecting object (Tml = l/c, l: measurement length, c: velocity of sound in the medium) and Nt is the number of transmitters. If a range on the order of metresis required, such as when inspecting a containment vessel underwater, the scanning time becomes prohibitively so long such that PIV processing becomes impractical. One way to reduce the scanning time required to generate an image without decreasing the field of view is developed in reference  . This method uses a linear sensor array instead of a convex sensor array. This array can generate steered acoustic beams by applying transmission delays to individually addressable elements (Figure 1(b)). The diverging beams are used to decrease the angular separation between adjacent measurement lines, increasing the field of view with fewer transmissions. This allows for faster scanning than with the sector scanning introduced previously. However, this advanced linear-array
Figure 1. Illustration of imaging methods: (a) Classical sector scanning with convex array, (b) Faster scanning with beam steered array and (c) Our proposed method.
method is still not well suited to inspections of submerged containment vessels due to the scanning time required. In addition to the potential imaging time required being impractically long, the sector imaging requires spatial coherence between acquired data-sets for one echo image. In order to maintain the spatial coherence, the time between acquisitions needs to be suitably short such that any physical drift of the array does not impact sector imaging quality but is limited by 2Tml (see Equation (1)).
To address this problem, this paper proposes a novel PIV imaging method. The proposed imaging method uses a divergent transmitted beam and the delay-and-sum (DAS) algorithm for image processing  . Our proposed method uses a linear sensor array but emits only a single diverging wave (DW) that covers a wide field of view, as shown in Figure 1(c). Following the transmission, echoes from waves scattered by suspended particles or stationary objects are received by all elements of the array at once. The received signals are then used to reconstruct the shape and the location of the scattering objects using the DAS algorithm. We refer to this imaging method as DW-DAS scanning, and expected that it would capture long-range and wide-view echo images much faster than the methods that use sector scanning.
This paper reports initial validation tests of the DW-DAS echo-PIV method for the long-range and wide-view measurements needed in the inspection of the containment vessel at FukushimaDaiichi.
2. Design Considerations
This section describes the mathematical details of the PIV imaging method and the prototype system of the DW-DAS echo-PIV method developed.
2.1. Formulisation of DW-DAS
Figure 2 illustrates the proposed imaging method. A virtual point source is assumed behind the sensor array. Assuming that sound waves propagate from the virtual point source in concentric waves, the travel time at which the wave reaches each element is calculated from the distance and the speed of sound in the medium. Practically, the differences in the travel times from the shortest one, Δτmin is considered as the transmission delay. Hence, the transmission delays, Δτ are calculated via the following equations:
Figure 2. Illustration of the proposed DW-DAS imaging method.
where x and y are the spatial coordinates, subscripts v and e distinguish the virtual source and the element respectively and c is the speed of sound in the medium.
The delayed transmitted wave described above should propagate in a diverging manner and should cover a wide field of view with a single transmission. After receiving scatter echoes with all sensor elements at the same time, the echo intensity at each measurement point (pixel) in the region of interest (ROI) is calculated by summing the amplitudes of all received signals with the corresponding travelling time from the virtual source to the elements via the measurement point. This process describes the delay-and-sum (DAS) imaging algorithm. Finally, considering the decrease in amplitude over distance due to wave dispersion, the DW-DAS is formulated as follows:
where I is the echo intensity, N is the number of elements, H is the Hilbert transform of the received echo signal, and d is the distance from the position (x, y).
2.2. Prototype Equipment
A laboratory-scale measurement system was developed for initial testing and optimisation of the system. The hardware of the developed system includes a sensor array (2K0.65 × 7.5I-8CH, Japan Probe Co. Ltd.), pulser/receiver (JPR-10C-5CH, Japan Probe Co. Ltd.), digitizer (NI PXI-5105, National Instruments) and a control computer. The specifications of the sensor array are as follows: the centre frequency is 2 MHz, which has the low attenuation required for long measurement lengths, the element width a = 0.65 mm, the element height b = 7.5 mm, and the pitch p = 0.7 mm (Figure 3(a)). Figure 3(b) plots the frequency response of the sensor array, which was measured experimentally. The peak response used is at 2.4 MHz, slightly higher than the designed centre frequency of 2 MHz. Hence, 2.4 MHz was employed as the transmission frequency to maximise the signal-to-noise ratio (SNR).
Figure 4 shows a flow chart for the image processing in the developed system. First, a DW transmission is repeated over the number of pulse repetitions (Nrep) at the pulse repetition frequency (PRF) and all received signals are recorded. Then, the received signals are converted to echo images, In ( ). The echo images are paired into combinations of #1-#2, #2-#3, #3-#4, etc. and those combinations are processed with a PIV algorithm to estimate the velocity vector field. Cross-correlation fast Fourier transforms (FFT) using a small window size with offsetting are used in the PIV algorithm  . Estimated velocity vectors, Vn are averaged over the number of image pairs (Nrep − 1) and a final vector flow map, is obtained. In addition, the averaged echo image, is also recorded for measuring the shape of surface reflectors (e.g. debris).
3. Performance Testing
The performance of the laboratory-scale prototype was evaluated and optimised. This section describes tests that determined the ideal position for a virtual point source, investigated the imaging resolution and optimised the density of tracer particles.
3.1. Locating a Virtual Point Source
As assumed in Section 2.1, the beam shape of the transmitted DW changes with the location of the virtual point source. This relationship was studied experimentally by measuring the sound field for a range of virtual point source positions to find the optimal virtual point source position that maximises the field of view.
Figure 5 shows schematics of the apparatus we used to measure the sound pressure distribution. This equipment consists of a water tank (1000 × 750 × 500 mm) with an automatic xyz-stage. The water tank was filled with tap water and the sensor array was deployed at x = 0 mm, y = 0 mm and z = 0 mm facing a hydrophone (5 K 0.5 × 0.5, Japan Probe Co. Ltd.). The sound pressure distribution of the wave transmitted from the sensor array was measured by moving the hydrophone within the tank (with the digitizer (NI PXI-5114, National Instruments)). The sound pressure at each point was defined as the peak-to-peak amplitude in the recorded signals. The measurement conditions were as follows: the water temperature was 30˚C, the speed of sound was 1509 m/s, 3 wave cycles were in each pulse, one measurement plane was in the xy-plane (0 ≤ x ≤ 150 mm, −40 ≤ y ≤ 40 mm, z = 0 mm) and three were parallel to the yz-plane (x = 40, 80, 120 mm, −40 ≤ y ≤ 40 mm, −20 ≤ z ≤ 20 mm) with a step size of 1.0 mm. Finally, the data from 128 pulse repetitions were averaged for each measurement point. The virtual point source was placed on the x axis (y = 0 mm, z = 0 mm) and was moved from x = −0.2 to x = −4.0 mm with a 0.2 mm pitch. These measurements were conducted for each position of the virtual point source.
Figure 6 shows the results for the virtual point source position of x = −1.8 mm. This virtual point source position maximises the beam width. The sound pressure distributions were normalised to the highest detected pressure.
Figure 6(a) indicates that the main lobe of the transmitted DW is within a range of ±20˚, so that angular area was assumed as the field of view. Figure 6(b) shows that the thickness of the transmitted divergent wave tends to increase with distance from the sensor array. This indicates that the resolution of 2-D echo imaging along the z axis will reduce with increasing x, so that the resolution may cause low detection accuracy on the xy plane. Thus, the optimal virtual point source position is x = −1.8 mm and the field of view was within ±20˚ at the measured length.
Figure 3. Sensor array used in this study. (a) Key array dimensions; (b) Frequency response.
Figure 4. Flow chart of image processing.
Figure 5. Schematic measurement equipment for sound pressure distribution.
Figure 6. The sound pressure distribution of the sensor array (step size: 1 mm). (a) xy-plane; (b) yz-plane.
3.2. Resolution Tests
Imaging resolution affects the system’s ability to detect the roughness of astationary object surface and to distinguish tracer particles, so we tested there solution of our prototype system.
The equipment described in the previous section was used for these tests (see Figure 5), except for the hydrophone, which was replaced with a tensevertical stainless wire (Ø287 μm) attached to the tip of the xyz-stage as a point-like scattering object. Imaging was conducted while moving the wire with the xyz stage, and the imaging resolution was investigated at each wire position we tested. Figure 7 illustrates our definition of imaging resolution using an example of the wire image. The obtained echo image was normalised to its maximum value. The image border was evaluated by following the −3 dB level. The maximum lengths of the image along the x and y axes were defined as the x and y resolutions, respectively. The experimental conditions are shown in Table 1.
The measurement results are shown in Figure 8(a) and Figure 8(b), which show the x and y resolutions. The resolutions are expressed as multiples of the centre-frequency wavelength, λ. The x resolution depends only on the distance from the y axis of the sensor array. In contrast, the y resolution decreased significantly as the distance from the sensor array increased. To conclude, the maximum x and y resolutions were 8λ and 25λ, respectively.
Figure 7. Definitions of x and y resolution when measuring a wire at x = 70 mm, y = 20 mm.
Figure 8. Imaging resolution expressed as multiples of the wavelength. (a) x resolution; (b) y resolution.
Table 1. Measurement condition and imaging settings.
HPF and LPF are high-pass filter and low-pass filter respectively.
3.3. Tracer Particle Density
Tracer particles are used in PIV. The tracer particle density must be optimised with consideration of the imaging resolution of the developed system to minimise miscorrelations during PIV processing. To determine the proper tracerparticle density, we conducted a flow measurement test.
Figure 9 shows schematics of the experimental apparatus. This apparatus consists of a magnetic stirrer (HS-30D, AS ONE Co. Ltd.), a stirrer bar, and a water tank (145 × 145 × 145 mm). The water tank was filled with tap water and tracer particles (Polyethylene, Ø300 μm, Kaneka Co. Ltd.) were mixed in. The stirrer bar was positioned at the centre of the water tank and was spun at 500 rpm by the magnetic stirrer. Measurements were recorded for 15 seconds after stopping the magnetic stirrer, using the prototype echo-PIV system and an optical-PIV system for the sake of comparison. These measurements were repeated 16 times to determine the average velocity field in the water tank. The tracer density was defined by the spatial volume of an individual tracer; specifically, the diameter of a sphere with the same volume of water. These density diameters were set to: 3λ, 6λ, 9λ, 12λ and 15λ. The other conditions and settings were as shown in Table 2. The window sizes in the echo-PIV processing were set as 128 × 128, 96 × 96 and 64 × 64 pixels with 50% overlap.
Figure 10 shows all measurement points, which were extracted within the ±20˚ range. The red and blue arrows show the time-averaged vector-flow maps as estimated by echo and optical PIV, respectively. At densities of 3λ, 6λ and 9λ,
Table 2. Measurement condition and imaging settings.
HPF and LPF are high-pass filter and low-pass filter respectively.
Figure 9. Experimental apparatus.
Figure 10. Results from echo and optical PIV (The grey-scale displays the first recorded echo image normalised to its maximum value). (a) Index of meas. point; (b) 3λ; (c) 6λ; (d) 9λ; (e) 12λ; (f) 15λ.
the echo vectors in the region x ≥ 100 mm agree better with the optical vectors as tracer density decreases. This trend is explained by how the transmitted DW could not reach and scatter from tracers far from the array, and how the multitude of echoes from dense tracers degrades the signal-to-noise ratio (SNR). Therefore, the tracer density should be at least 9λ.
The measurement results for 9λ, 12λ and 15λ densities are compared in Figure 11 using the measurement error of the velocity components v and u, which are along the x and y axes, respectively. The measurement error was calculated in comparison with the optical-PIV results. The horizontal axis is the index number of the measurement point as marked in Figure 10(a). As indicated in Figure 11, the three particle densities returned similar errors, and the average error with density of 9λ was the lowest. The error in the v component did not change considerably. In contrast, the u-component error dropped significantly for measurement points farther from the array by about 100 mm or more on the x axis. Tracer particles have more time to move between two subsequent echo images as the y resolution increases (see Figure 8(b)).
As a way to reduce the error caused by this small movement, we considered varying the PRF with the tracer density of 9λ. The images were paired to adjust the PRF to 25 or 50 Hz, and the pairs were then PIV-processed. Figure 12 shows that the error increased as the PRF decreased, likely due to miscorrelations introduced during PIV processing because of changes in the tracer shape between two paired images. The dynamic range of the velocity vector field is thus limited by the imaging resolution. To improve the imaging resolution, the sensor array can be made larger and the number of elements in the sensor array could be increased.
To conclude, the optimal density of the tracer particles is 9λ for the developed system. The velocity vector field can be estimated with an error less than 30% within 100 mm from the sensor array.
Figure 11. Velocity measurement errors using echo-PIV for three densities of tracer particles. (a) v component along the x axis. (b) u component along the y axis.
Figure 12. Measurement error of velocity using echo-PIV for different PRFs. (a) v component along the x axis. (b) u component along the y axis.
4. Feasibility Tests
4.1. Mock Containment Vessel Tests
Finally, to validate the feasibility of echo-PIV using DW-DAS for long-range and wide-view measurements, a mock up of a containment vessel containing simulated debris was created in the lab.
Figure 13(a) shows a schematic of the experimental apparatus. This apparatus consists of a water tank (290 × 290 × 390 mm), 1-D stage and a flow meter. The tank was filled with tap water and tracers were mixed in at the density of 9λ. A rough surfaced stone was placed at the bottom of the water tank to simulate debris. Figure 13(b) shows the details of the situation at the bottom of the water tank, as mapped by a 3D laser scanner (OPT MX, OPT Co. Ltd.), for the sake of comparison. The water tank had an outflow outlet (Ø24 mm) located at x = 100 mm, y = 0 mm and z = 0 mm. The sensor array was initially located at x = 0 mm, y = 0 mm and z = −40 mm. Measurements were recorded while moving the sensor array from −40 to 110 mm at 5 mm steps along the z axis while draining the stored water from the outlet. The measurement conditions and imaging settings are shown in Table 3. The window sizes in the echo PIV algorithm were set as 128 × 128, 96 × 96, 64 × 64 and 32 × 32 pixels with a 50% overlap.
4.2. Results & Discussion
The height corresponding to the maximum value of the averaged echo images at each measurement point was extracted and is displayed in grey-scale. The vectors displayed are only over 10 mm/s to facilitate visualizing the 3D flow behaviour.
High magnitude vectors were located around the y = 0 mm and z = 0 m and Figure 14(b) shows the velocity vector map obtained at z = 0 mm. The velocity vectors were convergent between −15 mm and 15 mm on the y axis, which indicates the location of a leak, closely matching the actual drain location.
Figure 13. Schematic of experimental apparatus for mapping flow and imaging simulated fuel debris (the outlet shows the location of the water flowing from the tank. The stone has a similar shape to the expected fuel debris). (a) Overview; (b) Detail of the water tank floor.
Figure 14. Experimental results of the combines flow mapping and object shape reconstruction. (a) All results (vector flow map & detected object shape). (b) Vector flow map (z = 0 mm).
Table 3. Measurement condition and imaging settings.
HPF and LPF are high-pass filter and low-pass filter respectively.
The absolute error in the map of the bottom surface of the water tank, with the 3D optical scan being the reference, was consistently low and the average was approximately 1.5 mm. On the other hand, errors at the surface of the stone were larger, with an average of around 9.7 mm. These errors are likely larger because the y resolution was less than that of the stone’s surface roughness, so a significant echo signal could not be received. The shape of the mock fuel debris was able to be determined to a reasonably high degree.
Consequently, we found that the prototype system could identify the leakage point and determine the shape of the mock fuel debris in the angular width within ±20˚ and at 100 mm away from the sensor array. The measurement range is expected to be increased with a larger sensor array to improve the imaging resolution as mentioned in Section 3.3, and therefore the feasibility of echo-PIV using DW-DAS for long-range and wide-view measurements was indicated. Furthermore, real-time measurements can be realized with the proposed method by speeding up the imaging  and PIV processing  steps for practical use.
We designed and tested an echo-PIV with DW-DAS system capable of visualising fluid flow and reconstructing simulated fuel debris shape in water. The system has direct applicability to the ongoing decommissioning of Fukushima Daiichi nuclear power plant.
Three initial optimisation tests allowed us to determine key parameters of the proposed imaging system. First, we optimised the location of the virtual point source in the PIV algorithm to be 1.8 mm behind the sensor array. This yields a field of view of ±20˚. We then demonstrated that the imaging resolution in the transverse direction depends greatly on the distance between the object of interest and the sensor array. Finally, we determined the appropriate tracer density that minimizes the measurement errors of velocity vector filed.
Having determined the optimal system parameters with these initial tests, we deployed the prototype system in a mock containment vessel. The prototype small-scale imaging system was located in a mock leak successfully and returned a suitably accurate reconstruction of mock fuel debris. These results suggest that the DW-DAS echo-PIV method is applicable for sensing leakage points and underwater debris at relatively long range and over a wide field of view.
a: Element width [m]
b: Element height [m]
c: Velocity of sound in a medium [m/s]
d: Distance from an element to a measurement position [m]
e: Index of elements
I: Echo intensity [-]
H: Hilbert transform of received echo signals
L: Measurement length [m]
n: index of obtained echo-intensity images
N: Number of elements
Nt: Number of transmitters
Nrep: Number of pulse repetitions
p: Element pitch [m]
Ts: Scanning time for one echo image [s]
Tml: Time a transmitted pulse needs to reach a reflecting object [s]
v: index of virtual point sources
V: Vector flow map [m/s]
λ: Centre-frequency wavelength [m]
Δτmin: The shortest travel time from a virtual point source to elements [s]
Δτ: Transmission delay [s]
 Tokyo Electric Power Company Holdings, Inc. (2018) Summary of Decommissioning and Contaminated Water Management.
 Hou, S.X., Luo, J.J., He, B., Li, R.S. and Shen, T. (2014) The Treatment of Radioactive Wastewater by Ultrasonic Standing Wave Method. Journal of Hazardous Materials, 274, 41-45.
 Hasegawa, H. and Kanai, H. (2011) High-Frame-Rate Echocardiography Using Diverging Transmit Beams and Parallel Receive Beamforming. Journal of Medical Ultrasonics, 38, 129-140.
 Langenberg, K.J., Berger, M., Kreutter, Th., Mayer, K. and Schmitz, V. (1986) Synthetic Aperture Focusing Technique Signal Processing. NDT International, 19, 177-189.
 Rougeron, G., Lambert, J., Iakovleva, E., Lacassagne, L. and Dominguez, N. (2014) Implementation of a GPU Accelerated Total Focusing Reconstruction Method within Civa Software. AIP Conference Proceedings, 1581, 1983-1990.
 Champagnat, F., Plyer, A., Le Besnerais, G., Leclaire, B., Davoust, S. and Le Sant, Y. (2011) Fast and Accurate PIV Computation Using Highly Parallel Iterative Correlation Maximization. Experiments in Fluids, 50, 1169-1182.