1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299
|
<!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 3.2 Final//EN">
<!--Converted with LaTeX2HTML 99.2beta8 (1.46)
original version by: Nikos Drakos, CBLU, University of Leeds
* revised and updated by: Marcus Hennecke, Ross Moore, Herb Swan
* with significant contributions from:
Jens Lippmann, Marek Rouchal, Martin Wilck and others -->
<HTML>
<HEAD>
<TITLE>6.3 Spectral estimation and xy-plots</TITLE>
<META NAME="description" CONTENT="6.3 Spectral estimation and xy-plots">
<META NAME="keywords" CONTENT="GMT_Docs">
<META NAME="resource-type" CONTENT="document">
<META NAME="distribution" CONTENT="global">
<META HTTP-EQUIV="Content-Type" CONTENT="text/html; charset=iso-8859-1">
<META NAME="Generator" CONTENT="LaTeX2HTML v99.2beta8">
<META HTTP-EQUIV="Content-Style-Type" CONTENT="text/css">
<LINK REL="STYLESHEET" HREF="GMT_Docs.css">
<LINK REL="next" HREF="node77.html">
<LINK REL="previous" HREF="node75.html">
<LINK REL="up" HREF="node73.html">
<LINK REL="next" HREF="node77.html">
</HEAD>
<BODY bgcolor="#ffffff">
<!--Navigation Panel-->
<A NAME="tex2html2094"
HREF="node77.html">
<IMG WIDTH="37" HEIGHT="24" ALIGN="BOTTOM" BORDER="0" ALT="next" SRC="next.gif"></A>
<A NAME="tex2html2088"
HREF="node73.html">
<IMG WIDTH="26" HEIGHT="24" ALIGN="BOTTOM" BORDER="0" ALT="up" SRC="up.gif"></A>
<A NAME="tex2html2082"
HREF="node75.html">
<IMG WIDTH="63" HEIGHT="24" ALIGN="BOTTOM" BORDER="0" ALT="previous" SRC="prev.gif"></A>
<A NAME="tex2html2090"
HREF="node1.html">
<IMG WIDTH="65" HEIGHT="24" ALIGN="BOTTOM" BORDER="0" ALT="contents" SRC="contents.gif"></A>
<A NAME="tex2html2092"
HREF="node149.html">
<IMG WIDTH="43" HEIGHT="24" ALIGN="BOTTOM" BORDER="0" ALT="index" SRC="index.gif"></A>
<BR>
<B> Next:</B> <A NAME="tex2html2095"
HREF="node77.html">6.4 A 3-D perspective</A>
<B> Up:</B> <A NAME="tex2html2089"
HREF="node73.html">6. Cook-book</A>
<B> Previous:</B> <A NAME="tex2html2083"
HREF="node75.html">6.2 Image presentations</A>
  <B> <A NAME="tex2html2091"
HREF="node1.html">Contents</A></B>
  <B> <A NAME="tex2html2093"
HREF="node149.html">Index</A></B>
<BR>
<BR>
<!--End of Navigation Panel-->
<H1><A NAME="SECTION001430000000000000000"></A>
<A NAME="10122"></A>
<A NAME="10123"></A>
<BR>
6.3 Spectral estimation and xy-plots
</H1>
<P>
In this example we will show how to use the <A NAME="tex2html370"
HREF="http://www.soest.hawaii.edu/gmt"><B>GMT</B></A> programs
<A NAME="tex2html371"
HREF="../fitcircle.html"><I><B>fitcircle</B></I></A><A NAME="10245"></A>, <A NAME="tex2html372"
HREF="../project.html"><I><B>project</B></I></A><A NAME="10254"></A>, <A NAME="tex2html373"
HREF="../sample1d.html"><I><B>sample1d</B></I></A><A NAME="10263"></A>,
<A NAME="tex2html374"
HREF="../spectrum1d.html"><I><B>spectrum1d</B></I></A><A NAME="10272"></A>, <A NAME="tex2html375"
HREF="../psxy.html"><I><B>psxy</B></I></A><A NAME="10281"></A>, and <A NAME="tex2html376"
HREF="../pstext.html"><I><B>pstext</B></I></A><A NAME="10290"></A>.
Suppose you have (lon, lat, gravity) along a satellite track
in a file called <U>sat.xyg</U>, and (lon, lat, gravity)
along a ship track in a file called <U>ship.xyg</U>.
You want to make a cross-spectral analysis of these data.
First, you will have to get the two data sets into equidistantly
sampled time-series form. To do this, it will be convenient to
project these along the great circle that best fits the sat track.
We must use <A NAME="tex2html377"
HREF="../fitcircle.html"><I><B>fitcircle</B></I></A><A NAME="10301"></A> to find this great circle and choose
the L<IMG
WIDTH="11" HEIGHT="29" ALIGN="MIDDLE" BORDER="0"
SRC="img12.gif"
ALT="$_2$"> estimates of best pole. We project the data using
<A NAME="tex2html378"
HREF="../project.html"><I><B>project</B></I></A><A NAME="10310"></A> to find out what their ranges are in the projected
coordinate. The <A NAME="tex2html379"
HREF="../minmax.html"><I><B>minmax</B></I></A><A NAME="10319"></A> utility will report the minimum and
maximum values for multi-column ASCII tables. Use this information
to select the range of the projected distance coordinate they have
in common. The script prompts you for that information after
reporting the values. We decide to make a file of equidistant
sampling points spaced 1 km apart from -1167 to +1169, and use
the <I>UNIX</I> utility <I>$AWK</I><A NAME="10328"></A> to accomplish this step. We can then
resample the projected data, and carry out the cross-spectral
calculations, assuming that the ship is the input and the satellite
is the output data. There are several intermediate steps that
produce helpful plots showing the effect of the various processing
steps (<U>example_3[a-f].ps</U>), while the final plot
<U>example_03.ps</U> shows
the ship and sat power in one diagram and the coherency on another
diagram, both on the same page. Note the extended use of
<A NAME="tex2html380"
HREF="../pstext.html"><I><B>pstext</B></I></A><A NAME="10336"></A> and <A NAME="tex2html381"
HREF="../psxy.html"><I><B>psxy</B></I></A><A NAME="10345"></A> to put labels and legends directly on
the plots. For that purpose we often use <B>-Jx</B>1i and specify
positions in inches directly. Thus, the complete automated
script reads:
<P>
<P>
<BR>
<P>
<BR CLEAR="ALL">
<HR>
<BR>
<PRE>
fitcircle sat.xyg -L2 >! report
set cpos = `grep "L2 Average Position" report`
set ppos = `grep "L2 N Hemisphere" report`
project sat.xyg -C$cpos[1]/$cpos[2] -T$ppos[1]/$ppos[2] -S -Fpz -Q >! sat.pg
project ship.xyg -C$cpos[1]/$cpos[2] -T$ppos[1]/$ppos[2] -S -Fpz -Q >! ship.pg
set plotr = `cat sat.pg ship.pg | minmax -I100/25 -C`
gmtset MEASURE_UNIT INCH
psxy -R$plotr[1]/$plotr[2]/$plotr[3]/$plotr[4] -U/-1.75i/-1.25i/"Example 3a in Cookbook" \
-JX8i/5i -X2i -Y1.5i -K -W1p sat.pg \
-Ba500f100:"Distance along great circle":/a100f25:"Gravity anomaly (mGal)":WeSn >! example_03a.ps
psxy -R -JX -O -Sp0.03i ship.pg >> example_03a.ps
$AWK '{ if (NR > 1) print $1 - last1; last1 = $1; }' ship.pg | pshistogram -W0.1 -G0 -JX3i -K \
-X2i -Y1.5i -B:."Ship": -U/-1.75i/-1.25i/"Example 3b in Cookbook" >! example_03b.ps
$AWK '{ if (NR > 1) print $1 - last1; last1 = $1; }' sat.pg | pshistogram -W0.1 -G0 -JX3i -O \
-X5i -B:."Sat": >> example_03b.ps
head -1 ship.pg >! ship.pg.extr
head -1 sat.pg >! sat.pg.extr
paste ship.pg.extr sat.pg.extr | $AWK '{ if ($1 > $3) print int($1); else print int($3); }' \
>! sampr1
tail -1 ship.pg >! ship.pg.extr
tail -1 sat.pg >! sat.pg.extr
paste ship.pg.extr sat.pg.extr | $AWK '{ if ($1 < $3) print int($1); else print int($3); }' \
>! sampr2
set sampr = `paste sampr1 sampr2`
$AWK 'BEGIN { for (i = '$sampr[1]'; i <= '$sampr[2]'; i++) print i }' /dev/null >! samp.x
sample1d sat.pg -Nsamp.x >! samp_sat.pg
filter1d ship.pg -Fm1 -T$sampr[1]/$sampr[2]/1 -E | sample1d -Nsamp.x >! samp_ship.pg
psxy -R$plotr[1]/$plotr[2]/$plotr[3]/$plotr[4] -JX8i/5i -X2i -Y1.5i -K -W1p samp_sat.pg \
-Ba500f100:"Distance along great circle":/a100f25:"Gravity anomaly (mGal)":WeSn \
-U/-1.75i/-1.25i/"Example 3c in Cookbook" >! example_03c.ps
psxy -R -JX -O -Sp0.03i samp_ship.pg >> example_03c.ps
paste samp_ship.pg samp_sat.pg | cut -f2,4 | spectrum1d -S256 -D1 -W -C >& /dev/null
psxy spectrum.coh -Ba1f3p:"Wavelength (km)":/a0.25f0.05:"Coherency@+2@+":WeSn -JX-4il/3.75i \
-R1/1000/0/1 -U/-2.25i/-1.25i/"Example 3d in Cookbook" -P -K -X2.5i -Sc0.07i -G0 \
-Ey/2 -Y1.5i >! example_03.ps
echo "3.85 3.6 18 0.0 1 11 Coherency@+2@+" | pstext -R0/4/0/3.75 -Jx1i -O -K >> example_03.ps
cat << END >! box.d
2.375 3.75
2.375 3.25
4 3.25
END
psxy -R -Jx -O -K -W1.5p box.d >> example_03.ps
psxy -Ba1f3p/a1f3p:"Power (mGal@+2@+km)"::."Ship and Satellite Gravity":WeSn spectrum.xpower \
-St0.07i -O -R1/1000/0.1/10000 -JX-4il/3.75il -Y4.2i -K -Ey/2 >> example_03.ps
psxy spectrum.ypower -R -JX -O -K -G0 -Sc0.07i -Ey/2 >> example_03.ps
echo "3.9 3.6 18 0.0 1 11 Input Power" | pstext -R0/4/0/3.75 -Jx -O -K >> example_03.ps
psxy -R -Jx -O -K -W1.5p box.d >> example_03.ps
psxy -R -Jx -O -K -G240 -L -W1.5p << END >> example_03.ps
0.25 0.25
1.4 0.25
1.4 0.9
0.25 0.9
END
echo "0.4 0.7" | psxy -R -Jx -O -K -St0.07i -G0 >> example_03.ps
echo "0.5 0.7 14 0.0 1 5 Ship" | pstext -R -Jx -O -K >> example_03.ps
echo "0.4 0.4" | psxy -R -Jx -O -K -Sc0.07i -G0 >> example_03.ps
echo "0.5 0.4 14 0.0 1 5 Satellite" | pstext -R -Jx -O >> example_03.ps
trend1d -Fxw -N2r samp_ship.pg >! samp_ship.xw
psxy -R$plotr[1]/$plotr[2]/$plotr[3]/$plotr[4] -JX8i/4i -X2i -Y1.5i -K -Sp0.03i \
-Ba500f100:"Distance along great circle":/a100f25:"Gravity anomaly (mGal)":WeSn \
-U/-1.75i/-1.25i/"Example 3e in Cookbook" samp_ship.pg >! example_03d.ps
psxy -R$plotr[1]/$plotr[2]/0/1.1 -JX8i/1.1i -O -Y4.25i -Bf100/a0.5f0.1:"Weight":Wesn -Sp0.03i \
samp_ship.xw >> example_03d.ps
trend1d -Fxrw -N2r samp_ship.pg | $AWK '{ if ($3 > 0.6) print $1, $2 }' | sample1d -Nsamp.x >! \
samp2_ship.pg
trend1d -Fxrw -N2r samp_sat.pg | $AWK '{ if ($3 > 0.6) print $1, $2 }' | sample1d -Nsamp.x >! \
samp2_sat.pg
set plotr = `cat samp2_sat.pg samp2_ship.pg | minmax -I100/25 -C`
psxy -R$plotr[1]/$plotr[2]/$plotr[3]/$plotr[4] -JX8i/5i -X2i -Y1.5i -K -W1p \
-Ba500f100:"Distance along great circle":/a50f25:"Gravity anomaly (mGal)":WeSn \
-U/-1.75i/-1.25i/"Example 3f in Cookbook" samp2_sat.pg >! example_03e.ps
psxy -R -JX -O -Sp0.03i samp2_ship.pg >> example_03e.ps
paste samp2_ship.pg samp2_sat.pg | cut -f2,4 | spectrum1d -S256 -D1 -W -C >& /dev/null
psxy spectrum.coh -Ba1f3p:"Wavelength (km)":/a0.25f0.05:"Coherency@+2@+":WeSn -JX-4il/3.75i \
-R1/1000/0/1 -U/-2.25i/-1.25i/"Example 3g in Cookbook" -P -K -X2.5i -Sc0.07i -G0 \
-Ey/2 -Y1.5i >! example_03f.ps
echo "3.85 3.6 18 0.0 1 11 Coherency@+2@+" | pstext -R0/4/0/3.75 -Jx -O -K >> example_03f.ps
cat << END >! box.d
2.375 3.75
2.375 3.25
4 3.25
END
psxy -R -Jx -O -K -W1.5p box.d >> example_03f.ps
psxy -Ba1f3p/a1f3p:"Power (mGal@+2@+km)"::."Ship and Satellite Gravity":WeSn spectrum.xpower \
-St0.07i -O -R1/1000/0.1/10000 -JX-4il/3.75il -Y4.2i -K -Ey/2 >> example_03f.ps
psxy spectrum.ypower -R -JX -O -K -G0 -Sc0.07i -Ey/2 >> example_03f.ps
echo "3.9 3.6 18 0.0 1 11 Input Power" | pstext -R0/4/0/3.75 -Jx -O -K >> example_03f.ps
psxy -R -Jx -O -K -W1.5p box.d >> example_03f.ps
psxy -R -Jx -O -K -G240 -L -W1.5p << END >> example_03f.ps
0.25 0.25
1.4 0.25
1.4 0.9
0.25 0.9
END
echo "0.4 0.7" | psxy -R -Jx -O -K -St0.07i -G0 >> example_03f.ps
echo "0.5 0.7 14 0.0 1 5 Ship" | pstext -R -Jx -O -K >> example_03f.ps
echo "0.4 0.4" | psxy -R -Jx -O -K -Sc0.07i -G0 >> example_03f.ps
echo "0.5 0.4 14 0.0 1 5 Satellite" | pstext -R -Jx -O >> example_03f.ps
\rm -f box.d report samp* *.pg *.extr spectrum.* .gmtcommands
</PRE>
<P>
<BR>
<BR>
<P>
<BR CLEAR="ALL">
<HR>
<P>
The final illustration (Figure <A HREF="node76.html#fig:GMT_example_03">6.3</A>) shows that the ship gravity
anomalies have more power than altimetry derived gravity for
short wavelengths and that the coherency between the two
signals improves dramatically for wavelengths <IMG
WIDTH="16" HEIGHT="29" ALIGN="MIDDLE" BORDER="0"
SRC="img36.gif"
ALT="$>$"> 20 km.
<P></P>
<DIV ALIGN="CENTER"><A NAME="fig:GMT_example_03"></A><A NAME="10404"></A>
<TABLE>
<CAPTION ALIGN="BOTTOM"><STRONG>Figure 6.3:</STRONG>
Spectral estimation and <IMG
WIDTH="26" HEIGHT="31" ALIGN="MIDDLE" BORDER="0"
SRC="img3.gif"
ALT="$x/y$">-plots</CAPTION>
<TR><TD><IMG
WIDTH="328" HEIGHT="544" BORDER="0"
SRC="img100.gif"
ALT="\begin{figure}\centering\epsfig{figure=eps/GMT_example_03.eps}\end{figure}"></TD></TR>
</TABLE>
</DIV><P></P>
<P>
<A NAME="10396"></A>
<A NAME="10397"></A>
<P>
<HR>
<!--Navigation Panel-->
<A NAME="tex2html2094"
HREF="node77.html">
<IMG WIDTH="37" HEIGHT="24" ALIGN="BOTTOM" BORDER="0" ALT="next" SRC="next.gif"></A>
<A NAME="tex2html2088"
HREF="node73.html">
<IMG WIDTH="26" HEIGHT="24" ALIGN="BOTTOM" BORDER="0" ALT="up" SRC="up.gif"></A>
<A NAME="tex2html2082"
HREF="node75.html">
<IMG WIDTH="63" HEIGHT="24" ALIGN="BOTTOM" BORDER="0" ALT="previous" SRC="prev.gif"></A>
<A NAME="tex2html2090"
HREF="node1.html">
<IMG WIDTH="65" HEIGHT="24" ALIGN="BOTTOM" BORDER="0" ALT="contents" SRC="contents.gif"></A>
<A NAME="tex2html2092"
HREF="node149.html">
<IMG WIDTH="43" HEIGHT="24" ALIGN="BOTTOM" BORDER="0" ALT="index" SRC="index.gif"></A>
<BR>
<B> Next:</B> <A NAME="tex2html2095"
HREF="node77.html">6.4 A 3-D perspective</A>
<B> Up:</B> <A NAME="tex2html2089"
HREF="node73.html">6. Cook-book</A>
<B> Previous:</B> <A NAME="tex2html2083"
HREF="node75.html">6.2 Image presentations</A>
  <B> <A NAME="tex2html2091"
HREF="node1.html">Contents</A></B>
  <B> <A NAME="tex2html2093"
HREF="node149.html">Index</A></B>
<!--End of Navigation Panel-->
<ADDRESS>
Paul Wessel
2001-04-18
</ADDRESS>
</BODY>
</HTML>
|