In-Class Exercise 10 [Geophysics]
In this exercise we’ll bring in data from two geophysical instruments with little pre-processing, georeference the data as 3d points and do some visualization of anomalies comparing the Magnetometer and Ground Penetrating Radar (GPR) datasets with a 2011 orthophoto as well as aerial LiDAR topographic data.
Part II. Interpolation and Visualization
Introduction
This exercise was developed from data collected at Fort Ross an early 19th century Russian merchant colony established on Kashaya Pomo land in what is now referred to as the Sonoma County coast. In the area west of the reconstructed Russian stockade an open field separates the Stockade from a later American Period house and ranch buildings known as Call Ranch.
UC Berkeley 2014 Geophysical Survey Blocks at Fort Ross
Professor Kent Lightfoot coordinated a study in June 2014 that included a number of researchers. The geophysical survey data we are using in this exercise were collected primarily by Rob Cuthrell using a G858 Geometrics cesium magnetometer/gradiometer and Scott Byram using a GSSI SIR 3000 Ground Penetrating Radar with a 400 mhz antenna. The geophysical survey occurred in a series of rectangles or "blocks" with one series of blocks for the magnetometer and another series for the GPR. The corners of each block were shot in with Total Station although differential GPS in this area (~15 cm error) would have provided adequate accuracy. While each instrument is quite different the data is gathered in regular transects that fall inside the polygons that define each block. Both Mag and GPR data is provided in a series of points commonly referred to as a "grid" but these grid points are in a local coordinate system they aren't georeferenced. Prior to analysis and interpretation today we will bring the grid of points into Arcmap as vector points and then we'll relocate them inside of the polygon blocks that provide their overall geographical location. [Note that GPR data is initially acquired as vertical profiles and these data for Grid 7 & 8 have had the initial processing step of being assembled into a horizontal "time-slice" for data from 0 to 8 nano-seconds depth]
The focus in this study area was to attempt to locate possible historical buildings or enclosures that align with the north-south alignment of the survey blocks (and therefore not aligned with the adjacent Call Ranch structures). Fieldwork data usually comes with Field notes and with geophysics projects and it’s especially important to have some field sketches of the placement of the grids. Two approaches are generally taken. Either (1) the study area has an existing local coordinate grid and geophysical transects comply with these or (2) survey plots are laid out and corners are mapped in as best suits the instrumentation and the needs of the fieldwork schedule. A GIS based approach should be able to handle either type.
Here we use the second type and in such cases it is helpful if surveyor create a field diagram showing the alignment of the survey blocks with respect to a few key local features as well as the walking direction of the survey work to ensure that we correctly align the data in Arcmap.
Download Ex 10 Geophysics data Links to an external site.in the Box directory (80mb)
Part I. Data Import and Visualization
A. IMPORT DATA
Steps
1. View the TXT file in a text editor like Notepad++ or Wordpad (often these files are too large for Notepad). Check the header row and how is the file fixed-width or delimited? For importing into ArcGIS comma-delimited is best, tab is OK, Space separated is worst. ArcGIS also reads Excel XLS format so today we'll import into Excel and bring the XLS into Arcmap. Leave the original TXT open for later comparison in Arcmap.
2. Excel:
Open each TXT file in MS Excel: browse to the folder but set the filename file to “All files *.*”
Next in the “Text Import Wizard” use the delimited wizard.
What line do you begin import? Check this by looking at the original TXT file. File Origin: 437 OEM United States text format is fine.
This may take some experimentation. GPR and Mag are slightly different.
Does your data have headers? (Answer: yes in one file, no in another).
Tips: scroll down halfway and check other values. In Step 3 use General data format except for the left-most column which we should choose Skip because it’s empty.
- the “Treat consecutive delimiters as one” in general should be checked (this is really handy).
- The first value of a column will often define the number type for the entire column (hence the “does my data have Headers” question, since Headers are text it throws off the format).
- Even when the “Headers” box is checked if the first value Excel encounters (say in Row 2) is an Integer then Excel may define the entire column as Integer and hide or round subsequent non-integer values. ArcGIS is even worse about this with bringing in CSV files. Watch out for this issue -- you may need to go in with a text editor and add ".0" to the first row of values to force Arcmap to treat that column as a real (non-integer) number. Today since we're just using XLS format we shouldn't have this problem.
- The end result should be data in three columns with a simplified Header (first row) that will import properly into ArcGIS.
Two further steps in Excel:
1. Create a header in Row 1 (no spaces or symbols allowed) such as X, Y, Value and let’s add a column called “Z” that contains a relative elevation value that represents a visualization depth for the dataset. This is not literally the depth of data points but rather a layering method for visualization. Fill the Z column with -0.4 in the Mag data and -0.5 in the GPR data.
To do that, type in the Header “Z” and then populate two of the fields with the depth value (-0.5 or -0.4), then with the next empty depth value cell selected click Data > Fill Down or “Flash Fill” in Excel 2013 so the entire column is populated with that depth.
2. Make sure to save the file to “Excel 97-2003 Workbook (*.xls)” format. Multiple Worksheets are allowed so, for example, multiple slices of GPR could be saved into one XLS as separate worksheets.
Repeat for all the geophysical TXT files.
B. Convert data to Arcmap point features
In both Mag and GPR data the instrument gathers higher precision data in one axis than another. That is, the instrument is collecting a data stream as you travel in parallel collection lines but you can only space your lines every 50 or 100cm. This is why the numbers jump up in 50 cm increments in the X column but the Y column has finer increments.
The other thing you will notice: there are no UTM or other coordinates. These files are currently in a local reference system and need to be moved to their proper position before we can view them together in a GIS.
1. Open a Blank Map in Arcmap.
2. Bring in the newly created excel tables. Scroll down the XLS table in Arcmap and compare it with the original TXT file still open in a simple text editor.
Inspect the data: How many rows in the Arcmap table versus the text file (a good text editor like Notepad++ can show line numbers). Do some spot checking (use Search if necessary to find corresponding values in the two programs). Inspect any unusual formatting or outlier values. Check the values of the final few rows in both files to make sure there was no truncating or misalignment of columns. Any issues at this stage are serious and need to be resolved before proceeding.
3. Right click Table layer and choose Display X, Y data.
Assign the X, Y, and Z values. For Coordinate System choose “NAD 1983 UTM 10N”. This is primarily to get the units defined as Meters (though the resulting points will be somewhere near the equator!).
Allow it to add an Object-ID field.
Click on the “Events” theme layer and choose Export > Export Data … save each to Shapefile in a new “Geophys” folder in the Ex10_Fort Ross folder and name the GPR file "GPR78" and Mag data as "Mag8".
Close map project, don’t save.
C. Georeference Mag vector points
Close Arcmap file and Open blank MXD. Add GIS files from the unzipped Ex10_FortRoss Links to an external site.folder downloaded previously. Add the shapefiles GPR Block 7 & 8 and Mag Block 8 as well as 2011ORTHO.TIF.
Geolocate: We will now shift the geophysical data points we brought in from Excel into their proper position in Sonoma County using the Block shapefiles as guides. We’ll do this using two slightly different methods: first, we’ll use Spatial Adjustment on it’s own with the Mag layer; second, we’ll move the points into the area and then use Spatial Adjustment.
1. Zoom in to Mag grid 8. Turn off GPR layer.2. Start Editing Mag grid 8. Turn on Snapping and choose “Use Snapping” from the menu. In the Snapping toolbar enable the Point and Vertex snapping tools (1st and 3rd from left).
3. Turn on Spatial Adjustment toolbar. This is like the Georeferencing toolbar but for adjusting vectors (including individual features and attributes).
4. Start with Mag grid 8 data. Right-click and “Zoom to Layer” Place links at one of the four corners and a connecting link to the mag block layer. Do each in turn by zooming in to each area (Z key) to place point carefully then return to previous view (blue arrow) and place. You may need to use “Zoom to Layer” to get back to mag grid 8 data (recall the georeferencing process when we connected points from one layer to features of another).
5. After all links have been placed, under Spatial Adjust choose “Set Adjust Data”. You can adjust individual vectors (objects) or entire files. Choose All Features and uncheck all but the Mag grid 8 box.
6. Spatial Adjustment > Adjust.
7. Happy with the results? Stop editing and save edits to Mag grid 8.
8. Save your Map MXD. Arcmap is liable to crash doing this work.
D. Georeference GPR vector points
Repeat with GPR layer but displace first to the neighborhood of Fort Ross and then use Spaital Adjustment. The GPR data needs to be rotated 90 degrees clockwise during Spatial Adjustment so it would be nice to have them on the same Display extent before using Spatial Adjust.
1. First turn on labeling on GPR Block layer because it helps with locating the data when zoomed far out. Right-click GPR Block layer > Label Features
2. Right click “GPR7_8” > Selection > Select All. Note the line at the bottom of Arcmap reads “Number of features selected = 1296”.
3. Click the globe (zooms out). Note the light blue blob on the bottom of the screen. The label “1” is all you can see of Fort Ross (make sure you have labeled the GPR Block 7-8 layer).
4. Start an editing session and select the black triangular tool just right of the “Editor” menu in the Editor toolbar.
5. Click on the light blue blob and drag it to the 1 until the snapping indicator appears and drop there. Note that Arcmap may crash at this point. If it didn’t crash Save edits and stop editing. Save your map.
6. At this point, you can unselect GPR layer, now that you have moved it. Now let’s use Spatial Adjustment as we did with Mag. Note that GPR came with the following comment “Grid needs to be rotated clockwise 90 degrees”. We can take that into account when we do the spatial adjustment For example the top-left corner in the original data becomes the top-right corner in the destination bounding box (GPR block layer).
a. Make sure to change “Set Adjust Data” to GPR because Mag will still be selected (and since they’re in the same folder / geodatabase they’re both editable when one is editable!).
7. Stop editing, save map.
8. Turn off labeling in GPR Block Layer.
E. Symbolize using Value fields.
These are 3D points (X,Y,Z position fields, plus a Value field).
Right click > Properties > Symbology
Quantities > Graduated Colors
Value: Value or Vrt_Grad
Select all 5 classes using the Shift key
Right-click > Properties for Selected Symbols
Choose Square1 (Black square) and make it Size: 10
Choose a Color Ramp (perhaps Red light to dark for Mag, Green light to dark for GPR
Part II. Interpolation and Visualization
A. Generate Raster layers
These data points can be converted directly to raster using the Point to Raster tool however the goal is to gain a generalized understanding of sub-surface features so we want the software to interpolate values in each cell using weighted values from adjacent cells in a process known as Kriging. Links to an external site.
Interpolate to Raster
1. Using Mag Data
Under ArcToolbox > Spatial Analyst > Interpolate > Kriging
- Input: Mag (point vectors)
- Z value field: Vert_Grd or Values
- Output Surface: Mag8_5cm.
- Semivariogram model: Linear
- Output Cell size: 0.05 m
- Leave the remaining using defaults.
The resulting data isn't displayed well. Go to Symbology > Classified and click Classify...
Change to Method: Quantiles, Classes: 18. Choose a palette with a wide spectrum of colors.
Clip Output
Rasters are generated to Extent of the layer and everything outside of Mag_Block8 polygon is junk (the blocks aren't aligned to axes). To clip to the polygon use
Arctoolbox > Data Management > Raster > Raster Processing > Clip
Use Input: MagKrig (or whatever you named it), Output Extent: Mag_Block8, Check "Use input features for clipping geometry". Resymbolize as above.
Generate contours using
Spatial Analyst > Surface > Contour.
Given the range of values (-320 to +1707) a contour interval of 20 is appropriate. Set the Base Contour to -100 (the extremes are noisy dipoles from metal)
2. Repeat using GPR Data.
Use a different spectrum of colors for comparison between the Mag and GPR layers.
3. Compare output
Set Transparency so you can see both. Which is more useful for interpreting features?
Use transparencies and layer order to upload. Label using only Title.
Save Map MXD.
B. Interpret and Create vector features
1. Create POLYLINE Shapefile, Add Title (30 chars), Source (50 chars), and Description (255 chars) text fields.
2. Trace Lines around features of interest and label with your observations.
Lines allow us to have incomplete features (they don't have to be closed like polygons), so they are more flexible for tracing features.
Save a PDF and upload as your Ex 10 assignment.
C. Export to Google Earth
It's easy to export data from Arcmap to Google Earth as a KMZ (a zipped KML). To demonstrate this and visit the texture-mapped reconstruction at Fort Ross virtually let's export the GPR point file (the Vector with 1296 points not the interpolated raster) to KMZ.
- Symbolize the GPR layer as you'd like it to appear in GoogleEarth (red squares size 11 showing Value field as Quantities as graduated colors, for example).
ArcToolbox > Conversion Tools > To KML > Layer to KML
- Save it to your project file.
- Leave "Clamp to ground" checked. Accept defaults.
- Quit Arcmap. Navigate to KMZ file and double-click. This should open Google Earth.
Note the GPR layer is in the Temporary Items in the Paces pane. Expand the layer (click +) and you can see every record in a sorted list.
- Double click to zoom in to layer at Fort Ross.
- Fly in and get low with Fort buildings visible. If you have a wheel-mouse try pressing the wheel like a button.
This process is just as simple with a raster/image such as a georeferenced map in Arcmap. We used a vector layer in this example to show that the attributes from the Shapefile's Attribute Table can be sent to GE as well. Note that Popup Infobox in the screenshot above shows the attribute table contents for a particular record.
The KMZ file can be saved to KML (unzipped) in GoogleEarth. Right click the layer in the Places pane and choose "Save Place As..." Type: KML format. That KML file can be opened in a Text editor and altered. The contents of the Pop-up info-box can be modified in the text editor. The Info box content is in HTML like a web browser. You can also include photos. There's great potential for a public outreach tour in GE or GoogleMaps (online) from your Arcmap project data.