Measurements, Simulations, Python Scripts, and Results Presented in "Estimating Soil Electrical Parameters in the Canadian High Arctic from Impedance Measurements of the MIST Antenna Above the Surface"

README This repository corresponds to measurements, simulations, Python scripts, and results presented in "Estimating Soil Electrical Parameters in the Canadian High Arctic from Impedance Measurements of the MIST Antenna Above the Surface" (Hendricksen et al. 2025), currently undergoing review for consideration of publication at the time of upload of this repository. If accepted, this repository will remain available such that the results are reproducible by anyone. If you wish to use any of the data found in this repository, please use that corresponding to Version 2 (V2). Ensure to extract all the zipped folders in the same directory since some scripts include relative paths. Data and Simulations The impedance measurements conducted with the MIST antenna from the study site in the Canadian High Arctic are contained in impedance_measurements.tar.gz The impedance measurements are stored in HDF5 files which can be opened manually with h5py, or through methods provided in soil.py and demonstrated in run_fits.py (see below). The Feko simulations are contained in feko_simulations.tar.gz which contains two subdirectories corresponding to the two soil models: The single-layer simulations are contained in feko_simulations/csa2022_2parameters_9samples_V3_20240128. The two-layer simulations are contained in feko_simulations/csa2022_5parameters_4samples_V3_20240108. It is likewise recommended to use the simulation loading methods demonstrated in run_fits.py to access the simulation data. Results Files The results are contained in results.tar.gz The nominal results are stored in results/single_layer/curve_fit/cubic_interp_absolute_sigma_True_mean_p0/results.hdf5: single-layer results/two_layer/curve_fit/cubic_interp_absolute_sigma_True_mean_p0/results.hdf5: two-layer One can access the different datasets contained in each file by checking the keys of the file with h5py. The included datasets are BIC: the Bayesian information criterion. bounds: the bounds of the spans for each soil parameter as displayed in the "Span" column of Table S2. chisq_red: the reduced chi-squared χred2χred2. data: the measured, calibrated antenna impedance as a function of time (axis 0) and frequency (axis 1). day_time: the measurement time referenced to the first measurement in units of days. frequency: the frequency bins spanning 25-125 MHz. initial_guess: the starting parameters for the fitting process. mist_temp: the temperature of the internal 50-ΩΩ load within the MIST receiver box. models: the best-fit models as a function of time (axis 0) and frequency (axis 1). p: the best-fit parameters. The parameter names can be accessed through the attribute param_names of this dataset. perr: the 68% confidence level of each parameter. s11_index: the concatenated indices of the impedance measurements for each dataset. Used primarily for debugging. Note that the first measurement (i.e., index 0) is typically conducted while humans are still operating the instrument, whose nearby presence biases the measured impedance; therefore, the first measurements are not considered in this analysis. total_uncertainty: the total uncertainty in ΩΩ used in the fitting process. The first (second) row corresponds to the real/resistance (imaginary/reactance) part of the total uncertainty. weatherstation_temp: the air temperature measured by the weather station, interpolated at the timestamps at which the MIST internal temperatures were recorded. Note that the full weather station data are made available (discussed in Python Codes below). Note that for each key where there is a frequency axis, the resistance and reactance are concatenated together, such that the first 101 indices correspond to the resistance from 25-125 MHz, and the second 101 indices correspond to the reactance from 25-125 MHz, for a total of 202 indices. The other files contained in the results directory correspond to tests of sensitivity of the estimated soil parameters to assumptions we have made in this analysis, whose results are summarized in Table S3. Specifically, they include: results/single_layer/curve_fit/cubic_interp_fix_L_absolute_sigma_True_mean_p0: the residual inductance L is fixed for the single-layer model. results/two_layer/curve_fit/cubic_interp_fix_L_absolute_sigma_True_mean_p0: the residual inductance L is fixed for the two-layer model. results/two_layer/curve_fit/cubic_interp_tighten_d_bounds_absolute_sigma_True_mean_p0: the bounds of the parameter search for the top-layer thickness tt are constrained to within ±1±1 sample standard deviation of the metal probe measurements. Python Codes The Python codes for fitting interpolated models to the data and generating plots are contained in python_codes.tar.gz The scripts require the following packages to be run (most of which are standard): numpy matplotlib scipy h5py datetime os glob argparse re json PyYAML pandas There are two main "helper" scripts which provide the software tools to interact with data and generate impedance models. These include soil.py, the main software tool which loads data, simulations, uncertainties, etc., and organizes fitting routines, and ampere.py, a generalized N-dimensional interpolation class called by soil.py when forming impedance models to fit to measured data. There are several analysis scripts included along with the software tools, including compare_nominal_and_constrained_rmsd.py: used to generate the values in Table S3. make_data_plot.py: used to generate a raw version of Figure 2. make_models_and_residuals_plot.py: used to generate Figure 3. make_parameter_plots.py: used to generate raw versions of Figures 4 and S5. make_simulation_plot.py: used to generate Figure S3. make_total_uncertainty_plot.py: used to generate Figure S4. run_fits.py: used to obtain the best-fit results presented in Hendricksen et al. (2025). It is recommended to start here if you would like to learn how to load data and simulations and fit models to data. NOTE: if this file is run with save_data = True (line 44), the results files contained in results.tar.gz (see the next section) will be overwritten. The uncertainties considered in this analysis are contained in the subdirectory python_codes/uncertainty_files, which includes python_codes/uncertainty_files/csa_2022_single_layer_all_uncertainties_20241130.txt: single-layer uncertainties. python_codes/uncertainty_files/csa_2022_two_layer_all_uncertainties_20241130.txt: two-layer uncertainties. Note that there are some quantities which are repeated between the two uncertainty files, such as δcalδcal, δmeasδmeas, and δFekoδFeko. Finally, the weather station data is contained in python_codes/csa_weather_station_data_2022-05-03_2023-04-25_original.xlsx which includes key metrics presented in the paper, including AirTemp_Avg: air temperature in ∘C∘C RH_Avg: relative humidity as a percentage Solar_W_Avg: solar radiation in Wm−2Wm−2 mean_wind_speed: mean wind speed in kmh−1kmh−1 The weather station data is loaded using the pandas package. An interpolated subset of the air temperature data corresponding to the same time that data was taken is included in each of the results files. Figures All the figures presented in the paper are contained in figures.tar.gz along with additional scripts and files to produce figures with broken axes, including Figures 2, 4, and S5, which require post-processing. Note that imagemagick and Inkscape are required to perform the post-processing described below. Figure 1 Figure 1 corresponds to mist_in_arctic_and_feko.pdf. Figure 2 Figure 2 corresponds to example_impedance_and_variation_crop.pdf. To reproduce Figure 2, run python_codes/make_data_plot.py, which will produce two temporary figures example_impedance_and_variation_1.jpg and example_impedance_and_variation_2.jpg. Run the shell script crop_example_impedance_and_variation.sh to produce example_impedance_and_variation_crop.jpg, an intermediate, cropped version of Figure 2, which needs to be polished after cropping. A scalable vector graphics (SVG) file example_impedance_and_variation_crop.svg is provided with example_impedance_and_variation_crop.jpg imported as its base layer. The SVG file can be opened with Inkscape and used to export the final PDF version of Figure 2, example_impedance_and_variation_crop.pdf. Figure 3 Figure 3 corresponds to example_fits_and_all_residuals.pdf. Figure 4 Figure 4 corresponds to best_fit_parameters_two_layer_crop.pdf. To reproduce Figure 3, run python_codes/make_parameter_plots.py, which will produce a temporary figure best_fit_parameters_two_layer.png. Run the shell script crop_best_fit_parameters_two_layer.sh to produce best_fit_parameters_two_layer_crop.png, an intermediate, cropped version of Figure 4, which needs to be polished after cropping. An SVG file best_fit_parameters_two_layer_crop.svg is provided with best_fit_parameters_two_layer_crop.png imported as its base layer. The SVG file can be opened with Inkscape to export the final PDF version of Figure 4, best_fit_parameters_two_layer_crop.pdf. Figure S1 Figure S1 corresponds to axel_heiberg_and_mars_region.pdf. Figure S2 Figure S2 corresponds to antenna.pdf, and is reused from Monsalve et al. (2024) (see the References.) Figure S3 Figure S3 corresponds to all_feko_sims.pdf. To reproduce Figure S3, run python_codes/make_simulation_plot.py. Figure S4 Figure S4 corresponds to total_uncertainty.pdf. To reproduce Figure S4, run python_codes/make_total_uncertainty_plot.py. Figure S5 Figure S5 corresponds to parameter_line_fits_crop.pdf. To reproduce Figure S5, run python_codes/make_parameter_plots.py for the two-layer model, which will produce a temporary figure parameter_line_fits.png. Run the shell script crop_parameter_line_fits.sh to produce parameter_line_fits_crop.png, an intermediate, cropped version of Figure S5, which needs to be polished after cropping. An SVG file parameter_line_fits_crop.svg is provided with parameter_line_fits_crop.png imported as its base layer. The SVG file can be opened with Inkscape to export the final PDF version of Figure S5, parameter_line_fits_crop.pdf.# README 本仓库对应发表于《基于地表MIST天线的阻抗测量估算加拿大高北极地区土壤电学参数》（Hendricksen等人，2025）的测量数据、仿真结果、Python脚本及研究成果，本仓库上传时该论文仍处于审稿待录用阶段。 若该论文最终被录用，本仓库将持续开放以确保所有研究结果均可被复现。若需使用本仓库中的数据，请选用对应版本2（V2）的数据集。请将所有压缩文件夹解压至同一目录，部分脚本依赖相对路径运行。 ## 数据与仿真 加拿大高北极地区研究站点使用MIST天线采集的阻抗测量数据，存放于`impedance_measurements.tar.gz`压缩包中。 该阻抗测量数据以HDF5（Hierarchical Data Format 5）格式存储，可通过h5py库手动打开，或通过soil.py中提供的方法进行读取，具体使用示例可参考run_fits.py（详见下文）。 Feko电磁仿真数据存放于`feko_simulations.tar.gz`压缩包中，该压缩包包含两个子目录，分别对应两种土壤模型： - 单层土壤模型仿真数据存放于`feko_simulations/csa2022_2parameters_9samples_V3_20240128`目录下。 - 双层土壤模型仿真数据存放于`feko_simulations/csa2022_5parameters_4samples_V3_20240108`目录下。 同样建议使用run_fits.py中演示的仿真数据加载方法来读取仿真数据。 ## 结果文件 研究结果存放于`results.tar.gz`压缩包中。 标准结果分别存放于以下路径： - 单层土壤模型：`results/single_layer/curve_fit/cubic_interp_absolute_sigma_True_mean_p0/results.hdf5` - 双层土壤模型：`results/two_layer/curve_fit/cubic_interp_absolute_sigma_True_mean_p0/results.hdf5` 可通过h5py库查看文件的键（keys），从而获取每个文件中包含的各类数据集。本仓库包含的数据集说明如下： 1. **BIC**：贝叶斯信息准则（Bayesian Information Criterion） 2. **bounds**：各土壤参数的搜索范围跨度，对应附表S2中“Span”列的取值 3. **chisq_red**：约化卡方统计量（reduced chi-squared, $chi^2_{red}$） 4. **data**：经过校准的实测天线阻抗数据，维度为时间（第0轴）× 频率（第1轴） 5. **day_time**：以首个测量时刻为基准的测量时间，单位为天 6. **frequency**：覆盖25~125 MHz频段的频率分档 7. **initial_guess**：拟合流程的初始参数猜测值 8. **mist_temp**：MIST接收机机箱内部50 Ω负载的温度 9. **models**：最优拟合阻抗模型数据，维度为时间（第0轴）× 频率（第1轴） 10. **p**：最优拟合参数。可通过该数据集的`param_names`属性获取参数名称 11. **perr**：各参数的68%置信区间误差 12. **s11_index**：各数据集对应的阻抗测量数据拼接索引，主要用于调试。注意：首个测量数据（即索引0）通常在操作人员未撤离仪器时采集，人员近距离干扰会导致阻抗测量结果出现偏差，因此本分析未纳入首个测量数据。 13. **total_uncertainty**：拟合过程中使用的总不确定度（单位：Ω）。第一行对应阻抗实部（电阻分量）的总不确定度，第二行对应阻抗虚部（电抗分量）的总不确定度。 14. **weatherstation_temp**：气象站测量的气温数据，经插值后与MIST内部温度记录的时间戳对齐。完整的气象站数据详见下述Python代码章节。 注意：对于包含频率轴的数据集，阻抗实部（电阻分量）与虚部（电抗分量）会被拼接存储：前101个索引对应25~125 MHz频段的电阻分量，后101个索引对应25~125 MHz频段的电抗分量，总索引数为202。 结果目录中包含的其他文件，用于测试分析假设对土壤参数估算结果的敏感性，其结果汇总于附表S3，具体包括： - `results/single_layer/curve_fit/cubic_interp_fix_L_absolute_sigma_True_mean_p0`：单层土壤模型分析中固定残余电感L的拟合结果 - `results/two_layer/curve_fit/cubic_interp_fix_L_absolute_sigma_True_mean_p0`：双层土壤模型分析中固定残余电感L的拟合结果 - `results/two_layer/curve_fit/cubic_interp_tighten_d_bounds_absolute_sigma_True_mean_p0`：双层土壤模型分析中将表层厚度$t$的参数搜索范围约束至金属探针测量样本标准差的±1倍区间内的拟合结果 ## Python代码 用于对数据拟合插值模型并生成绘图的Python代码，存放于`python_codes.tar.gz`压缩包中。 运行以下脚本需依赖以下Python包（多数为标准库）： - numpy - matplotlib - scipy - h5py - datetime - os - glob - argparse - re - json - PyYAML - pandas 本仓库包含两个核心“辅助”脚本，用于实现数据交互与阻抗模型生成： 1. **soil.py**：核心工具脚本，负责加载数据、仿真结果、不确定性参数等，并组织拟合流程 2. **ampere.py**：广义N维插值类，soil.py在构建拟合实测数据的阻抗模型时会调用该类 此外还包含多个分析脚本： 1. **compare_nominal_and_constrained_rmsd.py**：用于生成附表S3中的数值结果 2. **make_data_plot.py**：用于生成图2的原始版本 3. **make_models_and_residuals_plot.py**：用于生成图3 4. **make_parameter_plots.py**：用于生成图4及附表图S5的原始版本 5. **make_simulation_plot.py**：用于生成附表图S3 6. **make_total_uncertainty_plot.py**：用于生成附表图S4 7. **run_fits.py**：用于获取Hendricksen等人（2025）论文中展示的最优拟合结果 若希望学习如何加载数据与仿真结果、并对数据拟合模型，建议从该脚本开始学习。 注意：若将save_data参数设为True（第44行），则`results.tar.gz`压缩包中包含的结果文件将被覆盖（详见本章节前述内容）。 本分析中使用的不确定性参数存放于`python_codes/uncertainty_files`子目录中，包含以下文件： - `python_codes/uncertainty_files/csa_2022_single_layer_all_uncertainties_20241130.txt`：单层土壤模型不确定性参数文件 - `python_codes/uncertainty_files/csa_2022_two_layer_all_uncertainties_20241130.txt`：双层土壤模型不确定性参数文件 注意：两个不确定性参数文件中包含部分重复的量化指标，如$delta_{cal}$、$delta_{meas}$及$delta_{Feko}$。 最后，气象站数据存放于`python_codes/csa_weather_station_data_2022-05-03_2023-04-25_original.xlsx`文件中，该文件包含论文中提及的关键指标： - **AirTemp_Avg**：平均气温，单位为摄氏度（℃） - **RH_Avg**：平均相对湿度，单位为百分比 - **Solar_W_Avg**：平均太阳辐射，单位为$W/m^2$ - **mean_wind_speed**：平均风速，单位为$km/h$ 气象站数据可通过pandas库加载。各结果文件中已包含与MIST内部温度记录时间戳对齐的插值气温数据子集。 ## 配图 论文中所有配图均存放于`figures.tar.gz`压缩包中，同时附带用于生成带断裂轴绘图的脚本与文件，包括图2、图4及附表图S5，此类绘图需进行后期处理。注意：执行下述后期处理需依赖imagemagick与Inkscape工具。 ### 图1 图1对应文件为`mist_in_arctic_and_feko.pdf`。 ### 图2 图2对应文件为`example_impedance_and_variation_crop.pdf`。复现图2需执行以下步骤： 1. 运行`python_codes/make_data_plot.py`，将生成两个临时图片`example_impedance_and_variation_1.jpg`与`example_impedance_and_variation_2.jpg` 2. 执行shell脚本`crop_example_impedance_and_variation.sh`，生成裁剪后的中间文件`example_impedance_and_variation_crop.jpg`，该文件需进一步润色 3. 提供的SVG文件`example_impedance_and_variation_crop.svg`已将`example_impedance_and_variation_crop.jpg`作为基础图层导入，可通过Inkscape打开该SVG文件并导出为最终的PDF格式配图`example_impedance_and_variation_crop.pdf` ### 图3 图3对应文件为`example_fits_and_all_residuals.pdf`。 ### 图4 图4对应文件为`best_fit_parameters_two_layer_crop.pdf`。复现图4需执行以下步骤： 1. 运行`python_codes/make_parameter_plots.py`，将生成临时图片`best_fit_parameters_two_layer.png` 2. 执行shell脚本`crop_best_fit_parameters_two_layer.sh`，生成裁剪后的中间文件`best_fit_parameters_two_layer_crop.png`，该文件需进一步润色 3. 提供的SVG文件`best_fit_parameters_two_layer_crop.svg`已将`best_fit_parameters_two_layer_crop.png`作为基础图层导入，可通过Inkscape打开该SVG文件并导出为最终的PDF格式配图`best_fit_parameters_two_layer_crop.pdf` ### 附表图S1 附表图S1对应文件为`axel_heiberg_and_mars_region.pdf`。 ### 附表图S2 附表图S2对应文件为`antenna.pdf`，该图复用自Monsalve等人（2024）的研究（详见参考文献）。 ### 附表图S3 附表图S3对应文件为`all_feko_sims.pdf`。复现该图需运行`python_codes/make_simulation_plot.py`。 ### 附表图S4 附表图S4对应文件为`total_uncertainty.pdf`。复现该图需运行`python_codes/make_total_uncertainty_plot.py`。 ### 附表图S5 附表图S5对应文件为`parameter_line_fits_crop.pdf`。复现该图需执行以下步骤： 1. 针对双层土壤模型运行`python_codes/make_parameter_plots.py`，将生成临时图片`parameter_line_fits.png` 2. 执行shell脚本`crop_parameter_line_fits.sh`，生成裁剪后的中间文件`parameter_line_fits_crop.png`，该文件需进一步润色 3. 提供的SVG文件`parameter_line_fits_crop.svg`已将`parameter_line_fits_crop.png`作为基础图层导入，可通过Inkscape打开该SVG文件并导出为最终的PDF格式配图`parameter_line_fits_crop.pdf`