Water Temperature Observation by Coastal Acoustic Tomography in Artificial Upwelling Area

Artificial upwelling is a geoengineering method to repair and improve marine ecosystems, and its operation requires long-term and continuous temperature field observation. However, existing methods are rarely seen to accomplish such observation. In this study, we investigate the coastal acoustic tomography (CAT) to obtain the long-term horizontal temperature field of an artificial upwelling area in an anechoic tank. We conduct four sets of experiments with different CAT station numbers and compare their data with those collected from temperature sensors. By analyzing the travel time from the CAT experiments, the horizontal temperature field of the upwelling area could be mapped. The CAT results and the comparison results show that the surface temperature of the observed area decreases by approximately 3 °C after upwelling, while the temperature of where the CAT is deployed decreases by about 1 °C; the temperature is lowest at the center of the upwelling area. Increasing the number of stations and station spacing would improve the temperature mapping accuracy. Therefore, the feasibility of using the CAT system to observe artificial upwelling is proved valid. This study indicates the potential application of CAT in temperature field observation in artificial upwelling area in the sea.


Introduction
Currently, more than half of the world's population lives within 60 km from a coast. The destruction of the marine ecosystem continues to worsen owing to increased activities, such as fisheries, mineral and energy extraction, coastal development, and marine waste emissions [1]. Consequently, the use of artificial methods to repair and improve marine ecosystems has become a crucial and necessary remedy.
Artificial upwelling is a geoengineering method that artificially lifts deep seawater to the surface using pump equipment in the operating sea area. Besides changing nutrient distribution and the nitrogen and phosphorus levels, it can enhance fish stocks [2] and improve regional carbon sink capacity through efficient absorption of CO 2 in the atmosphere [3]. It may also help mitigate local environmental extremes [4]. The optimal operation of an artificial upwelling system requires a supporting observation system; research to establish this is under way [5]. Koutaro et al. from Saga University and Tohoku University, Japan, used tracer detection systems for continuous measurement of upwelling [6]. This method could only observe the upwelling velocity of a fixed point and could not measure temperature changes at an arbitrary position of upwelling area. Xiamen University, Zhejiang University, China Ocean University, and several other universities built an artificial upwelling demonstration area in Ao Shan Bay, Qingdao, China in September 2017, with the goal of establishing a system that allows Subsequently, the temperature field is mapped using the inversion method. Note that the propagation of sound in a small area (no more than 100 m of station spacing in this study) can be approximated as linear propagation. Hence, the length of the acoustic propagation path can be assumed to be equal to the station spacing, and the length of the direct path need not be calculated by ray-tracing simulation.
Sensors 2019, 19, x FOR PEER REVIEW 3 of 16 in this study) can be approximated as linear propagation. Hence, the length of the acoustic propagation path can be assumed to be equal to the station spacing, and the length of the direct path need not be calculated by ray-tracing simulation.

Equations for Sound Speed Deviation Fields
When sound travels between two points in water, it is affected by the flow velocity and the temperature of the travel path, resulting in different travel times in different directions. As shown in

Equations for Sound Speed Deviation Fields
When sound travels between two points in water, it is affected by the flow velocity and the temperature of the travel path, resulting in different travel times in different directions. As shown in Figure 2, points A and B are two observation sites, while u is the flow component along the direction of the straight line (extending from A to B). The travel time can be obtained as follows: in this study) can be approximated as linear propagation. Hence, the length of the acoustic propagation path can be assumed to be equal to the station spacing, and the length of the direct path need not be calculated by ray-tracing simulation.

Equations for Sound Speed Deviation Fields
When sound travels between two points in water, it is affected by the flow velocity and the temperature of the travel path, resulting in different travel times in different directions. As shown in Figure 2, points A and B are two observation sites, while u is the flow component along the direction of the straight line (extending from A to B). The travel time can be obtained as follows: The equation depends on the underwater reference sound speed (C 0 ), the flow velocity (v), the unit vector along the path direction (n), and the change of sound velocity caused by temperature variation ∆C. With the constant temperature and the flow velocity of 0, we can simplify the formula as follows: When ∆C C 0 , u = v * n C 0 , the Taylor expansion of Equation (1) can be expressed as follows: Ignoring the higher-order terms, i.e., exponents of 2 or higher, we obtain the following: Combining Equations (2) and (4), we obtain the following: Since this study requiresthe temperature calculation, it is necessary to eliminate parameters that are irrelevantto the temperature. The resulting equation is as follows: The Fourier expansion formula of ∆C can be expressed as [21]: where: D = D j = a, A 00 , B 00 , A 01 , B 01 , · · · · · · , A NxNy Q(x, y) = Q j = 1, 1, 0, cos 2πy L y , sin 2πy L y , · · · · · · , cos 2π N x x L x + N y y L y , sin 2π N x x L x + N y y L y where x, y represents the two-dimensional position coordinates, and L x , L y represents the length and width of the compute domain, respectively. Substituting Equation (7) into Equation (6), we obtain the following: Equation (8) establishes the relationship between δτ and the sound velocity deviation ∆C. By solving this equation, the sound velocity in the calculated domain can be obtained. The corresponding temperature field can be obtained by substituting it into the empirical formula of sound velocity which was proposed by Medwin in 1975 [22]. This empirical formula is as follows: In this equation, C indicates the speed of sound, while T, S, and Z represent temperature, salinity, and depth, respectively. Owing to the small size of the experimental area in this study, the acoustic signal path is approximated as a straight line, and therefore the change of depth item can be ignored. Additionally, the salinity is considered to be zero.

Inversion Process
The matrix form of Equation (7) is as follows: y=Ex+n (10) In this equation, y = { δτ i } represents the time information obtained from the experiment, x = D j represents the unknown coefficient matrix to be calculated, n = n j represents the error, and represents the known coefficient matrix. The solution to Equation (10) involves an inversion process, and it is necessary to obtain the minimum error of x under the condition where a stochastic error of n is present. The conical least squares method can be used to solve this problem. The cost function is represented as: By setting the derivative of x to 0, we can obtain the value of variable x that minimizes the cost function: Errors can be expressed as:n The coefficient α can be calculated using the L-curve method [23]. This makes a trade-off between solution and error.

Experiments
In this section, the experiment setup is introduced first, including the parameters of CAT instruments, the station layouts, and the experiment environment. After that, the method of generating upwelling is elaborated upon and the experiment processes are detailed. Figure 3 is the experimental layout in the anechoic tank. The upwelling device is positioned in the center of the observed area. Four CAT instruments are placed at each side of the tank. Each is attached with a GPS receiver for time synchronization and a dual-use transducer with a frequency of 50 kHz. The CAT instrument and transducer were powered by 12 V and 24 V DC power supplies respectively. The transducer was placed at a water depth of 1 m. We believed that the station numbers and station distances possibly affected the CAT observation results. Right triangle, isosceles triangle, rectangle and lozenge are considered four typical shapes, therefore, four groups of experiments with different station layouts were designed and carried out, as displayed in Figure 4. The station distances were mainly determined by the dimensions of the experimental area in the anechoic tank.
The experiments were conducted in the anechoic tank at Zhejiang University, as shown in Figure 5. The anechoic tank, with a length of 50 m, width of 15 m, and depth of 10 m, was used to simulate a marine sound environment, and was free from disturbances to measure a variety of acoustic parameters. Its band coverage is between 1-100 kHz. A sliding experimental platform was installed above the anechoic tank. Thus, it allowed the release of the transducer, positioned the station accurately, and helped in acquiring the experimental data.   The experiments were conducted in the anechoic tank at Zhejiang University, as shown in Figure  5. The anechoic tank, with a length of 50 m, width of 15 m, and depth of 10 m, was used to simulate a marine sound environment, and was free from disturbances to measure a variety of acoustic parameters. Its band coverage is between 1-100 kHz. A sliding experimental platform was installed above the anechoic tank. Thus, it allowed the release of the transducer, positioned the station accurately, and helped in acquiring the experimental data.  The experiments were conducted in the anechoic tank at Zhejiang University, as shown in Figure  5. The anechoic tank, with a length of 50 m, width of 15 m, and depth of 10 m, was used to simulate a marine sound environment, and was free from disturbances to measure a variety of acoustic parameters. Its band coverage is between 1-100 kHz. A sliding experimental platform was installed above the anechoic tank. Thus, it allowed the release of the transducer, positioned the station accurately, and helped in acquiring the experimental data.   The experiments were conducted in the anechoic tank at Zhejiang University, as shown in Figure  5. The anechoic tank, with a length of 50 m, width of 15 m, and depth of 10 m, was used to simulate a marine sound environment, and was free from disturbances to measure a variety of acoustic parameters. Its band coverage is between 1-100 kHz. A sliding experimental platform was installed above the anechoic tank. Thus, it allowed the release of the transducer, positioned the station accurately, and helped in acquiring the experimental data. In this experiment, the pneumatic lifting method was used to generate the upwelling. The same method was used in the upwelling experiments conducted by Zhejiang University in QiandaoLake, China in October 2011 and December 2012 [24,25]. As displayed in Figure 6a, we designed a gas injection frame consisting of three concentric water pipe rings with a maximum diameter of 3 m. The pipe rings were made of irregular copolymer polypropylene. The gas injection frame had 90 holes that were evenly distributed with a 4.8 mm diameter. The two inlets of the gas injection frame were connected with two air compressors, to enable underwater gas injection. Before the experiment, a rectangular area of 4 × 5 m is prepared on the surface of the water, and the gas injection frame is placed 4 m underneath with the suspended lifting arm installed on the experimental platform. The upwelling is shown in Figure 6b. In this experiment, the pneumatic lifting method was used to generate the upwelling. The same method was used in the upwelling experiments conducted by Zhejiang University in QiandaoLake, China in October 2011 and December 2012 [24,25]. As displayed in Figure 6a, we designed a gas injection frame consisting of three concentric water pipe rings with a maximum diameter of 3 m. The pipe rings were made of irregular copolymer polypropylene. The gas injection frame had 90 holes that were evenly distributed with a 4.8 mm diameter. The two inlets of the gas injection frame were connected with two air compressors, to enable underwater gas injection. Before the experiment, a rectangular area of 4×5m is prepared on the surface of the water, and the gas injection frame is placed 4m underneath with the suspended lifting arm installed on the experimental platform. The upwelling is shown in Figure 6b. The experiment lasted for 28 days, from 17 June 2018 to 14 July 2018. Each time, the CAT instruments were first placed in the still water and worked for 4 h. Then, the air compressors were switched on and the upwelling was created and lasted for 4 h. The transducer could not receive and send a signal at the same time; therefore, each instrument was only allowed to send a signal once per cycle and the time interval between signal-sending was 1 min. For example, for the stations with the right triangle layout, K1, K2 and K3 took turns sending signals with 1 min intervals, and the whole cycle took 3 min. As the reflected signal was absorbed in the anechoic tank, any multipath effect was eliminated, resulting in a high quality of the received signals. The 7-order M sequence as the acoustic signal was used and this M sequence was repetitively sent 8 times every one minute. These 8 repeated signals were averaged in data processing to improve the signal-to-noise ratio (SNR). Each code element of the sequence utilized two carrier cycles. After collecting the received data from instruments, we could obtain the travel time between each pair of stations by cross-correlation. Then the method described in Section 2 could be used for data processing to obtain the temperature field.

Results and Discussions
In this section, we first perform the signal cross-correlation to acquire the travel time. Based on this, Section 4.2 maps the horizontal temperature field with the inversion method. In Section 4.3, we compare the CAT results with the data collected by temperature sensors that are employed at different water depths. Section 4.4 discusses the accuracy of the CAT observation results

The Signal Correlation
The signal correlation is performed to obtain travel time from the collected data in the experiment. Figure 7 displays the travel time for different layouts after cross-correlation. We take a three-dimensional plot of travel time for the right triangle layout as an example (shown in Figure 7a) to explain. The red circle is the point with the highest SNR for each set of acoustic signals. The SNR of the correlation peak is significantly higher than the other positions, and the travel time can be The experiment lasted for 28 days, from 17 June 2018 to 14 July 2018. Each time, the CAT instruments were first placed in the still water and worked for 4 h. Then, the air compressors were switched on and the upwelling was created and lasted for 4 h. The transducer could not receive and send a signal at the same time; therefore, each instrument was only allowed to send a signal once per cycle and the time interval between signal-sending was 1 min. For example, for the stations with the right triangle layout, K1, K2 and K3 took turns sending signals with 1 min intervals, and the whole cycle took 3 min. As the reflected signal was absorbed in the anechoic tank, any multipath effect was eliminated, resulting in a high quality of the received signals. The 7-order M sequence as the acoustic signal was used and this M sequence was repetitively sent 8 times every one minute. These 8 repeated signals were averaged in data processing to improve the signal-to-noise ratio (SNR). Each code element of the sequence utilized two carrier cycles. After collecting the received data from instruments, we could obtain the travel time between each pair of stations by cross-correlation. Then the method described in Section 2 could be used for data processing to obtain the temperature field.

Results and Discussions
In this section, we first perform the signal cross-correlation to acquire the travel time. Based on this, Section 4.2 maps the horizontal temperature field with the inversion method. In Section 4.3, we compare the CAT results with the data collected by temperature sensors that are employed at different water depths. Section 4.4 discusses the accuracy of the CAT observation results

The Signal Correlation
The signal correlation is performed to obtain travel time from the collected data in the experiment. Figure 7 displays the travel time for different layouts after cross-correlation. We take a three-dimensional plot of travel time for the right triangle layout as an example (shown in Figure 7a) to explain. The red circle is the point with the highest SNR for each set of acoustic signals. The SNR of the correlation peak is significantly higher than the other positions, and the travel time can be accurately identified. As shown in Figure 7b-e, the travel time for different layouts could be obtained after eliminating the wrong circles. In addition, the number of signals in Figure 7e is significantly lower, possibly due to a weak indoor GPS signal. This might lead to the CAT instruments becoming desynchronized and further the loss of data. Overall, the quality of the collected data was regarded satisfactory and sufficient for subsequent data processing. accurately identified. As shown in Figure 7b-e, the travel time for different layouts could be obtained after eliminating the wrong circles. In addition, the number of signals in Figure 7e is significantly lower, possibly due to a weak indoor GPS signal. This might lead to the CAT instruments becoming desynchronized and further the loss of data. Overall, the quality of the collected data was regarded satisfactory and sufficient for subsequent data processing.

Horizontal Temperature Field Mapping
After obtaining the travel time, we use the inversion method to map the temperature field. The resulting horizontal temperature field for each station layout is displayed in Figure 8. When the upwelling is in progress, its position is marked with a dashed circle in the figure. In each station layout, the results before starting the upwelling and during the upwelling are presented with 2-h intervals. For example, for the right triangle layout, the horizontal temperature field is mapped at 14:00, and 16:00 with no upwelling, and at 18:00 and 20:00 with upwelling occurring. In addition, to prevent the effect of regularities in the station layout, the calculated domain was rotated by 45°. Therefore, the coordinates in Figure 8 represent those after rotating the original width and length direction of the tank.

Horizontal Temperature Field Mapping
After obtaining the travel time, we use the inversion method to map the temperature field. The resulting horizontal temperature field for each station layout is displayed in Figure 8. When the upwelling is in progress, its position is marked with a dashed circle in the figure. In each station layout, the results before starting the upwelling and during the upwelling are presented with 2-h intervals. For example, for the right triangle layout, the horizontal temperature field is mapped at 14:00, and 16:00 with no upwelling, and at 18:00 and 20:00 with upwelling occurring. In addition, to prevent the effect of regularities in the station layout, the calculated domain was rotated by 45 • . Therefore, the coordinates in Figure 8 represent those after rotating the original width and length direction of the tank.

Horizontal Temperature Field Mapping
After obtaining the travel time, we use the inversion method to map the temperature field. The resulting horizontal temperature field for each station layout is displayed in Figure 8. When the upwelling is in progress, its position is marked with a dashed circle in the figure. In each station layout, the results before starting the upwelling and during the upwelling are presented with 2-h intervals. For example, for the right triangle layout, the horizontal temperature field is mapped at 14:00, and 16:00 with no upwelling, and at 18:00 and 20:00 with upwelling occurring. In addition, to prevent the effect of regularities in the station layout, the calculated domain was rotated by 45°. Therefore, the coordinates in Figure 8 represent those after rotating the original width and length direction of the tank. Based on the results of horizontal temperature field, we can see that: (1). The horizontal temperature distribution of a certain depth could be obtained and displayed. This depth is decided by the transducer placement. To obtain accurate temperature mapping, we suggest each transducer be placed at the same depth. (2). The temperature of the upwelling area is significantly lower than that of its surrounding regions; this is obvious in Figure 8a,b,d. This is because the cooler water at the depth of 4 m is brought to the surface by upwelling, reducing the water temperature at the surface of the upwelling area. The lowest water temperature is often at the upwelling center. (3). It is most obvious in Figure8d that the temperature in the upwelling area is much lower than its surrounding area. Compared with Figure 8a,b, the number of stations in Figure 8d is 4; compared with Figure 8c, the station spacing of Figure 8d is larger. We conclude that increasing Based on the results of horizontal temperature field, we can see that: (1) The horizontal temperature distribution of a certain depth could be obtained and displayed. This depth is decided by the transducer placement. To obtain accurate temperature mapping, we suggest each transducer be placed at the same depth. (2) The temperature of the upwelling area is significantly lower than that of its surrounding regions; this is obvious in Figure 8a,b,d. This is because the cooler water at the depth of 4 m is brought to the surface by upwelling, reducing the water temperature at the surface of the upwelling area. The lowest water temperature is often at the upwelling center. (3) It is most obvious in Figure 8d that the temperature in the upwelling area is much lower than its surrounding area. Compared with Figure 8a,b, the number of stations in Figure 8d is 4; compared with Figure 8c, the station spacing of Figure 8d is larger. We conclude that increasing the number of stations and station spacing may help to improve the accuracy of the temperature mapping. (4) The temperature distribution on Figure 8c did not change much before and after the upwelling, possibly because stations K1 and K2, K3 and K4 are positioned too close to each other, resulting in a non-ideal inversion process. However, due to the noticeable change in average temperature, upwelling phenomenon could still be detected.

Comparison of CAT Results with Temperature Sensor Results
To verify the CAT observation results, we also measured the water temperature in the vertical direction. First, five temperature sensors named Minilog-II-T were deployed vertically at the depth of 0, 2, 4, 6, and 8 m respectively, to form a temperature chain. The temperature chain was attached to the top of the gas injection frame. Second, the sensors transmit temperature data to the computer via Bluetooth. The observed results are compared with the former obtained CAT results as shown in Figure 9. the number of stations and station spacing may help to improve the accuracy of the temperature mapping. (4). The temperature distribution on Figure 8c did not change much before and after the upwelling, possibly because stations K1 and K2, K3 and K4 are positioned too close to each other, resulting in a non-ideal inversion process. However, due to the noticeable change in average temperature, upwelling phenomenon could still be detected.

Comparison of CAT Results with Temperature Sensor Results
To verify the CAT observation results, we also measured the water temperature in the vertical direction. First, five temperature sensors named Minilog-II-T were deployed vertically at the depth of 0, 2, 4, 6, and 8 m respectively, to form a temperature chain. The temperature chain was attached to the top of the gas injection frame. Second, the sensors transmit temperature data to the computer via Bluetooth. The observed results are compared with the former obtained CAT results as shown in Figure 9. Based on the comparison results in Figure 9, we can see that: 1. All figures in Figure 9 show that the CAT results are distributed between the depth of 0 and 2 m. This is consistent with the deployed depth of transducer, i.e., approximate 1 m. 2. The water temperature of the anechoic tank decreases from the surface to the bottom. Especially between 0 and 2 m, the temperature drops sharply. This phenomenon indicates that the closer the transducer is placed on the surface, the more accurate temperature changes it could capture. 3. Figure 9 also shows that the surface temperature is between 20.5 and 21 °C. After upwelling the water at depth of 4m, the surface temperature drops by about 3 °C. 4. After the upwelling, the CAT observation results are approximate to the water temperature at the surface. Before the upwelling begins, temperature from CAT observation fluctuates and remains within a certain range. After the upwelling, the temperature begins to drop, but stabilizes as the water flow with different depths is mixed. The temperature change captured by CAT before and after the upwelling is approximately 1 °C. 5. Figure 9d produces the best results in term of the least temperature fluctuating range with 0.5 °C in CAT observation results. This is also consistent with the conclusion from the inversion results in Section 4.2.

The Accuracy Analysis of CAT Observation Results
The accuracy of CAT observation results is affected by the resolution of travel time. The resolution of travel time r t using M sequence is calculated by: where Q is the cycle-per-digit of M sequence and f is the frequency. The resolution of sound speed deviation r C δ is expressed by: where L is the distance of two stations and 0 C is the reference sound speed. To further increase observation results accuracy, r C δ could be recalculated by: Based on the comparison results in Figure 9, we can see that: (1) All figures in Figure 9 show that the CAT results are distributed between the depth of 0 and 2 m. This is consistent with the deployed depth of transducer, i.e., approximate 1 m. (2) The water temperature of the anechoic tank decreases from the surface to the bottom. Especially between 0 and 2 m, the temperature drops sharply. This phenomenon indicates that the closer the transducer is placed on the surface, the more accurate temperature changes it could capture. (3) Figure 9 also shows that the surface temperature is between 20.5 and 21 • C. After upwelling the water at depth of 4m, the surface temperature drops by about 3 • C. (4) After the upwelling, the CAT observation results are approximate to the water temperature at the surface. Before the upwelling begins, temperature from CAT observation fluctuates and remains within a certain range. After the upwelling, the temperature begins to drop, but stabilizes as the water flow with different depths is mixed. The temperature change captured by CAT before and after the upwelling is approximately 1 • C. (5) Figure 9d produces the best results in term of the least temperature fluctuating range with 0.5 • C in CAT observation results. This is also consistent with the conclusion from the inversion results in Section 4.2.

The Accuracy Analysis of CAT Observation Results
The accuracy of CAT observation results is affected by the resolution of travel time. The resolution of travel time t r using M sequence is calculated by: where Q is the cycle-per-digit of M sequence and f is the frequency. The resolution of sound speed deviation δC r is expressed by: where L is the distance of two stations and C 0 is the reference sound speed. To further increase observation results accuracy, δC r could be recalculated by: where N is decided by the signal repeat number and is 8 as mentioned in Section 3. Within the temperature range of 17-21 • C, sound speed and temperature has an approximately linear relationship, i.e., per 1 m/s increase of sound speed means a 0.31 • C increase of water temperature. Therefore, the temperature deviation δT caused by the change of sound speedcould be calculated. The parameters of resolution calculation mentioned in Equations (14)- (16) are listed in the Table 1. δT is 0.33-1.32 • C when station-to-station distance is 7.5-30 m. This value agrees well with the overall temperature fluctuation range of CAT displayed in Figure 9.

Conclusions
In this study, we investigate a CAT method to carry out long-term and continuous temperature field observation in the artificial upwelling area in the anechoic tank. Four different CAT layouts were designed and used to observe the temperature in the experiments. The collected temperature data were processed with cross-correlation to obtain travel time. From travel time, we applied an inversion method to map the horizontal temperature field. From the CAT results and the comparison results, we concluded that: • The change in temperature field could be successfully observed by comparing the horizontal temperature field before and after the upwelling.

•
The CAT observation results are basically consistent with the temperature sensor results. This proves the feasibility of CAT for artificial upwelling observation.

•
The temperature field in the experimental area has the following characteristics before and after the upwelling. Before the upwelling, the water temperature in the vertical direction changes obviously. After upwelling, the water temperature on the surface decreases rapidly by approximately 3 • C, and stabilizes due tothe exchange of water flow from different depths.
These characteristics are considered an important reference for subsequent upwelling observation in sea experiments. Moreover, there is still space to improve. For example, the intervals for sending acoustic signals is too long currently. This may lead to the temperature field changes remaining not fully captured at the beginning of the upwelling. Furthermore, each station is required to take turns sending signals because simultaneous transmission, and reception was not possible. This may influence the observation accuracy. In the future, we will conduct offshore experiments with larger station spacing and more CAT stations to obtain long-term and continuous temperature field observations in the artificial upwelling area.
Author Contributions: The work described in this article is the collaborative development of all authors. H.H. contributed to the idea of experiments setup and designed the methods. Y.G. and Z.W. carried out the tank experiments. Z.W. and Y.S. made contributions to data collection, processing and analysis. Y.G. and Y.W. participated in the drafting and revising of the paper.