Note
Go to the end to download the full example code.
CSV & Rasters: Multi-dataset Survey with Derivative Products
This example demonstrates the typical workflow for creating a GS file for an AEM survey in its entirety, i.e., the NetCDF file contains all related datasets together, e.g., raw data, processed data, inverted models, and derivative products. Specifically, this survey contains:
Minimally processed (raw) AEM data and raw/processed magnetic data provided by SkyTEM
Fully processed AEM data used as input to inversion
Laterally constrained inverted resistivity models
Point-data estimates of bedrock depth derived from the AEM models
Interpolated magnetic and bedrock depth grids
Note: To make the size of this example more managable, some of the input datasets have been downsampled relative to the source files in the data release referenced below.
Source Reference: Minsley, B.J, Bloss, B.R., Hart, D.J., Fitzpatrick, W., Muldoon, M.A., Stewart, E.K., Hunt, R.J., James, S.R., Foks, N.L., and Komiskey, M.J., 2022, Airborne electromagnetic and magnetic survey data, northeast Wisconsin (ver. 1.1, June 2022): U.S. Geological Survey data release, https://doi.org/10.5066/P93SY9LI.
import matplotlib.pyplot as plt
from os.path import join
import numpy as np
import gspy
from gspy import Survey, Metadata
import xarray as xr
from pprint import pprint
import warnings
warnings.filterwarnings('ignore')
Initialize the Survey
# Path to example files
data_path = '..//data_files//skytem_csv'
# Survey metadata file
metadata = join(data_path, "data//skytem_survey.yml")
# Establish the Survey
survey = Survey.from_dict(metadata)
Survey YAML file
1dataset_attrs:
2 title: SkyTEM Airborne Electromagnetic (AEM) Survey, Northeast Wisconsin Bedrock Mapping
3 institution: USGS Geology, Geophysics, and Geochemistry Science Center
4 source: SkyTEM raw data, USGS processed data and inverted resistivity models, and depth to bedrock surface
5 history: (1) Data acquisition 01/2021 - 02/2021 by SkyTEM Canada Inc.; (2) AEM and magnetic data processing by SkyTEM Canada Inc. 02/2021 - 03/2021; raw and minimally processed AEM data, and processed magnetic data, received by USGS from SkyTEM Canada Inc 03/2021; Minimally processed AEM data exported to netCDF group /survey/data/raw_data group; (3) Minimally processed binary data and system response information received from the contractor were imported into the Aarhus Workbench software (v 6.0.1.0) where data were processed by USGS 03/2021 - 06/2021. Processed AEM data exported to netCDF group /survey/data/processed_data; (4) Processed data were inverted in Aarhus Workbench software using laterally constrained inversion to recover 40-layer fixed depth blocky resistivity models by USGS 03/2021 - 06/2021; Inverted resistivity models exported to netCDF group /survey/models/inverted_models. (5) Resistivity models were imported into the Geoscene3D software (v. 12.0.0.680) and points were generated at the first depth where resistivity exceeded 325 ohm-meters. These points were visually inspected and manually adjusted in selected areas to produce an AEM-derived estmiate of the elevation of the top of bedrock by USGS together with WGNHS 06/2021 - 07/2021. Points were exported to netCDF group /survey/data/depth_to_bedrock. (6) Bedrock elevation points were interpolated using kriging in Geoscene3D software to produce a regular bedrock elevation grid 07/2021. (7) A bedrock depth grid was calculated in QGIS software (v. 3.14.1-Pi) by subtracting the bedrock elevation from land surface elevation. (8) Bedrock elevation, bedrock depth, and SkyTEM-provided magnetic grids were aligned to a common 100m x 100m grid and exported to netCDF group /survey/derived_products/maps.
6 references: Minsley, Burke J., B.R. Bloss, D.J. Hart, W. Fitzpatrick, M.A. Muldoon, E.K. Stewart, R.J. Hunt, S.R. James, N.L. Foks, and M.J. Komiskey, 2021, Airborne electromagnetic and magnetic survey data, northeast Wisconsin, 2021, U.S. Geological Survey data release, https://doi.org/10.5066/P93SY9LI.
7 comment: This dataset includes minimally processed (raw) AEM and raw/processed magnetic data provided by SkyTEM, fully processed data used as input to inversion, laterally constrained inverted resistivity models, and derived estimates of bedrock depth.
8 summary: Airborne electromagnetic (AEM) and magnetic survey data were collected during January and February 2021 over a distance of 3,170 line kilometers in northeast Wisconsin. These data were collected in support of an effort to improve estimates of depth to bedrock through a collaborative project between the U.S. Geological Survey (USGS), Wisconsin Department of Agriculture, Trade, and Consumer Protection (DATCP), and Wisconsin Geological and Natural History Survey (WGNHS). Data were acquired by SkyTEM Canada Inc. with the SkyTEM 304M time-domain helicopter-borne electromagnetic system together with a Geometrics G822A cesium vapor magnetometer. The survey was acquired at a nominal flight height of 30 - 40 m above terrain along parallel flight lines oriented northwest-southeast with nominal line spacing of 0.5 miles (800 m). AEM data were inverted to produce models of electrical resistivity along flight paths, with typical depth of investigation up to about 300 m and 1 - 2 m near-surface resolution. Shallow resistivity transitions were used to estimate depth to bedrock across the survey area.
9 content: Wisconsin SkyTEM survey information
10
11survey_information:
12 contractor_project_number: 20022
13 contractor: SkyTEM Canada Inc
14 client: U.S. Geological Survey
15 survey_type: EM/Mag
16 survey_area_name: Northeast Wisconsin Bedrock Mapping
17 state: WI
18 country: USA
19 acquisition_start: "20210117"
20 acquisition_end: "20210207"
21 survey_attributes_units: SI
22
23spatial_ref:
24 wkid: 3071
25 authority: EPSG
26 vertical_crs: NAVD88
27
28flightline_information:
29 traverse_line_spacing: 800 m
30 traverse_line_direction: north-west to south-east, 114 degrees / 294 degrees
31 tie_line_spacing: No tie lines were flown in this survey
32 nominal_terrain_clearance: 30 m
33 final_line_kilometers: 3170 km
34 traverse_line_numbers: 100101 - 115201
35 repeat_line_numbers: 920001 - 920006
36 zero_line_numbers: No high altitude lines were flown in this survey
37
38survey_equipment:
39 aircraft: Eurocopter Astar 350 B3
40 magnetometer: Geometrics G822A, Kroum KMAG4 counter
41 magnetometer_installation: front of EM transmitter frame
42 electromagnetic_system: SkyTEM304M time‑domain electromagnetic system
43 electromagnetic_installation: Rigid transmitter frame 40 m beneath helicopter, Receiver coils at rear of transmitter frame 2 m vertical offset
44 spectrometer_system: n/a
45 radar_altimeter_system: n/a
46 laser_altimeter_system: Two MDL ILM 300R laser altimeters
47 laser_altimeter_installation: both laser altimeters were mounted on the suspended EM frame
48 laser_altimeter_information: Sampling rate 30 Hz, range 0.2 to 200 m and uncertainty 10 to 30 cm.
49 inclinometer_system: Two Bjerre Technology inclinometers
50 inclinometer_installation: Both inclinometers were installed on the rear of the EM frame near the Z coil.
51 inclinometer_information: Measured x and y tilt at 2 Hz sampling rate. The x angle is parallel to flight direction (positive when front is above horizontal); y angle is perpendicular to flight direction (negative when right side is above horizontal).
52 gsp_system: Differential GPS using OEMV1-L1 chipset and Trimble Bullet III antennas. Two DGPS units recorded by EM system; a third GPS recorded by the magnetic system.
53 gps_installation: DGPS antenna mounted on top of the boom in front of the frame
54 gps_information: Sampling rate 1 Hz. The uncertainty in the xyz-directions is 1 m after processing.
55 acquisition_system: SkyTEM 304M. The EM system records full waveforms and gate stacks; magnetic DAS synchronized with TEM pulses.
56 gps_information: Sampling rate 1 Hz. The uncertainty in the xyz-directions is 1 m after processing.
57 acquisition_system: SkyTEM 304M. The EM system records full waveforms and gate stacks; magnetic DAS synchronized with TEM pulses.
Create a Data Branch
data_container = survey.gs.add_container('data', **dict(content = "raw and processed data",
comment = "<extra info goes here>"))
Attach leaves to the data branch
Raw Data
# Import raw AEM data from CSV-format.
# Define input data file and associated metadata file
d_data1 = join(data_path, 'data//skytem_contractor_data.csv')
d_supp1 = join(data_path, 'data//skytem_contractor_data.yml')
raw_systems = Metadata.read(join(data_path, "data//skytem_system.yml"))
# raw_systems = {"skytem_system" : survey["nominal_system"],
# "magnetic_system" : survey["magnetic_system"]}
# Add the raw AEM data as a tabular dataset,
# pass the EM system from the survey
rd = data_container.gs.add(key='raw_data', data=d_data1,
metadata_file=d_supp1, system=raw_systems)
Processed Data
# Import processed AEM data from CSV-format.
# Define input data file and associated metadata file
d_data2 = join(data_path, 'data//skytem_processed_data.csv')
d_supp2 = join(data_path, 'data//skytem_processed_data.yml')
print(rd['skytem_system'])
<xarray.DataTree 'skytem_system'>
Group: /survey/data/raw_data/skytem_system
Dimensions: (index: 2000,
hm_gate_times: 32,
lm_gate_times: 28,
gate_times: 22, nv: 2,
n_loop_vertices: 8, xyz: 3,
n_transmitter: 2,
transmitter_lm_waveform_time: 21,
transmitter_hm_waveform_time: 36,
n_receiver: 2, n_couplet: 4,
dim_0: 1)
Coordinates:
* gate_times (gate_times) float64 176B 5....
* nv (nv) int64 16B 0 1
* n_loop_vertices (n_loop_vertices) int64 64B ...
* xyz (xyz) int64 24B 0 1 2
* n_transmitter (n_transmitter) int64 16B 0 1
* transmitter_lm_waveform_time (transmitter_lm_waveform_time) float64 168B ...
* transmitter_hm_waveform_time (transmitter_hm_waveform_time) float64 288B ...
* n_receiver (n_receiver) int64 16B 0 1
* n_couplet (n_couplet) int64 32B 0 1 2 3
Inherited coordinates:
* index (index) int32 8kB 0 1 ... 1999
* hm_gate_times (hm_gate_times) float64 256B ...
* lm_gate_times (lm_gate_times) float64 224B ...
Dimensions without coordinates: dim_0
Data variables: (12/35)
gate_times_bnds (gate_times, nv) float64 352B ...
lm_gate_times_bnds (lm_gate_times, nv) float64 448B ...
hm_gate_times_bnds (hm_gate_times, nv) float64 512B ...
n_loop_vertices_bnds (n_loop_vertices, nv) float64 128B ...
xyz_bnds (xyz, nv) float64 48B -0.5 ....
transmitter_label (n_transmitter) <U2 16B 'LM'...
... ...
couplet_data_type (n_couplet) <U4 64B 'dBdt' ....
couplet_gate_times (n_couplet) <U13 208B 'lm_ga...
data_normalized (dim_0) bool 1B True
skytem_skb_gex_available (dim_0) bool 1B True
reference_frame <U26 104B 'right-handed posi...
coil_orientations <U4 16B 'X, Z'
Attributes:
type: system
mode: airborne
method: electromagnetic, time domain
instrument: SkyTEM 304M
name: skytem_system
Example of how systems can be selected and modified to accurately match the processed data
proc_systems = {"skytem_system" : rd["skytem_system"].isel(lm_gate_times=np.s_[1:],
hm_gate_times=np.s_[10:]),
"magnetic_system" : rd["magnetic_system"]}
Add the processed AEM data as a tabular dataset, passing the updated systems
pd = data_container.gs.add(key='processed_data', data=d_data2,
metadata_file=d_supp2, system=proc_systems)
Processed Data YAML file
1dataset_attrs:
2 content: processed data
3 comment: This dataset includes processed AEM data produced by USGS
4 type: data
5 structure: tabular
6 mode: airborne
7 method: electromagnetic, time domain
8 instrument: SkyTEM 304M
9
10coordinates:
11 x: E_N83WTM
12 y: N_N83WTM
13 z: ELEVATION
14 t: TIMESTAMP
15
16
17variables:
18 pINDEX:
19 standard_name: processing_index
20 long_name: Unique index number for processing
21 units: not_defined
22 missing_value: not_defined
23
24 sLINE_NO:
25 standard_name: master_line
26 long_name: Master line number
27 units: not_defined
28 missing_value: not_defined
29
30 E_N83WTM:
31 standard_name: easting_nad83
32 long_name: Easting, Wisconsin Transverse Mercator (WTM), North American Datum of 1983 (NAD83)
33 units: meter
34 missing_value: not_defined
35
36 N_N83WTM:
37 standard_name: northing_nad83
38 long_name: Northing, Wisconsin Transverse Mercator (WTM), North American Datum of 1983 (NAD83)
39 units: meter
40 missing_value: not_defined
41
42 TIMESTAMP:
43 standard_name: timestamp
44 long_name: Time, decimal days
45 units: days since 1900-1-1 0:0:0
46 calendar: gregorian
47 missing_value: not_defined
48 datum: January 1, 1900
49
50 RECORD:
51 standard_name: record
52 long_name: Workbench record number
53 units: not_defined
54 missing_value: not_defined
55
56 ELEVATION:
57 standard_name: elevation
58 long_name: Digital elevation model
59 units: meter
60 missing_value: not_defined
61 positive: up
62 datum: North American Vertical Datum of 1988 (NAVD88)
63
64 ALT:
65 standard_name: altitude
66 long_name: DGPS instrument altitude
67 units: meter
68 missing_value: not_defined
69
70 NUMDATA:
71 standard_name: number_of_data
72 long_name: Number of active time gates
73 units: not_defined
74 missing_value: not_defined
75
76 LM_Data:
77 standard_name: em_data_lmz
78 long_name: EM data, low moment z-component
79 units: picoVolt per Ampere per meter^4
80 missing_value: -9999.99
81 system_couplet: lm_z
82 dimensions: [index, lm_gate_times]
83
84 LM_DataSTD:
85 standard_name: em_data_error_lmz
86 long_name: EM data error standard deviation, low moment z-component
87 units: picoVolt per Ampere per meter^4
88 missing_value: -9999.99
89 system_couplet: lm_z
90 dimensions: [index, lm_gate_times]
91
92 HM_Data:
93 standard_name: em_data_hmz
94 long_name: EM data, high moment z-component
95 units: picoVolt per Ampere per meter^4
96 missing_value: -9999.99
97 system_couplet: hm_z
98 dimensions: [index, hm_gate_times]
99
100 HM_DataSTD:
101 standard_name: em_data_error_hmz
102 long_name: EM data error standard deviation, high moment z-component
103 units: picoVolt per Ampere per meter^4
104 missing_value: -9999.99
105 system_couplet: hm_z
106 dimensions: [index, hm_gate_times]
107
108 TX_ALTITUDE:
109 standard_name: transmitter_altitude
110 long_name: Processed transmitter altitude
111 units: meter
112 missing_value: -9999.99
113
114 TX_ALTITUDE_STD:
115 standard_name: transmitter_altitude_error
116 long_name: Standard deviation for transmitter altitude
117 units: meter
118 missing_value: -9999.99
119
120 RX_ALTITUDE:
121 standard_name: receiver_altitude
122 long_name: Processed receiver altitude
123 units: meter
124 missing_value: -9999.99
125
126 RX_ALTITUDE_STD:
127 standard_name: receiver_altitude_error
128 long_name: Standard deviation for receiver altitude
129 units: meter
130 missing_value: -9999.99
131
132 txrx_dx:
133 standard_name: txrx_dx
134 long_name: Nominal inline transmitter-receiver offset
135 units: meter
136 missing_value: -9999.99
137
138 txrx_dy:
139 standard_name: txrx_dy
140 long_name: Nominal transverse transmitter-receiver offset
141 units: meter
142 missing_value: -9999.99
143
144 txrx_dz:
145 standard_name: txrx_dz
146 long_name: Calculated vertical transmitter-receiver offset
147 units: meter
148 missing_value: -9999.99
149
150 LINE_NO:
151 standard_name: line_number
152 long_name: Line number
153 units: not_defined
154 missing_value: not_defined
Create a Models Branch
# Create a new container for models
model_container = survey.gs.add_container('models', **dict(content = "Inverted models",
comment = "This is a test"))
Inverted Models
# Import inverted AEM models from CSV-format.
# Define input data file and associated metadata file
m_data3 = join(data_path, 'model//skytem_inverted_models.csv')
m_supp3 = join(data_path, 'model//skytem_inverted_models.yml')
# Add the inverted AEM models as a tabular dataset
mods = model_container.gs.add(key='inverted_models', data=m_data3,
metadata_file=m_supp3)
Inverted Models YAML file
1dataset_attrs:
2 content: inverted resistivity models
3 comment: This dataset includes inverted resistivity models derived from processed AEM data produced by USGS
4 type: model
5 structure: tabular
6 mode: airborne
7 method: electromagnetic, time domain
8 instrument: SkyTEM 304M
9 property: electrical resistivity
10
11inversion_parameters:
12 dataset_attrs:
13 type: parameters
14 method: electromagnetic, time domain
15 instrument: SkyTEM 304M
16 mode: airborne
17 property: electrical resistivity
18
19 variables:
20 model_file: WI_SkyTEM_2021_InvertedModels.csv
21 inversion_software: Aarhus Workbench
22 software_version: "v 6.0.1.0"
23 software_reference: "placeholder"
24 date: "03/2021 - 06/2021"
25 description: "Processed data were inverted in Aarhus Workbench software (v 6.0.1.0) using laterally constrained inversion to recover 40-layer fixed depth blocky resistivity models by USGS 03/2021 - 06/2021; Inverted resistivity models were exported to netCDF 11/2021."
26 data_file: WI_SkyTEM_2021_ProcessedData.csv
27
28coordinates:
29 x: E_N83WTM
30 y: N_N83WTM
31 z: ELEVATION
32 t: TIMESTAMP
33
34dimensions:
35 layer_depth:
36 standard_name: layer_depth
37 long_name: Depth to model layer
38 units: meters
39 missing_value: not_defined
40 centers: [0.375, 1.16 , 2.02 ,
41 2.965, 4.005, 5.145,
42 6.39 , 7.755, 9.255,
43 10.9 , 12.7 , 14.675,
44 16.845, 19.22 , 21.825,
45 24.685, 27.815, 31.25 ,
46 35.02 , 39.15 , 43.68 ,
47 48.65 , 54.095, 60.065,
48 66.615, 73.795, 81.67 ,
49 90.31 , 99.78 , 110.16 ,
50 121.545, 134.03 , 147.72 ,
51 162.73 , 179.19 , 197.24 ,
52 217.035, 238.745, 262.55 , 343.75]
53 bounds: [[ 0.0 , 0.75],
54 [ 0.75, 1.57],
55 [ 1.57, 2.47],
56 [ 2.47, 3.46],
57 [ 3.46, 4.55],
58 [ 4.55, 5.74],
59 [ 5.74, 7.04],
60 [ 7.04, 8.47],
61 [ 8.47, 10.04],
62 [ 10.04, 11.76],
63 [ 11.76, 13.64],
64 [ 13.64, 15.71],
65 [ 15.71, 17.98],
66 [ 17.98, 20.46],
67 [ 20.46, 23.19],
68 [ 23.19, 26.18],
69 [ 26.18, 29.45],
70 [ 29.45, 33.05],
71 [ 33.05, 36.99],
72 [ 36.99, 41.31],
73 [ 41.31, 46.05],
74 [ 46.05, 51.25],
75 [ 51.25, 56.94],
76 [ 56.94, 63.19],
77 [ 63.19, 70.04],
78 [ 70.04, 77.55],
79 [ 77.55, 85.79],
80 [ 85.79, 94.83],
81 [ 94.83, 104.73],
82 [104.73, 115.59],
83 [115.59, 127.5 ],
84 [127.5 , 140.56],
85 [140.56, 154.88],
86 [154.88, 170.58],
87 [170.58, 187.8 ],
88 [187.8 , 206.68],
89 [206.68, 227.39],
90 [227.39, 250.1 ],
91 [250.1 , 275.0 ],
92 [275.0 , 412.5 ]]
93
94variables:
95 pINDEX:
96 standard_name: processing_index
97 long_name: Unique index number for processing
98 units: not_defined
99 missing_value: not_defined
100
101 sLINE_NO:
102 standard_name: master_line
103 long_name: Master line number
104 units: not_defined
105 missing_value: not_defined
106
107 E_N83WTM:
108 standard_name: easting_nad83
109 long_name: Easting, Wisconsin Transverse Mercator (WTM), North American Datum of 1983 (NAD83)
110 units: meter
111 missing_value: not_defined
112 axis: x
113
114 N_N83WTM:
115 standard_name: northing_nad83
116 long_name: Northing, Wisconsin Transverse Mercator (WTM), North American Datum of 1983 (NAD83)
117 units: meter
118 missing_value: not_defined
119 axis: y
120
121 TIMESTAMP:
122 standard_name: timestamp
123 long_name: Time, decimal days since January 1, 1900
124 units: day
125 missing_value: not_defined
126 axis: t
127 datum: January 1, 1900
128
129 RECORD:
130 standard_name: record
131 long_name: Workbench record number
132 units: not_defined
133 missing_value: not_defined
134
135 ELEVATION:
136 standard_name: elevation
137 long_name: Digital elevation model
138 units: meter
139 missing_value: not_defined
140 axis: z
141 positive: up
142 datum: North American Vertical Datum of 1988 (NAVD88)
143
144 ALT:
145 standard_name: altitude
146 long_name: DGPS instrument altitude
147 units: meter
148 missing_value: not_defined
149
150 INVALT:
151 standard_name: inverted_altitude
152 long_name: Inverted instrument altitude
153 units: meter
154 missing_value: not_defined
155
156 INVALTSTD:
157 standard_name: inverted_altitude_uncertainty
158 long_name: Standard deviation of inverted instrument altitude
159 units: meter
160 missing_value: not_defined
161
162 DELTAALT:
163 standard_name: inverted_altitude_difference
164 long_name: Measured minus inverted altitude
165 units: meter
166 missing_value: not_defined
167
168 NUMDATA:
169 standard_name: number_of_data
170 long_name: Number of active time gates
171 units: not_defined
172 missing_value: not_defined
173
174 RESDATA:
175 standard_name: data_residual
176 long_name: Error-weighted inversion data misfit (target = 1.0)
177 units: not_defined
178 missing_value: not_defined
179
180 RESTOTAL:
181 standard_name: total_residual
182 long_name: Total inversion residual (data and model regularization)
183 units: not_defined
184 missing_value: not_defined
185
186 RHO_I:
187 standard_name: layer_resistivity
188 long_name: Inverted layer resistivity
189 units: Ohm*meter
190 missing_value: not_defined
191 dimensions: [index, layer_depth]
192
193 RHO_I_STD:
194 standard_name: layer_resistivity_uncertainty
195 long_name: Uncertainty in inverted layer resistivity
196 units: not_defined
197 missing_value: not_defined
198 dimensions: [index, layer_depth]
199
200 DOI_CONSERVATIVE:
201 standard_name: depth_of_investigation_conservative
202 long_name: Conservative estimate of depth of investigation (DOI)
203 units: meter
204 missing_value: not_defined
205
206 DOI_STANDARD:
207 standard_name: depth_of_investigation_standard
208 long_name: Standard estimate of depth of investigation (DOI)
209 units: meter
210 missing_value: not_defined
211
212 DEP_TOP:
213 standard_name: depth_top
214 long_name: Top of model layers
215 units: meter
216 missing_value: not_defined
217 dimensions: [index, layer_depth]
218
219 DEP_BOT:
220 standard_name: depth_bottom
221 long_name: Bottom of model layers
222 units: meter
223 missing_value: not_defined
224 dimensions: [index, layer_depth]
225
226 LINE_NO:
227 standard_name: line_number
228 long_name: Line number
229 units: not_defined
230 missing_value: not_defined
Derivative Products
Bedrock Picks
Adding bedrock picks to the ‘data’ branch
# Import AEM-based estimated of depth to bedrock from CSV-format.
# Define input data file and associated metadata file
d_data4 = join(data_path, 'data//top_dolomite_blocky_lidar.csv')
d_supp4 = join(data_path, 'data//bedrock_picks.yml')
# Add the AEM-based estimated of depth to bedrock as a tabular dataset
bedrock = data_container.gs.add(key='depth_to_bedrock', data=d_data4,
metadata_file=d_supp4)
Bedrock Picks YAML file
1dataset_attrs:
2 content: bedrock elevation points
3 comment: This dataset includes AEM-derived point estimates of the elevation of the top of bedrock produced by USGS
4 type: data
5 structure: tabular
6 mode: airborne
7 method: electromagnetic, time domain
8 instrument: SkyTEM 304M
9
10coordinates:
11 x: E_N83WTM
12 y: N_N83WTM
13
14variables:
15 ID:
16 standard_name: identifier
17 long_name: Unique identifier
18 units: not_defined
19 missing_value: not_defined
20
21 E_N83WTM:
22 standard_name: easting
23 long_name: Easting, Wisconsin Transverse Mercator (WTM), North American Datum of 1983 (NAD83)
24 units: meter
25 missing_value: not_defined
26
27 N_N83WTM:
28 standard_name: northing
29 long_name: Northing, Wisconsin Transverse Mercator (WTM), North American Datum of 1983 (NAD83)
30 units: meter
31 missing_value: not_defined
32
33 BR_ELEVATION:
34 standard_name: top_bedrock_elevation
35 long_name: Elevation, top of dolomite bedrock, North American Vertical Datum of 1988 (NAVD88)
36 units: meter
37 missing_value: not_defined
38
39 ZSTD:
40 standard_name: elevation_uncertainty
41 long_name: Standard devation of top bedrock elevation
42 units: meter
43 missing_value: not_defined
44
45 OriginType:
46 standard_name: point_origin
47 long_name: Point origin type; 3 = automated pick from resistivity value; 0 = manual pick
48 units: not_defined
49 missing_value: not_defined
50
51 EditDate:
52 standard_name: edit_date
53 long_name: Date of interpretation point
54 units: not_defined
55 format: M/D/YYYY H:MM
56 missing_value: not_defined
57 dtype: str
Raster Maps
Create a 3rd container for the derived raaster maps
derived_maps = survey.gs.add_container('derived_maps', **dict(content = "raster products derived from airborne data and models"))
# Import interpolated bedrock and magnetic maps from TIF-format.
# Define input metadata file (which contains the TIF filenames linked to variable names)
m_supp5 = join(data_path, 'data//magnetics_bedrock_picks.yml')
# Add the interpolated maps as a raster dataset
maps = derived_maps.gs.add(key='maps', metadata_file=m_supp5)
Gridded Maps YAML file
1dataset_attrs :
2 content : gridded magnetic and bedrock maps
3 comment: This dataset includes AEM-derived estimates of the elevation of the top of bedrock produced by USGS
4 type: data
5 structure: raster
6 mode: airborne
7 method: ["electromagnetic, time domain", "magnetic, total field"]
8 instrument: skytem
9 property: [total magnetic intensity, depth to bedrock]
10
11coordinates:
12 x: E_Nad83
13 y: N_Nad83
14
15dimensions:
16 x: E_Nad83
17 y: N_Nad83
18
19variables:
20 magnetic_tmi:
21 standard_name: total_magnetic_intensity
22 long_name: Total magnetic intensity, diurnally corrected and filtered
23 units: nanoTesla
24 missing_value: -9999.99
25 files : [mag_tmi.tif]
26 dimensions: [x, y]
27
28 magnetic_rmf:
29 standard_name: residual_magnetic_field
30 long_name: Residual magnetic field, IGRF corrected from 2015 model
31 units: nanoTesla
32 missing_value: -9999.99
33 files : [mag_rmf.tif]
34 dimensions: [x, y]
35
36 bedrock_top_elevation:
37 standard_name: bedrock_top_elevation
38 long_name: Elevation, top of dolomite bedrock, North American Vertical Datum of 1988 (NAVD88)
39 units: foot
40 missing_value: -9999.99
41 files : [top_bedrock.tif]
42 dimensions: [x, y]
43
44 bedrock_depth:
45 standard_name: bedrock_depth
46 long_name: Depth to bedrock
47 units: foot
48 missing_value: -9999.9
49 files : [bedrock_depth.tif]
50 dimensions: [x, y]
51
52 E_Nad83:
53 standard_name: easting_nad83
54 long_name: Easting, Wisconsin Transverse Mercator (WTM), North American Datum of 1983 (NAD83)
55 units: meter
56 missing_value: not_defined
57 axis : x
58
59 N_Nad83:
60 standard_name: northing_nad83
61 long_name: Northing, Wisconsin Transverse Mercator (WTM), North American Datum of 1983 (NAD83)
62 units: meter
63 missing_value: not_defined
64 axis : y
Save to NetCDF file
d_out = join(data_path, 'skytem.nc')
survey.gs.to_netcdf(d_out)
Export just one branch to file
The gspy goal is to have the complete survey in a single file. However, we can also save containers or datasets separately.
data_container.gs.to_netcdf(join(data_path, 'test_datacontainer.nc'))
Opening a GS NetCDF
new_survey = gspy.open_datatree(d_out)['survey']
View the Data Tree
print(new_survey.gs.tree)
/survey
/survey/data
/survey/models
/survey/derived_maps
/survey/data/raw_data
/survey/data/processed_data
/survey/data/depth_to_bedrock
/survey/models/inverted_models
/survey/derived_maps/maps
/survey/data/raw_data/skytem_system
/survey/data/raw_data/magnetic_system
/survey/data/processed_data/skytem_system
/survey/data/processed_data/magnetic_system
/survey/models/inverted_models/inversion_parameters
print(new_survey)
<xarray.DataTree 'survey'>
Group: /survey
│ Dimensions: ()
│ Coordinates:
│ spatial_ref float64 8B ...
│ Data variables:
│ survey_information float64 8B ...
│ flightline_information float64 8B ...
│ survey_equipment float64 8B ...
│ Attributes:
│ type: survey
│ title: SkyTEM Airborne Electromagnetic (AEM) Survey, Northeast Wi...
│ institution: USGS Geology, Geophysics, and Geochemistry Science Center
│ source: SkyTEM raw data, USGS processed data and inverted resistiv...
│ history: (1) Data acquisition 01/2021 - 02/2021 by SkyTEM Canada In...
│ references: Minsley, Burke J., B.R. Bloss, D.J. Hart, W. Fitzpatrick, ...
│ comment: This dataset includes minimally processed (raw) AEM and ra...
│ summary: Airborne electromagnetic (AEM) and magnetic survey data we...
│ content: Wisconsin SkyTEM survey information /survey; raw and proce...
│ gspy_version: 2.2.8
│ conventions: GS-2.0, CF-1.13
├── Group: /survey/data
│ │ Dimensions: ()
│ │ Data variables:
│ │ spatial_ref float64 8B ...
│ │ Attributes:
│ │ content: raw and processed data
│ │ comment: <extra info goes here>
│ │ type: container
│ ├── Group: /survey/data/raw_data
│ │ │ Dimensions: (index: 2000, hm_gate_times: 32, lm_gate_times: 28)
│ │ │ Coordinates:
│ │ │ * index (index) float64 16kB 0.0 1.0 2.0 ... 1.998e+03 1.999e+03
│ │ │ * hm_gate_times (hm_gate_times) float64 256B 2.886e-05 ... 0.003544
│ │ │ * lm_gate_times (lm_gate_times) float64 224B -1.135e-06 ... 0.001394
│ │ │ spatial_ref float64 8B ...
│ │ │ x (index) float64 16kB ...
│ │ │ y (index) float64 16kB ...
│ │ │ z (index) float64 16kB ...
│ │ │ t (index) float64 16kB ...
│ │ │ Data variables: (12/30)
│ │ │ _60hz_intensity (index) float64 16kB ...
│ │ │ alt (index) float64 16kB ...
│ │ │ anglex (index) float64 16kB ...
│ │ │ angley (index) float64 16kB ...
│ │ │ base_mag (index) float64 16kB ...
│ │ │ curr_hm (index) float64 16kB ...
│ │ │ ... ...
│ │ │ mag_filt (index) float64 16kB ...
│ │ │ mag_raw (index) float64 16kB ...
│ │ │ n_wgs84 (index) float64 16kB ...
│ │ │ rmf (index) float64 16kB ...
│ │ │ time (index) object 16kB ...
│ │ │ tmi (index) float64 16kB ...
│ │ │ Attributes:
│ │ │ content: raw data
│ │ │ comment: This dataset includes minimally processed (raw) AEM and raw/...
│ │ │ type: data
│ │ │ structure: tabular
│ │ │ mode: airborne
│ │ │ method: ['electromagnetic, time domain', 'magnetic']
│ │ │ instrument: SkyTEM 304M
│ │ ├── Group: /survey/data/raw_data/skytem_system
│ │ │ Dimensions: (gate_times: 22, nv: 2,
│ │ │ lm_gate_times: 28,
│ │ │ hm_gate_times: 32,
│ │ │ n_loop_vertices: 8, xyz: 3,
│ │ │ n_transmitter: 2,
│ │ │ transmitter_lm_waveform_time: 21,
│ │ │ transmitter_hm_waveform_time: 36,
│ │ │ n_receiver: 2, n_couplet: 4,
│ │ │ dim_0: 1)
│ │ │ Coordinates:
│ │ │ * gate_times (gate_times) float64 176B 5....
│ │ │ * nv (nv) float64 16B 0.0 1.0
│ │ │ * n_loop_vertices (n_loop_vertices) float64 64B ...
│ │ │ * xyz (xyz) float64 24B 0.0 1.0 2.0
│ │ │ * n_transmitter (n_transmitter) float64 16B ...
│ │ │ * transmitter_lm_waveform_time (transmitter_lm_waveform_time) float64 168B ...
│ │ │ * transmitter_hm_waveform_time (transmitter_hm_waveform_time) float64 288B ...
│ │ │ * n_receiver (n_receiver) float64 16B 0.0...
│ │ │ * n_couplet (n_couplet) float64 32B 0.0 ...
│ │ │ Dimensions without coordinates: dim_0
│ │ │ Data variables: (12/35)
│ │ │ gate_times_bnds (gate_times, nv) float64 352B ...
│ │ │ lm_gate_times_bnds (lm_gate_times, nv) float64 448B ...
│ │ │ hm_gate_times_bnds (hm_gate_times, nv) float64 512B ...
│ │ │ n_loop_vertices_bnds (n_loop_vertices, nv) float64 128B ...
│ │ │ xyz_bnds (xyz, nv) float64 48B ...
│ │ │ transmitter_label (n_transmitter) <U2 16B ...
│ │ │ ... ...
│ │ │ couplet_data_type (n_couplet) <U4 64B ...
│ │ │ couplet_gate_times (n_couplet) <U13 208B ...
│ │ │ data_normalized (dim_0) bool 1B ...
│ │ │ skytem_skb_gex_available (dim_0) bool 1B ...
│ │ │ reference_frame <U26 104B ...
│ │ │ coil_orientations <U4 16B ...
│ │ │ Attributes:
│ │ │ type: system
│ │ │ mode: airborne
│ │ │ method: electromagnetic, time domain
│ │ │ instrument: SkyTEM 304M
│ │ │ name: skytem_system
│ │ └── Group: /survey/data/raw_data/magnetic_system
│ │ Dimensions: (n_transmitter: 1, n_receiver: 1,
│ │ n_couplet: 1, n_base_magnetometer: 1,
│ │ base_mag_locations: 2)
│ │ Coordinates:
│ │ * n_transmitter (n_transmitter) float64 8B 0.0
│ │ * n_receiver (n_receiver) float64 8B 0.0
│ │ * n_couplet (n_couplet) float64 8B 0.0
│ │ * n_base_magnetometer (n_base_magnetometer) float64 8B 0.0
│ │ * base_mag_locations (base_mag_locations) float64 16B 1.0 2.0
│ │ Data variables: (12/28)
│ │ transmitter_label (n_transmitter) <U7 28B ...
│ │ transmitter_description (n_transmitter) <U142 568B ...
│ │ receiver_label (n_receiver) <U19 76B ...
│ │ receiver_sensor_type (n_receiver) <U23 92B ...
│ │ receiver_sensor_model (n_receiver) <U6 24B ...
│ │ receiver_sensor_manufacturer (n_receiver) <U10 40B ...
│ │ ... ...
│ │ diurnal_correction <U110 440B ...
│ │ tieline_levelling <U33 132B ...
│ │ microlevelling <U31 124B ...
│ │ igrf_model_date <U21 84B ...
│ │ igrf_model_location <U54 216B ...
│ │ igrf_model_height <U69 276B ...
│ │ Attributes:
│ │ type: system
│ │ mode: airborne
│ │ method: magnetic
│ │ instrument: Geometrics G-822A cesium‑vapor magnetometer
│ │ name: magnetic_system
│ ├── Group: /survey/data/processed_data
│ │ │ Dimensions: (index: 2000, lm_gate_times: 27, hm_gate_times: 22)
│ │ │ Coordinates:
│ │ │ * index (index) float64 16kB 0.0 1.0 2.0 ... 1.998e+03 1.999e+03
│ │ │ * lm_gate_times (lm_gate_times) float64 216B 3.65e-07 ... 0.001394
│ │ │ * hm_gate_times (hm_gate_times) float64 176B 5.636e-05 ... 0.003544
│ │ │ spatial_ref float64 8B ...
│ │ │ x (index) float64 16kB ...
│ │ │ y (index) float64 16kB ...
│ │ │ z (index) float64 16kB ...
│ │ │ t (index) float64 16kB ...
│ │ │ Data variables: (12/17)
│ │ │ pindex (index) float64 16kB ...
│ │ │ sline_no (index) float64 16kB ...
│ │ │ record (index) float64 16kB ...
│ │ │ alt (index) float64 16kB ...
│ │ │ numdata (index) float64 16kB ...
│ │ │ lm_data (index, lm_gate_times) float64 432kB ...
│ │ │ ... ...
│ │ │ rx_altitude (index) float64 16kB ...
│ │ │ rx_altitude_std (index) float64 16kB ...
│ │ │ txrx_dx (index) float64 16kB ...
│ │ │ txrx_dy (index) float64 16kB ...
│ │ │ txrx_dz (index) float64 16kB ...
│ │ │ line_no (index) float64 16kB ...
│ │ │ Attributes:
│ │ │ content: processed data
│ │ │ comment: This dataset includes processed AEM data produced by USGS
│ │ │ type: data
│ │ │ structure: tabular
│ │ │ mode: airborne
│ │ │ method: electromagnetic, time domain
│ │ │ instrument: SkyTEM 304M
│ │ ├── Group: /survey/data/processed_data/skytem_system
│ │ │ Dimensions: (gate_times: 22, nv: 2,
│ │ │ lm_gate_times: 27,
│ │ │ hm_gate_times: 22,
│ │ │ n_loop_vertices: 8, xyz: 3,
│ │ │ n_transmitter: 2,
│ │ │ transmitter_lm_waveform_time: 21,
│ │ │ transmitter_hm_waveform_time: 36,
│ │ │ n_receiver: 2, n_couplet: 4,
│ │ │ dim_0: 1)
│ │ │ Coordinates:
│ │ │ * gate_times (gate_times) float64 176B 5....
│ │ │ * nv (nv) float64 16B 0.0 1.0
│ │ │ * n_loop_vertices (n_loop_vertices) float64 64B ...
│ │ │ * xyz (xyz) float64 24B 0.0 1.0 2.0
│ │ │ * n_transmitter (n_transmitter) float64 16B ...
│ │ │ * transmitter_lm_waveform_time (transmitter_lm_waveform_time) float64 168B ...
│ │ │ * transmitter_hm_waveform_time (transmitter_hm_waveform_time) float64 288B ...
│ │ │ * n_receiver (n_receiver) float64 16B 0.0...
│ │ │ * n_couplet (n_couplet) float64 32B 0.0 ...
│ │ │ Dimensions without coordinates: dim_0
│ │ │ Data variables: (12/35)
│ │ │ gate_times_bnds (gate_times, nv) float64 352B ...
│ │ │ lm_gate_times_bnds (lm_gate_times, nv) float64 432B ...
│ │ │ hm_gate_times_bnds (hm_gate_times, nv) float64 352B ...
│ │ │ n_loop_vertices_bnds (n_loop_vertices, nv) float64 128B ...
│ │ │ xyz_bnds (xyz, nv) float64 48B ...
│ │ │ transmitter_label (n_transmitter) <U2 16B ...
│ │ │ ... ...
│ │ │ couplet_data_type (n_couplet) <U4 64B ...
│ │ │ couplet_gate_times (n_couplet) <U13 208B ...
│ │ │ data_normalized (dim_0) bool 1B ...
│ │ │ skytem_skb_gex_available (dim_0) bool 1B ...
│ │ │ reference_frame <U26 104B ...
│ │ │ coil_orientations <U4 16B ...
│ │ │ Attributes:
│ │ │ type: system
│ │ │ mode: airborne
│ │ │ method: electromagnetic, time domain
│ │ │ instrument: SkyTEM 304M
│ │ │ name: skytem_system
│ │ └── Group: /survey/data/processed_data/magnetic_system
│ │ Dimensions: (n_transmitter: 1, n_receiver: 1,
│ │ n_couplet: 1, n_base_magnetometer: 1,
│ │ base_mag_locations: 2)
│ │ Coordinates:
│ │ * n_transmitter (n_transmitter) float64 8B 0.0
│ │ * n_receiver (n_receiver) float64 8B 0.0
│ │ * n_couplet (n_couplet) float64 8B 0.0
│ │ * n_base_magnetometer (n_base_magnetometer) float64 8B 0.0
│ │ * base_mag_locations (base_mag_locations) float64 16B 1.0 2.0
│ │ Data variables: (12/28)
│ │ transmitter_label (n_transmitter) <U7 28B ...
│ │ transmitter_description (n_transmitter) <U142 568B ...
│ │ receiver_label (n_receiver) <U19 76B ...
│ │ receiver_sensor_type (n_receiver) <U23 92B ...
│ │ receiver_sensor_model (n_receiver) <U6 24B ...
│ │ receiver_sensor_manufacturer (n_receiver) <U10 40B ...
│ │ ... ...
│ │ diurnal_correction <U110 440B ...
│ │ tieline_levelling <U33 132B ...
│ │ microlevelling <U31 124B ...
│ │ igrf_model_date <U21 84B ...
│ │ igrf_model_location <U54 216B ...
│ │ igrf_model_height <U69 276B ...
│ │ Attributes:
│ │ type: system
│ │ mode: airborne
│ │ method: magnetic
│ │ instrument: Geometrics G-822A cesium‑vapor magnetometer
│ │ name: magnetic_system
│ └── Group: /survey/data/depth_to_bedrock
│ Dimensions: (index: 82864)
│ Coordinates:
│ * index (index) float64 663kB 0.0 1.0 2.0 ... 8.286e+04 8.286e+04
│ spatial_ref float64 8B ...
│ x (index) float64 663kB ...
│ y (index) float64 663kB ...
│ Data variables:
│ id (index) float64 663kB ...
│ br_elevation (index) float64 663kB ...
│ zstd (index) float64 663kB ...
│ origintype (index) float64 663kB ...
│ editdate (index) object 663kB ...
│ Attributes:
│ content: bedrock elevation points
│ comment: This dataset includes AEM-derived point estimates of the ele...
│ type: data
│ structure: tabular
│ mode: airborne
│ method: electromagnetic, time domain
│ instrument: SkyTEM 304M
├── Group: /survey/models
│ │ Dimensions: ()
│ │ Data variables:
│ │ spatial_ref float64 8B ...
│ │ Attributes:
│ │ content: Inverted models
│ │ comment: This is a test
│ │ type: container
│ └── Group: /survey/models/inverted_models
│ │ Dimensions: (layer_depth: 40, nv: 2, index: 2000)
│ │ Coordinates:
│ │ * layer_depth (layer_depth) float64 320B 0.375 1.16 2.02 ... 262.6 343.8
│ │ * nv (nv) float64 16B 0.0 1.0
│ │ * index (index) float64 16kB 0.0 1.0 2.0 ... 1.998e+03 1.999e+03
│ │ spatial_ref float64 8B ...
│ │ x (index) float64 16kB ...
│ │ y (index) float64 16kB ...
│ │ z (index) float64 16kB ...
│ │ t (index) float64 16kB ...
│ │ Data variables: (12/18)
│ │ layer_depth_bnds (layer_depth, nv) float64 640B ...
│ │ pindex (index) float64 16kB ...
│ │ sline_no (index) float64 16kB ...
│ │ record (index) float64 16kB ...
│ │ alt (index) float64 16kB ...
│ │ invalt (index) float64 16kB ...
│ │ ... ...
│ │ rho_i_std (index, layer_depth) float64 640kB ...
│ │ dep_top (index, layer_depth) float64 640kB ...
│ │ dep_bot (index, layer_depth) float64 640kB ...
│ │ doi_conservative (index) float64 16kB ...
│ │ doi_standard (index) float64 16kB ...
│ │ line_no (index) float64 16kB ...
│ │ Attributes:
│ │ content: inverted resistivity models
│ │ comment: This dataset includes inverted resistivity models derived fr...
│ │ type: model
│ │ structure: tabular
│ │ mode: airborne
│ │ method: electromagnetic, time domain
│ │ instrument: SkyTEM 304M
│ │ property: electrical resistivity
│ └── Group: /survey/models/inverted_models/inversion_parameters
│ Dimensions: ()
│ Data variables:
│ model_file <U33 132B ...
│ inversion_software <U16 64B ...
│ software_version <U9 36B ...
│ software_reference <U11 44B ...
│ date <U17 68B ...
│ description <U253 1kB ...
│ data_file <U32 128B ...
│ Attributes:
│ type: parameters
│ method: electromagnetic, time domain
│ instrument: SkyTEM 304M
│ mode: airborne
│ property: electrical resistivity
│ name: inversion_parameters
└── Group: /survey/derived_maps
│ Dimensions: ()
│ Data variables:
│ spatial_ref float64 8B ...
│ Attributes:
│ content: raster products derived from airborne data and models
│ type: container
└── Group: /survey/derived_maps/maps
Dimensions: (x: 799, nv: 2, y: 1155)
Coordinates:
* x (x) float64 6kB 6.551e+05 6.552e+05 ... 7.349e+05
* nv (nv) float64 16B 0.0 1.0
* y (y) float64 9kB 4.953e+05 4.952e+05 ... 3.799e+05
spatial_ref float64 8B ...
Data variables:
x_bnds (x, nv) float64 13kB ...
y_bnds (y, nv) float64 18kB ...
magnetic_tmi (y, x) float64 7MB ...
magnetic_rmf (y, x) float64 7MB ...
bedrock_top_elevation (y, x) float32 4MB ...
bedrock_depth (y, x) float32 4MB ...
Attributes:
content: gridded magnetic and bedrock maps
comment: This dataset includes AEM-derived estimates of the elevation...
type: data
structure: raster
mode: airborne
method: ['electromagnetic, time domain', 'magnetic, total field']
instrument: skytem
property: ['total magnetic intensity', 'depth to bedrock']
Plotting Examples
plt.figure()
new_survey['data']['raw_data']['height'].plot()
plt.tight_layout()
pcd = new_survey['data']['processed_data']
plt.figure()
pcd['tx_altitude'].plot()
plt.tight_layout()
m = new_survey['derived_maps']['maps']
plt.figure()
m['magnetic_tmi'].plot(cmap='jet')
plt.tight_layout()
plt.show()
Total running time of the script: (0 minutes 1.565 seconds)