next up previous contents index
Next: 2.1.3 More exercises Up: 2.1 General Information Previous: 2.1.1 Examples   Contents   Index


2.1.2 Exercises

  1. Plot the data as a green-blue polygon instead.

  2. Try using a predefined pattern.

A common question is : ``How can I plot symbols connected by a line with psxy?''. The surprising answer is that we must call psxy twice. While this sounds cumbersome there is a reason for this: Basically, polygons need to be kept in memory since they may need to be clipped, hence computer memory places a limit on how large polygons we may plot. Symbols, on the other hand, can be plotted one at the time so there is no limit to how many symbols one may plot. Therefore, to connect symbols with a line we must use the overlay approach:

 
psxy data -R -JX -B -P -K -W0.5p >! plot.ps
psxy data -R -JX -O -W -Si0.2i >> plot.ps

Our final psxy example involves a more complicated scenario in which we want to plot the epicenters of several earthquakes over the background of a coastline basemap. We want the symbols to have a size that reflects the magnitude of the earthquakes, and that their color should reflect the depth of the hypocenter. You will find the two files quakes.ngdc and quakes.cpt in your directory. The first few lines in the quakes.ngdc looks like this:

Historical Tsunami Earthquakes from the NGDC Database
Year  Mo  Da  Lat+N  Long+E  Dep  Mag
1987  01  04  49.77  149.29  489  4.1
1987  01  09  39.90  141.68  067  6.8

Thus the file has three header records (including the blank line), but we are only interested in columns 5, 4, 6, and 7. In addition to extract those columns we must also scale the magnitudes into symbols sizes in inches. Given their range it looks like multiplying the magnitude by 0.02 will work well. Reformatting this file to comply with the psxy input format can be done in a number of ways, including manual editing, using MATLAB, a spreadsheet program, or UNIX tools. Here, without further elaboration, we simply use the UNIX tool awk to do the job ($5 refers to the 5'th column etc., and NR is the current record number):

awk '{if (NR > 3) print $5, $4, $6, 0.02*$7}' quakes.ngdc >! quakes.d

The awk statement is automatically applied to each record, hence the output file quakes.d should now look like this (try it!):

149.29  49.77  489  0.082
141.68  39.90  067  0.136
...etc etc

We will follow conventional color schemes for seismicity and assign red to shallow quakes (depth 0-100 km), green to intermediate quakes (100-300 km), and blue to deep earthquakes (depth $>$ 300 km). The quakes.cpt file establishes the relationship between depth and color:

# color palette for seismicity
#z0  red  green blue    z1   red green blue
0    255      0    0   100   255     0    0
100    0    255    0   300     0   255    0
300    0      0  255  1000     0     0  255

Apart from comment lines (starting with #), each record in the cpt file governs the color of a symbol whose z value falls in the range between $z_0$ and $z_1$. If the lower and upper red/green/blue triplets differ then an intermediate color will be linearly interpolated given the $z$ value. Here, we have chosen constant color intervals.

We may now complete our example using the Mercator projection; we throw in a map scale out of pure generosity:

 
pscoast -R130/150/35/50 -JM6i -B5 -P -G200 -Lf134/49/42.5/500 -K >! map.ps
psxy -R -JM -O -Cquakes.cpt quakes.d -Sci -W0.25p >> map.ps

where the i appended to the -Sc option ensures that symbols sizes are interpreted to be in inches.


next up previous contents index
Next: 2.1.3 More exercises Up: 2.1 General Information Previous: 2.1.1 Examples   Contents   Index
Paul Wessel 2006-05-31