Guth 1987 GSA Poster

Department of Geophysical Sciences

When I participated in the NSF short course Teaching Physics Using the World Wide Web in 1996, it became obvious to me that web pages are identical in form to poster sessions common at professional meetings in the geophysical sciences. It seemed to me that the WWW would be the ideal medium where undergraduate students could present the results of their senior research projects, capstone projects or independent studies and have them be seen by someone other than their advisors.

Since I throw nothing away (as evidenced by my office), I still had the poster materials I presented at the annual meeting of the Geological Society of America in 1987. This page tests the feasibility of using the web for student poster presentations.

REFERENCE

Guth, Lawrence R., 1987, Microcomputer contouring of stereographic projections: variations on a theme in 128k: Geological Society of America Abstracts with Programs, v. 19, No. 7, p. 688.

ABSTRACT

MICROCOMPUTER CONTOURING OF STEREOGRAPHIC PROJECTIONS: VARIATIONS ON A THEME IN 128k
            GUTH, Lawrence R., Department of Geology and Geophysics, Rice University Houston, TX 77251
Contoured lower hemisphere, equal area projections are useful to display many types of geologic data. Once a computer code is developed to produce these diagrams, small programming changes can lead to a large family of stand alone programs for various applications. I have written such a set of Turbo Pascal programs for the IBM-PC. All are based on using the Kalsbeek Counting Net as a finite element grid containing 331 nodes and 600 triangular elements. The nodes are treated as vectors in space during the counting process to assign scalar values to the nodes, while during contouring the nodes are considered to be locations in the plane of the projection. Data are stored in ASCII files with each datum read from file and acted on sequentially, so only the array of 331 scalar values need be stored in RAM with the number of data points limited solely by the size of mass storage.
            The base program generates contoured diagrams of linear or planar data from either field or microstructural observations. The user may select Schmidt, Kamb, or spherical tessellation counting methods or chose an arbitrary counting area. Output includes a scatter diagram of the input data and the contoured stereonet with a shading option.
            One variation of this program accommodates a different data structure and allows the user to do a dynamic analysis of a population of faults. The ASCII file is composed of the fault plane, the slip direction and the type of fault. For each fault in the population, the Schmid factor for each of the 331 nodes is calculated and the counting variable at the node is incremented if the Schmid factor is positive. Contouring then gives the principal stress axes by the method of Pressure and Tension dihedra (Angelier and Mechler, 1977).

ACKNOWLEDGMENTS

This programming was done to aid with research on the kinematics of the deformational structures found on Isla de Margarita, Venezuela. Field work for that research was supported by GSA Research Grants 3305-84 and 3437-85. Additional support was provided through a NSF grant to my advisor, Dr. Hans G. Avé Lallemant. Finally, this work would not have been possible without the Kary Data Scholarship awarded to me by the Geophysical Society of Houston. This scholarship enabled me to upgrade the IBM-PC I bought back in 1982 (when I worked for Big Oil and had money).

INTRODUCTION

KALSBEEK GRID
The Kalsbeek counting net is shown as it is used in this set of programs, a triangular finite element grid. All 331 nodes are shown, with some of them labeled. Half of the 600 triangular elements have been drawn in as well.

As will be seen in the flow chart, what I want to do with this grid is very simple. The major programming headaches are in the bookkeeping for all the nodes and triangles. Once the bookkeeping problems have been solved, you might as well get all the mileage out of the code as possible. That is the point of this poster.

The Kalsbeek grid is not the only grid used for the contouring of stereonets, and it may not be the best. The nodes are not quite uniformly distributed over the surface of the hemisphere and the spacing of the grid is rather coarse. The biggest hole in the net occurs just inside the primitive where locations exist that are almost 6 degrees away from any node. Fortunately, these shortcomings are only significant in rather extreme circumstances. The advantages of direct comparison to hand contoured diagrams and general familiarity with the grid I think outweigh its deficiencies.

HEMISPHERICAL CAPS
The illustration shows that a cone with its apex at the center of a sphere will intersect that sphere and cut out a circular spherical cap. This is generally referred to as a hemispherical cap, since we normally only use the lower hemisphere in contouring structural data and the area of this cap is easily related to the area of a hemisphere as shown. The angle q is called the search radius.
As will be seen, these caps are used with the SCHMIDT, USER and KAMB counting methods and with the MagPlot.PAS program in the following manner:

1) The direction cosines are calculated for the position of the data point and each node in turn.

2) Since these are unit vectors, the dot product between the data point and the node can be interpreted as the cosine of the angle between these two vectors.

3) All computations using hemispherical caps are done in the lower hemisphere, so we've got to remember to count points that are on the opposite side of the net. Due to the symmetry of the cosine function, IF WE USE THE ABSOLUTE VALUE OF THE DOT PRODUCT, WE TEST ONLY THE ACUTE ANGLE BETWEEN THE TWO VECTORS.

4) The absolute value of the dot product is tested against the cosine of the search radius as defined in the illustration.

FLOW CHART INTRODUCTION
I program to help me with other research. I don't use the programming as an end unto itself. The result is that the programs I write are unsophisticated and probably do not use the full power of the computer and the programming language. The earliest version of the contouring program was written in IBM-PC BASICA 1.1 and later translated into Turbo Pascal 3.0. Although the programs presented here are greatly improved over the original BASICA version, the approach to the programming has remained the same.
Turbo Pascal 3.0 has a maximum code size of 64k and a maximum data size of 64k. For my purposes, it is easier to live within these restrictions than to learn how to use overlays and stacks and heaps that I didn't have to use in BASIC and FORTRAN. As a result, each of the programs presented requires no more than 128k of computer memory.

THE PROGRAM GUTS
The base driver is simplified to show the features common for all of these programs. It consists of a nested loop, with the outer loop running through the input data and the inner loop testing each node against a data point. Changing the test criterion fundamentally changes the output, and this forms the foundation of this entire poster.

COMMON MODULES
Modules composed of two or more related Pascal procedures (subroutines) have been written to be general (to my CGA graphics, IBM Proprinter system) so they can pass unmodified between the several different programs. These do the following operations:
--- standard arithmetic functions;
--- standard procedures for bit-image graphics;
--- standard procedures to do stereonet rotations;
--- standard procedures to generate finite element grid in space and stereographic projection; and
--- standard procedures for finite element contouring.

I/O PROCEDURES
The I/O procedures cannot be made perfectly general, because some programs require different input data structures and the output needs new titles for proper identification. However, the base program can often serve as a template, so that these changes are easily made.

STEREONET CONTOURING

The program SContour.PAS is the base program from which all the other programs were derived. In fact, the idea of changing the test criterion to alter the counting methods led to the realization that other data structures and test criterion could be incorporated into the basic framework of the program to spin off other programs that did very different analyses.
The four different counting methods supported by the program are described separately. The output includes scatter diagrams showing the input data, unshaded contour diagrams that can then be colored for slides, contour diagrams with the source data superimposed, and shaded contour diagrams.

THE CONTOURED DATA SET Within the contouring program, the test that relates a data position to the nodes of the finite element grid can be set to 4 different counting methods. The results from these 4 methods are shown so you can see how changing the test of the counting process generates sometimes markedly different contoured diagrams. To make this a valid comparison, the same data set was used in each of the examples. The following just describes the data set for anyone who might be interested. These are quartz c-axis petrofabric data from the Matasiete Trondhjemite outcropping at Puerto Fermín on Isla de Margarita, Venezuela. Approximately 100 c-axes were measured in each of three mutually perpendicular thin sections on a universal stage. These U-stage data were reduced using a series of programs I've developed to help with a research project on the kinematics of the deformational structures exposed in the rocks of Margarita. Here the data set has been rotated from the geographic reference frame so that the mesoscopic foliation is vertical, striking EW and the lineation is horizontal in the foliation plane. The dotted great circle girdle on these diagrams is the rotated geographic horizon. A scatter diagram of the input data is shown on the right and several other diagrams relating these data to Bingham and Fisher distributions are shown below.

The Bingham distribution is a simple eigen vector analysis of the sample population to find the maximum, intermediate and minimum moment axes. The ellipses about each axis is the a₉₅ "cone" of confidence for the location of each axis. The bottom two displays were generated by SferStat.PAS, a separate program I wrote to examine the spherical statistics of my structural and petrofabric data.	These histograms show the angular relationship between the sample data and the moment axes (top three) and the Fisher axis (bottom). The Fisher axis is the vector sum of the sample population when each is treated as a unit vector and all vector additions are done using the ends of the vectors forming acute angles. The class width is 3 degrees and the "+" indicates the expected class percentage for a uniform distribution.

SPHERICAL TESSELLATION
     The dictionary defines a tessellation to be a mosaic composed of small square tiles. This definition is expanded in mathematics to include any shape or set of shapes that can cover a surface with no gaps or overlaps. It follows that a spherical tessellation must use a shape or set of shapes to completely cover the surface of a sphere with no gaps or overlaps.
     It's readily apparent that the circular hemispherical caps used with the other counting methods will not tessellate -- that some other shape must be used. Furthermore, it turns out that the Kalsbeek counting net which forms the finite element grid for this entire family of programs has nodal points that are not quire uniformly distributed over the surface of a sphere. Working with actual shapes would be cumbersome and likely each spherical polygon of the tessellation would have to be individually defined. The following visualization was used to define another property of the tessellation which could be more easily tested.
     Imagine that the hemisphere is made of double-pane insulated glass. The 331 nodes are seen as straws that go through the inner pane of glass and end in the dead air space between the panes. If you now started blowing soap bubbles through these straws so that each bubble nucleated at the same time and grew at the same rate, the bubbles would form circular disks centered at the node and grow outward until they started to impinge on each other. The final result would be a tessellation formed by a set of spherical polygons, each centered over its associated node. Because these spherical polygons were generated from circular disks all growing at the same rate, each point within the polygon is closer to the central node than to any other node of the grid. THIS PROPERTY OF THE TESSELLATION WAS USED AS THE TEST IN THE COUNTING PROCESS SO THAT ONLY THE BIN ASSOCIATED WITH THE NODE CLOSEST TO THE DATA POINT IS INCREMENTED.

CONSIDERATIONS: Viewed in the above manner, this particular spherical tessellation counting method in fact quantizes the data so that they exist only at the nodal locations and no where else. The test for the spherical tessellation counting used here specifies that each datum be shifted to the closest mode possible before contouring. For comparison with the other counting methods, the program keeps track of the largest angular shift needed to place a datum coincident with a node, and the area of a circular hemispherical cap with that radius is reported as a "counting area".

SCHMIDT COUNTING METHOD
The Schmidt counting method increments all the nodes within a 1% area centered on each data point. From the previous discussion under the heading "HEMISPHERICAL CAPS" cos q = 1 - 0.01
cos q = 0.99
q = 8.1 degrees

and the program's test becomes IF THE ABSOLUTE VALUE OF THE DOT PRODUCT BETWEEN THE DATA POINT AND THE NODE IS GREATER THAN OR EQUAL TO 0.99, THEN INCREMENT THE BIN ASSOCIATED WITH THAT NODE.
After I got out the plot with a contour interval of 1% per 1% area, I thought that the girdle could be better defined using a different contour interval. The program allows the user to chose new contour intervals without recounting. The second plot at 0.5% per 1% area I think shows the pattern much better (in the original printout).



Contour Interval = 1% per 1% area Unshaded area is 0%	Contour Interval = 0.5% per 1% area Unshaded area is 0% (lower right stereonet only)

USER DEFINED COUNTING AREA
The user defined counting area method increments all the nodes within some defined area centered on each data point. From the previous discussion under the heading "HEMISPHERICAL CAPS" cos q = 1 - (Percentage area/100)

and the program's test becomes IF THE ABSOLUTE VALUE OF THE DOT PRODUCT BETWEEN THE DATA POINT AND THE NODE IS GREATER THAN OR EQUAL TO cos q THEN INCREMENT THE BIN ASSOCIATED WITH THAT NODE. In this example, the counting area was chosen to be 2% of the area of the hemisphere.

From the notes about the Kalsbeek counting net used as a finite element grid, it was mentioned that if the counting cap became too small, it would be possible to have a data point slip between the nodes uncounted. For the Kalsbeek grid, the user defined counting area can be no smaller than 0.55% of the area of the hemisphere.

KAMB COUNTING METHOD
Kamb (1959) proposed using a counting area which is related to the size of the sample population (N). He defined this counting area so that for a uniform distribution sampled randomly, the expected number of points falling in the counting area would be three times the standard deviation for the sample population. Rearranging his relationship, this fractional area of the hemisphere is Area = 9/(N + 9)

From the previous discussion under the heading "HEMISPHERICAL CAPS"

cos q = 1 - 9/(N + 9)

and the program's test becomes IF THE ABSOLUTE VALUE OF THE DOT PRODUCT BETWEEN THE DATA POINT AND THE NODE IS GRATER THAN OR EQUAL TO cos q, THEN INCREMENT THE BIN ASSOCIATED WITH THAT NODE.
From the notes about the Kalsbeek counting net used as a finite element grid, it was mentioned that if the counting cap became too small, it would be possible to have a data point slip between the nodes uncounted. For the Kalsbeek grid, the smallest area is 0.55%. This then puts an upper limit on the size of the sample population that can be contoured using a Kamb counting cap with the Kalsbeek grid -- specifically, N can be no greater than 1,627.

Kamb Shading
The Kamb counting method is actually comparing your data to a uniform population and I've chosen a shading scheme to emphasize this. Rearranging the equation of Kamb and put into the program's terminology:

E = expected nodal (bin) value = 9N/(N + 9)
s = standard deviation = 3N/(N + 9)
O = observed nodal (bin) value

I contour the values of (O - E)/s at each node, or the number of standard deviations above or below the mean for a uniform distribution sampled randomly.
My eye tends to be attracted to dark shading in black and white illustrations, so I've constructed a gray scale where the shading becomes darker going in both negative and positive directions away from the mean. Above the gray scale shown is the standard normal distribution which gives a probability of 95.4% that a node will have a value within 2s of the mean when a uniform population is randomly sampled. I use this as approximating a 95% confidence level, that contour levels greater than 2s and less than -2s indicate the fabric is significantly different from a uniform distribution. I find this particularly useful in interpreting quartz c-axis data, where significant holes in the fabric are just as important as the girdle of high concentrations. To show the effects of this shading, I've contoured the data with the user defined area option, setting the counting area to be the same as that used in the Kamb counting.



Contour Interval = 1s from the mean Unshaded area is within 1s of the mean for a uniform distribution	Contour Interval = 1% per 2.98% area Unshaded area is pole free (0%)

STEREONET CONTOURING SUMMARY DIAGRAM
For this web poster, I've added this section to make comparisons of the contouring methods easier to see. With the three large panels available for posters at the GSA meeting, this synoptic display was not needed. Due to the marginal results of the scans, I've displayed unshaded contours with the data points superimposed.

TESSELLATION CONTOURING Contour Interval = 0.5% per 0.54% area	SCHMIDT CONTOURING Contour Interval = 1.00% per 1.00% area


USER DEFINED CONTOURING Contour Interval = 1.00% per 2.00% area	KAMB CONTOURING Contour Interval = ±1.00 standard deviation from the mean

ALPHA-95 CONTOURING

MAGPLOT.PAS
      As I've mentioned, this programming was done to help with our research on Isla de Margarita, Venezuela. Margarita lies about 20 km off the Venezuelan mainland in the Caribbean Sea and is therefore caught up in the dextral deformation associated with the Caribbean-South American plate boundary. Paleomagnetics from many islands along this boundary suggest that major rotations have occurred due to the Caribbean-South American plate interaction. In relating our kinematic information obtained on Margarita to possible causative plate geometries, it will be important to remove any post-deformational rotations.
     Hargraves and Skerlec (1980) published paleomagnetic results from Margarita which other authors have used to argue that Margarita is unrotated. Plotting these data up, I could not see it. To try to simplify the pattern, I modified the contouring program SContour.PAS
     I wanted to look at the regions of overlapping alpha-95 (a₉₅) cones of confidence. Instead of using a single search radius for all data points as was done with the contouring program, I change the search radius for each paleomagnetic vector, setting the search radius to be equal to a₉₅. The output is displayed, and I still cannot convince myself that the island is unrotated from these data. This is included to show the first and easiest modification to the contouring program.

FRY PLOT CONTOURING

FRY STRAIN ANALYSIS
     When doing Fray strain analysis, my eye is often fooled because I tend to see the hole in the pattern and not the high concentration of points defining the centers of the nearest neighbors. Contouring seemed to be the best way of emphasizing patterns within the Fry plot.
     Going from spherical data to a stereographic projection is a one-to-one map, so there exists an inverse map to go from the projection back to spherical data. Alternatively, planar data can be projected onto the lower hemisphere using the same inverse map. Since I already had written the code to contour lower hemisphere stereographic projections, I wrote my Fry program to write out to file all points within a user defined radius from the center of the Fry plot. These points are stored on file as plunge and trend data, having been mapped from the plane onto the lower hemisphere. This file is then directly readable to the stereonet contouring program I had already written.
     The output of the stereonet contouring program is shown for test data digitized from Ramsay and Huber's (1983) Figure 7-7. Because of the large number of data points created in the analysis (here 2138), the Kamb counting method could not be used. Spherical tessellation was used to show the data with the least smoothing possible. The central peak of the crater results when you lose track of where you are on the digitizing tablet and input the same point twice.

CONSIDERATIONS
As already mentioned, the large number of points may make the search radius of the Kamb counting method too small for the Kalsbeek grid. Other special considerations come from the fact that these are planar data that have been treaded as spherical data. The most important of these will be that data from near the primitive but on opposite sides of the stereonet will be counted as being close together, when in fact these data are the farthest apart possible. There will therefore be a rim on the outer edge of the diagram that will contain false data, with the width of this rim depending on the counting scheme used. When the projected file is created, the user just defines the Fry plot field-of-view radius large enough so that the edge data can be later disregarded.
The second consideration is the well known fact that the equal area projection turns counting circles on the sphere into ellipses in the planar projection. This is shown in the diagram on the right that I generated with the program MagPlot.PAS. The result of this is a slight spatial filtering, emphasizing tangential patterns over radial. This last consideration benefits our goal of defining an ellipse centered on the stereonet.

DYNAMIC FAULT ANALYSIS

The above illustration shows the coordinate frame used to change the stereonet contouring program into a program used to perform a dynamic analysis for a population of faults based on the method of Pressure and Tension dihedra of Angelier and Mechler (1977). The data file read by the program contains for each fault the orientation of the fault plane, the plunge and trend of the slip direction and the type of fault (normal, reverse, dextral or sinistral). A coordinate frame is hung on the hanging wall block so that its outward normal (the lower hemisphere pole to the fault plane) is assigned to be +N_k and the direction of motion of the hanging wall block is +L_k. Unlike any other quantity in any other program, +L_k can be in either the upper or lower hemisphere.
Against this coordinate frame, each of the 331 nodes of the Kalsbeek finite element grid is tested. The test criterion is the Schmid factor, defined as Schmid factor = cos(+N_k ^ node) X cos(+L_k^ node)
or, in terms of vector dot products
Schmid factor = (+N_k • node)(+L_k • node) Comparing this equation with the illustration shows that if the Schmid factor is positive, the node must exist in the Pressure Quadrant. A maximum compressive stress axis (s₁) parallel to that node would produce a shear stress on the fault plane in the sip direction that would give the observed sense of shear. If the Schmid factor is negative, the maximum compressive stress parallel to that node would produce the wrong sense of motion on the fault. Incrementing all nodal bins with a positive Schmid factor defines the region in which the maximum compressive stress had to exist to get that fault to move in the observed direction. Looping through all faults of the sample population and contouring the nodal values superimposes the Pressure Quadrants. Hopefully, a region in space will be defined in which the maximum compressive stress had to exist to be compatible with all observed fault motions. Field data are used to infer that all faults in the population are coeval and therefore represent brittle deformation under the same stress field.

The output presented used test data from Angelier (1979). The fault plane and slickenside striation data are shown, and all of the faults were observed to have normal motions. The dynamic analysis using this program agrees well with the results obtained with the inversions done by Angelier (1979) and Michael (1984).

Slickenside striations (N = 38) Contour Interval 2% per 1% area	Poles to fault planes (N = 38) Contour Interval 2% per 1% area


Faults compatible with a maximum compressive stress axis parallel to nodal orientation Contour Interval = 10% of sample population	Faults compatible with a maximum compressive stress axis parallel to nodal orientation Contour Interval = 10% of sample population (Unshaded version showing the inferred principal stress axes based on the program output)

REFERENCES CITED

TEST RESULTS

     It took me 15 hours to put this poster session together using Page Composer in Netscape Communicator. Setting up the HTML document took essentially no time at all. I spent 5 hours just scanning in figures. Because the scanner software is so old in the Cartographic Lab, it did not support any of the file formats currently acceptable to Web browsers. I had to scan in the image, save it as a Mac PICT file, then use GIFConverter to translate it into a JPG file. Current graphics programs that save directly in JPG or GIF format have got to produce better figures.
     I retyped the entire poster into Page Composer, since the original documents were WordStar files on 5.25-inch floppies. I thought retyping would be faster than finding a machine with a 5.25-inch drive that could also translate WordStar into ASCII. The retyping, however, took the remaining 10 hours -- with some of that time trying to remember the organization of the poster presented 11 years ago.
     Based on this test, I would expect undergraduate students should be able to convert their final project into an HTML poster in under a day. This assumes that the project is still fresh in their heads and they generated the graphics and text for the normal paper version of their final report using the latest software available in the Fitchburg State College open computer labs. In that case, the same graphics and text files could be directly imported into the HTML authoring software also available in the open computer labs on campus.