USING ELECTROMAGNETIC FIELD ANALYSIS FOR SEA ICE THICKNESS MEASUREMENT

Observation of sea ice thickness by remote sensing is one of key issues to understand an effect of global warming. However, ground truth must be necessary to discuss this kind of approach. Although there are several methods to acquire ice thickness, Ground Penetrating Radar (GPR) can be good solution because it can discriminate snow-ice and ice-sea water interface thanks to comparative higher spatial resolution than the other methods. In this paper, we carried out GPR measurement in brackish lake and an electromagnetic field analysis in order to interpret the GPR data. The simulation model was assumed considering the actual snow and ice thickness acquired in field measurement. From the simulation results, although it seems difficult to identify the reflection at snow and ice interface due to a thin layer thickness and a low dielectric constant, snow and ice thickness may be estimated by using multiple reflection components.


INTRODUCTION
Sea ice thickness is one of important parameters to understand an effect of global warming.The remote sensing approach such as Synthetic Aperture Radar (SAR) has possibility to estimate ice thickness (Wakabayashi et al., 2004a).However, still ground truth must be necessary to discuss the remote sensing approach.
There are several methods to obtain the ground truth of ice and snow thickness.Simplest but most correct way is making a hole and measuring the thickness directly by a ruler.However this method is destructive, time consuming and gives only point information.Although, there are other methods in order to cover larger area such as electromagnetic induction, sonar and so on, ground penetrating radar (GPR) can be good solution because it can discriminate snow-ice and ice-sea water interface thanks to comparative higher spatial resolution than the other methods (Liu et al., 2014a).GPR is a tool in order to image a subsurface structure (shallower than several meters) such as buried object, layer structure and so on by using a penetration and reflection of microwave.The GPR can also be used in ice thickness measurement (Tebaldini et al., 2016a).
On Feb. 2018, we carried out GPR measurement in Lake Saroma, which is a brackish lake, for relatively large area (200m by 300m, approximately) aiming to obtain grand truth for remote sensing data (Matsumoto et al., 2018).Although several Bscan images could be obtained, it was difficult to identify where reflection signals come from.In this paper, we carried out an electromagnetic filed analysis in order to interpret Bscan image acquired in field measurement.Firstly, methodology for ice thickness estimation and field experiment in Lake Saroma will be explained.Considering the result of the direct measurement, the simulation model was assumed and the electromagnetic field analysis was implemented.The interpretation of simulated result will be presented.

METHODOROGY
Figure 1 shows skematic view of ice thickness measurement by GPR.A three layer structure, which consists of snow, ice and water, is assumed.GPR transmits electromagnetic wave into subsurface and received reflected wave caused by snow-ice and ice-water interface.GPR records amplitude and arrival time of reflection signal.Therefore if we identify the each reflection signal, layer thickness H can be determined by where ∆t = time difference of reflection signal v = propagation velocity in a medium and given by where c = speed of electromagnetic wave in vacuum  r = dielectric constant of a medium Note that the propagation velocity, which is determined by dielectric constant of material, is still unknown.
Our methodology to estimate snow and ice thickness is as following.Firstly, we carry out GPR survey and acquire cross section image, which is called Bscan image, by collecting a reflection signal along survey line.The direct measurement on the survey line should be conducted just after GPR survey at several points for calibration and validation.If the reflection signal from interface can be identified in Bscan image, time difference of the signal also can be determined.By considering the actual thickness of snow and ice, the propagation velocity at which direct measurement carried out can be determined.Therefore snow and ice thickness along survey line can be estimated by assuming that the propagation velocity is constant.
In this method, the important points are 1) Whether the reflection signal can be identified, and 2) How the arrival time of the signal can be determined.
An electromagnetic field analysis was conducted in order to confirm these points.The gprMax which is a free numerical GPR simulator based on FDTD (Finite-Difference Time-Domain) method (Warren et al., 2016a) was used for this purpose.We assumed a three layer model which consists of snow, ice and sea water, having layer thickness of 20cm, 20cm and 40cm, respectively.The dielectric constant of them is assumed as 2, 3 and 81, and the conductivity is assumed as 10 - 5 S/m, 10 -3 S/m, and 100S/m, respectively by referring to a list of electric conductivity (Daniels, 2004).A case of dry snow was assumed in this simulation.This case may be more difficult because travel time of snow layer is shorter due to lower dielectric constant resulting in a difficulty to discriminate the reflection from snow-ice interface and direct antenna coupling.A 600MHz GPR system was modelled by using a first derivative of Gaussian pulse with central frequency of 600MHz as a feeding source.By calculating the time difference of peak of reflection signal, the thickness of snow and ice were estimated as 18.95cm and 19.23cm, respectively.The results were less than 5% error and almost agreed with actual thickness in the simulation model.Therefore this method can be work if we identify the peak of reflection signal from layer interface.

Test Field
Field experiment was carried out in Lake Saroma, north east of Hokkaido, Japan, on 18 th to 22 of February 2018.Lake Saroma is facing to Sea of Okhotsk and connected by two inlets resulting in brackish lake.Figure 3

Result
According to the direct measurement, the overall average and standard deviation of snow were 10.1cm and 2.2cm, and those of ice were 26.0cm and 3.0cm, respectively.We also observed the condition between snow and ice layer.As for No.1, there was a mixture of snow and water above ice layer at distance of 9m to 18m.As for No.2 and No.4, we found that snow just above ice layer is bit wet condition at all points.Although these signal may be a reflection which is related to snow-ice and ice-water interface, it is not obvious where the reflections come from.Besides, since the average snow thickness is thinner than the simulation model in chapter 2, the method to estimate thickness may be different to apply.Therefore the simulation model must be adjusted by considering the actual thickness acquired in the direct measurement in order to interpret the GPR data.

NUMERICAL ANALYSIS
Aiming to interpret the Bscan images, we carried out an electromagnetic field analysis.We assumed layer thickness of snow, ice and sea water is 10cm, 26cm and 44cm, respectively, considering the result of direct measurement.The dielectric constant and the conductivity are same value mentioned in chapter 2. We call this model as Model 1. Besides, we assumed other two models, Model 2 and Model 3 which have different layer thickness in order to identify reflection signal from layer interface.Figure 6 shows these simulation models.Compared to Model 1, Model 2 has same thickness of snow layer but 10cm thicker ice, on the other hand, Model 3 has 10cm thicker snow but same thickness of ice.
Figure 7 shows path of reflection signal in simulated models.We assumed six reflection paths such as P 0 to P 5 as shown in this figure.Three interfaces such as air-snow, snow-ice and icesea water are represented as B 1 , B 2 and B 3 , respectively.Twoway travel time in the mediums is represented as T 1 , T 2 and T 3 , respectively.P 0 and T 0 indicate antenna direct coupling and its travel time.P 1 , P 2 and P 3 are paths reflected at B 1 , B 2 and B 3 , respectively.P 4 and P 5 indicate a path of multiple reflections.As for P 4 , firstly, the electromagnetic wave is reflected at B 3 .After that the electromagnetic wave is reflected at B 2 and B 3 then back to the antenna.As for P 5 , the electromagnetic wave is reflected at B 3 , B 1 and B 3 then back to the receiver.The reflection coefficient of layer interface of medium1 and medium 2 can be calculated by following formula.
Figure 8 shows simulation results of assumed models.Analytic signal representation is used in order to discuss the arrival time of the peak signal.The four main reflections were observed in Figure 8 (a).Although the peak signal which may come from snow-ice interface can be observed in Figure 2, this peak cannot be observed in Figure7 (a) due to thinner snow thickness.Therefore it is not simple to calculate layer thickness compared to the case of Figure 2.
Let us focus on an analysis of Figure 8 (a), i.e., Model 1.First peak can be mainly a direct antenna coupling because the direct antenna coupling has much larger amplitude than the other reflection signals.The reflection from air-snow interface is much weaker and distance between antenna and snow is too close, the reflection path P 1 cannot be identified.The reflection path P 2 seems to be also included at the first peak due to thin snow thickness.If snow thickness is thicker as Model 3, the reflection path P 2 can be observed as shown in Figure 8 (c).Therefore we concluded that reflection signal at snow-ice interface may not be observed in measurement data due to thin snow thickness.Summarizing the interpretation above, the first peak includes the reflection path P 0 , P 1 and P 2 .
Second peak of signal which has second strongest signal amplitude can be reflection path P 3 .Observing the result of Model 2 and Model 3 which are the deeper water level from top of snow as shown in Figure 8 (b) and (c), arrival time of this peak are later than the case of Model 1. Besides, absolute value of reflection coefficient at B 3 is larger than the other interface.This can be the reason why the signal amplitude is stronger than the others.Therefore we concluded that second peak of the signal comes from ice-sea water interface, i.e., reflection path P 3 .
Third and fourth peak of signal can be considered that they correspond to reflection path P 4 and P 5 , i.e., multiple reflections.If snow is thicker as Model 3, time difference of P 4 and P 5 is larger as shown in Figure 8 (c).The time difference of them is 0.85nsec for Model 1, 0.90nsec for Model 2 and 1.93nsec for Model 3, respectively.If ice is thicker as Model 2, time difference of P 3 and P 4 is larger as shown in Figure 8 (b).The time difference of them is 3.03nsec for Model 1, 4.13nsec for Model 2 and 2.96nsec for Model 3, respectively.From this observation, third and fourth peak of signal can be considered as reflection path P 4 and P 5 .
The travel time of each reflection path is shown in Figure 7. From this figure, the snow thickness can be estimated by the time difference of P 4 and P 5 .The ice thickness can be estimated by that of P 3 and P 4 .As for Model 1, the snow and ice thickness was estimated as 9.0cm and 26.22cm, respectively.The errors of them are 10% and 0.85%, respectively.Even though the reflection at snow-ice cannot be identified, we confirmed that the snow and ice thickness can be estimated using time difference of second, third and fourth peak of signal.Therefore we concluded that multiple reflections are useful information for the thickness estimation if the snow layer is thin and the reflection from snow-ice interface cannot be identified directly.The comparison of simulated data and measured GPR data will be a future work.

CONCLUSION
In this paper, an electromagnetic field analysis was carried out in order to interpret GPR data which was acquired as a ground truth of snow and ice thickness for remote sensing data.Simulation model was assumed by considering the result of direct measurement.From the simulation results, although it seems difficult to identify the reflection at snow and ice interface due to thin layer thickness and low dielectric constant, snow and ice thickness may be estimated by using multiple reflection components.The comparison of simulated data and measured GPR data will be a future work.
Figure 2 shows simulated signal.Absolute value of analytic signal representation is used in order to detect a single peak of the reflection signal.The reflection signals at  1 ,  2 and  3 correspond to Path 1, Path 2 and Path 3 in Figure 1, respectively.
(a)  shows Satellite image around the test field located in west side of the lake.The area indicated by an arrow can be discriminated by the crack as shown in Figure3(b).Figure3(c) shows location of GPR survey lines in the right side of the crack area indicated by rectangular in Figure3(a).GPR survey was conducted at five locations in the area.

Figure 4
Figure4shows a GPR system used in the experiment.The radar unit is "RIS Hi-Mod", product of IDS GeoRadar s.r.l.The Bscan image of 600MHz was obtained by collecting the reflection signal at every 3cm using hip chain encoder along the scan direction.A direct measurement along survey line was also conducted in No.1, No.2 and No.4 with interval of 3m as shown in Figure3 (c).We drilled and measured thickness of snow and ice by a ruler directly.

Figure 3 .
Figure 3. Test field located in the west side of the lake.(a) Sattellite image taken by Sentinel-2 (Jan.18 th , 2018).(b) Crack arround the area.(c) GPR survey line arround the rectanguler in (a).The arrows indicate scan direction.

Figure 4 .
Figure 4.The GPR system used in the experiment.

Figure 5
Figure 5 shows GPR Bscan images acquired in survey line No.1, No.2 and No.4, respectively.Horizontal axis indicates GPR position acquired by hip chain encoder.Vertical axis indicates arrival time.The layer like structure can be seen in part of Figure 5 (a) and (b), and whole of (c) as indicated by arrows.Although these signal may be a reflection which is related to snow-ice and ice-water interface, it is not obvious where the reflections come from.Besides, since the average snow thickness is thinner than the simulation model in chapter 2, the method to estimate thickness may be different to apply.Therefore the simulation model must be adjusted by considering the actual thickness acquired in the direct measurement in order to interpret the GPR data.
Bscan image obtained from 600MHz antenna.

Figure
Figure 6.Simulation model