Subversion Repositories lagranto.arpege

Rev

Details | Last modification | View Log | RSS feed

Rev Author Line No. Line
3 michaesp 1
\documentclass[a4paper,10pt]{article}
2
\usepackage[utf8x]{inputenc}
3
 
4
\usepackage{graphicx}
5
\graphicspath{{/home/sprenger/lagranto/docu/tutorial/}}
6
 
7
\textwidth16cm
8
\textheight22.5cm
9
\oddsidemargin0.5cm
10
\evensidemargin0cm
11
\topmargin0.cm
12
\headsep0cm
13
\topskip-0.5cm
14
 
15
\title{Lagranto - Tutorial}
16
\author{Michael Sprenger and Heini Wernli}
17
 
18
\begin{document}
19
 
20
\maketitle
21
 
22
\section{Definition of the Problem}
23
 
24
In this tutorial we consider a simple way how to find airstreams which transport air from the surface into the upper troposphere and lower stratosphere, i.e. the UTLS region. More specifically, we intend
25
\begin{itemize}
26
\item[1)] to initialize starting positions over the North Atlantic and Europe (80 W to 20 E, 40 N to 80 N) for 00 UTC, 20 October 1989. The starting positions are horizontally equidistant with 80 km horizontal spacing and are set 100 hPa above ground level.
27
\item[2)] to calculate trajectories 72 hours forward in time.
28
\item[3)] to select trajectories which ascend to levels above 400 hPa within 48 hours and are found at starting time below 700 hPa.
29
\item[4)] to trace several meteorological fields along the trajectories: potential vorticity (PV), potential temperature (TH), relative (RH) and specific humidity (Q).
30
\item[4] to select subsamples of trajectories: a) those reaching the stratosphere; b) those travelling at most 2000 km; ... 
31
\item[5)] to show the densities of the trajectories on a geographical map.
32
\end{itemize}
33
 
34
\section{Meteorological Data}
35
 
36
The following dates and netCDF files are needed for a Lagranto calculation covering the time period from 00\,UTC 20 October 1989 to 00\,UTC 23 October 1989. 
37
\begin{verbatim}
38
> datelist stdout -create 19891020_00 19891023_00
39
 19891020_00
40
 19891020_06
41
 19891020_12
42
 ...
43
 19891022_06
44
 19891022_12
45
 19891022_18
46
 19891023_00
47
\end{verbatim}
48
 
49
\noindent
50
There are two different files involved, the P and the S files. Lagranto expects them to be in the running directory:
51
 
52
\begin{verbatim}
53
> ls -1
54
 P19891020_00 P19891020_06 P19891020_12 P19891020_18 
55
 P19891021_00 P19891021_06 P19891021_12 P19891021_18
56
 ... 
57
 S19891020_00 S19891020_06 S19891020_12 S19891020_18 
58
 S19891021_00 S19891021_06 S19891021_12 S19891021_18
59
 ...
60
\end{verbatim}
61
 
62
\noindent
63
The meteorological fields on the primary P files are at least: zonal wind (U, in m/s); meridional wind (V, in m/s); vertical wind (OMEGA, in Pa/s); 
64
surface pressure (PS, in hPa). Additional fields might be available on the P files, e.g. temperature (T), specific
65
humidity (Q),... Secondary fields can be saved in the S files, which must have the same grid structure as the P files.
66
In the example the following fields are saved on the S files: potential temperature (TH, in K); Ertel potential vorticity (PV, in pvu); relative humidity (in \%). 
67
Furthermore, the surface pressure (PS) is also saved on the S files; it must be exactly identical to the one in the P files.\\
68
 
69
\noindent
70
If the P and S files are stored at another place, they might be linked with a simple Shell script (in csh), usimg the command {\em datelist}:
71
\begin{verbatim}
72
> foreach date ( `datelist stdout -create 19891020_00 19891023_00` )
73
>   ln -s {SOURCE DIR}/P${date} {DEST DIR}/P${date}
74
>   ln -s {SOURCE DIR}/S${date} {DEST DIR}/S${date}
75
> end
76
\end{verbatim}
77
Note that the command {\em datelist} offers several options how to work with date list - creating, stepping through, comparing.
78
 
79
\section{Starting Positions}
80
 
81
In a first step the starting positions must be specified. To this aim a file {\em regionf} must be created with the definition of the region. The file is in the same directory as the meteorological data (section 1):
82
 
83
\begin{verbatim}
84
> more regionf 
85
"1 -80 20 40 80"
86
\end{verbatim}
87
 
88
\noindent
89
The first number specifies a region ID (here 1) and the other values are: west boundary (80 W), east boundary (20 E), southern boundary (20 N) and northern boundary (80 N). The starting positions are then created with
90
 
91
\begin{verbatim}
92
> create_startf 19891020_00 startf.2 'region.eqd(1,80) ...
93
                             ... @ level(100) @ hPa,agl' -changet
94
\end{verbatim}
95
 
96
\noindent
97
The starting positions are written to {\em startf.2}, i.e. in format 2, and cover the region 1 specified in {\em regionf}. The horizontal start points are equidistantly distributed with 80 km spacing, and they are all at 100 hPa above ground. All points refer to the starting date 00 UTC, 20 October 1989.
98
In total, 3750 starting positions are written to {\em startf.2}. The first few lines of the file look as follows:
99
 
100
\begin{verbatim}
101
> head -10 startf.2
102
Reference date 19891020_0000 / Time range       0 min
103
 
104
  time      lon     lat     p     level
105
---------------------------------------
106
 
107
   0.00   -79.61   40.45   862   100.000
108
   0.00   -78.66   40.45   860   100.000
109
   0.00   -77.71   40.45   873   100.000
110
   0.00   -76.76   40.45   886   100.000
111
   0.00   -75.82   40.45   893   100.000
112
\end{verbatim}
113
 
114
\noindent
115
The different columns are: time, longitude, latitude, pressure (in hPa) and level (in hPa,agl). Note that the 'same' starting file could have been created without a region file. In this case, the command would have been:
116
 
117
\begin{verbatim}
118
> create_startf 19891020_00 startf.2 'box.eqd(-80,20,40,80,80) ...
119
                             ... @ level(100) @ hPa,agl' -changet
120
\end{verbatim}
121
 
122
\noindent
123
However, note that the starting positions are not exactly the same as with the previous command: it is only guaranteed that the starting points are equidistantly distributed within the region.
124
 
125
 
126
\section{Trajectory Calculation}
127
 
128
In a next step, the trajectories are calculated 72 h forward in time, with starting date 00 UTC, 20 October 1989. The command is:
129
 
130
\begin{verbatim}
131
> caltra 19891020_00 19891023_00 startf.2 traj.4 -j
132
\end{verbatim}
133
 
134
\noindent
135
The taring positions are taken from{\em startf.2} and the output is written to {\em traj.4}. Furthermore, the jumping flag {\em -j} is set, i.e. if trajectories run into the ground they are lifted a little and allowed to move on.\\
136
 
137
\noindent
138
Note that the output file {\em traj.4} is not in ASCII format. To look at the file, use the command {\em trainfo}, for instance:
139
 
140
\begin{verbatim}
141
> trainfo traj.4 vars
142
time lon lat p 
143
 
144
> trainfo traj.4 dim
145
3750           13            4
146
 
147
> trainfo trai.4 startdate
148
19891020_0000
149
\end{verbatim}
150
 
151
\noindent
152
It is also possible to convert the trajectory file into ASCII format with the command {\em reformat}:
153
 
154
\begin{verbatim}
155
> reformat traj.4 traj.1
156
> more traj.1
157
Reference date 19891020_0000 / Time range    4320 min
158
 
159
  time      lon     lat     p
160
-----------------------------
161
 
162
   0.00   -79.61   40.45   862
163
   6.00   -80.57   43.23   791
164
  12.00   -82.23   45.89   782
165
  18.00   -84.94   47.07   744
166
\end{verbatim}
167
 
168
\noindent
169
The command {\em trafinfo} cna also be used to look at the trajectory tables, i.e. without a conversion to ASCII format. To this aim, use
170
 
171
\begin{verbatim}
172
> trainfo traj.4 list
173
\end{verbatim}
174
 
175
\noindent
176
Whereas the ASCII format is most convenient for visual inspection, it is the least compact format. In particular, if the output of {\em caltra} should be further processed, e.g. with {\em trace} or {\em select}, the binary formats should be used (see documentation for {\em reformat}).
177
 
178
\section{Pre-Selection of Trajectories}
179
 
180
Quite often, the position of the air parcels is sufficient to select trajectories. It is then most efficient to pre-select these trajectories and do the tracing of additional fields along the trajectories only on the pre-selected ones. In the example, airstreams should be identified which ascend from below 700 hPa at initial time to levels above 300 hPa. The ascent has to take place within 48 h. This selection can be achieved with the command:
181
 
182
\begin{verbatim}
183
> select traj.4 wcb.1 'GT:p:700:FIRST & LT:p(MIN):400:0 to 48'
184
\end{verbatim}
185
 
186
\noindent
187
The first criterion selects all trajectories for which the pressure (p) at the initial time (FIRST) is greater than (GT) 700 hPa. This criterison is logically AND-combined with the second criteron: consider all times between 0 and 48 h and take the minimum pressure, i.e. p(MIN), over this time interval; if the minium pressure is less than (LT) 400 hPa, the trajectory is selected. A sample trajectory looks like:
188
 
189
\begin{verbatim}
190
> trainfo wcb.1 list 
191
Reference date 19891020_0000 / Time range    4320 min
192
 
193
  time      lon     lat     p
194
-----------------------------
195
 
196
   0.00   -72.98   40.45   918
197
   6.00   -76.45   43.14   879
198
  12.00   -78.53   46.69   808
199
  18.00   -80.08   48.70   770
200
  24.00   -84.49   48.71   563
201
  30.00   -87.89   43.32   377
202
  36.00   -80.69   37.24   396
203
  42.00   -73.05   39.00   477
204
  48.00   -67.62   47.21   488
205
  54.00   -63.53   54.61   455
206
  60.00   -53.79   58.53   447
207
  66.00   -38.79   59.08   452
208
  72.00   -27.72   55.51   493
209
\end{verbatim}
210
 
211
\noindent
212
In total, 99 trajectories are selected. Further Lagrangian selection criteria might be reasonable. For instance, it could be of interest whether the air parcels are far away from their initial position after 48 h:
213
 
214
\begin{verbatim}
215
> select traj.4 wcb.1 'GT:p:700:FIRST & LT:p(MIN):400:0 to 48 & GT:DIST0:5000:48'
216
\end{verbatim}
217
 
218
\noindent
219
The last criterion test whether the spherical distance from the initial position exceeds 5000 km at 48 h (only met by 4 trajectories). Note that the field {\em DIST0} is not available on the trajectory file {\em traj.4}, but is implicitely calculated. \\
220
 
221
\noindent
222
Similarly, it can be tested whether a trajectory passes through a target region (e.g. 20E-30E,50N-60N). Such a region might be defined in the region file {\em regionf}:
223
 
224
\begin{verbatim}
225
> more regionf 
226
# Starting positions
227
"1 -80  20 40 80"
228
# Target region
229
"2 20 30 50 60"
230
\end{verbatim}
231
 
232
\noindent
233
The call to {\em select} now looks as follows:
234
 
235
\begin{verbatim}
236
> select traj.4 wcb.1 'GT:p:700:FIRST & LT:p(MIN):400:0 to 48 & ...
237
                       ... TRUE:INREGION:2:42 to 54(ANY)'
238
\end{verbatim}
239
 
240
\noindent
241
The last criterion is interpreted in the following way: consider the times from 42 h to 48 h (42 to 54) and check whether a trajectory is at any of these times (ANY) in the traget region 2, as specified in {\em regionf}. A sample trajectory is given below:
242
 
243
\begin{verbatim}
244
> more wcb.1
245
Reference date 19891020_0000 / Time range    4320 min
246
 
247
  time      lon     lat     p
248
-----------------------------
249
 
250
   0.00   -44.56   40.45   914
251
   6.00   -41.22   38.95   885
252
  12.00   -37.50   37.59   849
253
  18.00   -33.54   36.76   823
254
  24.00   -29.45   36.45   770
255
  30.00   -24.56   37.44   622
256
  36.00   -18.89   41.17   429
257
  42.00   -10.48   49.48   349
258
  48.00     8.75   57.02   355
259
  54.00    28.74   56.97   354
260
  60.00    37.14   53.71   320
261
  66.00    38.99   49.75   334
262
  72.00    37.57   45.94   349
263
\end{verbatim}
264
 
265
\section{Tracing Meteorological Fields}
266
 
267
Meteorological fields can be traced along the trajectories with the command {\em trace}. Most often, a list of fields to trace will be listed in a file {\em tracevars}:
268
 
269
\begin{verbatim}
270
> more tracevars 
271
PS    1.  0 P
272
Q  1000.  0 P
273
TH    1.  0 S
274
RH    1.  1 *
275
\end{verbatim}
276
 
277
\noindent
278
The following fields are to be traced: surface pressure (PS), specific humidity (Q), potential temperature (TH) and relative humidity (RH). PS and Q are available on the P files and need not to be calculated; Q will be scaled by a factor 1000 to convert from Kg/kg to g/kg. TH is found on the S file and need not to be calculated. Finally, RH is found neither on the P nor on the S file and must becomputed - the flag 1 in the third column.\\
279
 
280
\noindent
281
With the {\em tracevars} file ready, the tracing is started with:
282
 
283
\begin{verbatim}
284
> trace wcb.1 wcb.1
285
\end{verbatim}
286
 
287
\noindent
288
Note that the input and output file are allowed to have the same name. The following table shows the first few lines of the new trajectory file:
289
 
290
\begin{verbatim}
291
> more wcb.1 
292
Reference date 19891020_0000 / Time range    4320 min
293
 
294
  time      lon     lat     p        PS         Q        TH        RH
295
---------------------------------------------------------------------
296
 
297
   0.00   -44.56   40.45   914  1014.093     9.664   294.921    86.667
298
   6.00   -41.22   38.95   885  1012.965     8.380   296.549    77.992
299
  12.00   -37.50   37.59   849  1014.056     8.565   299.122    81.321
300
  18.00   -33.54   36.76   823  1012.074     8.720   300.798    85.407
301
  24.00   -29.45   36.45   770  1012.254     7.666   303.487    85.262
302
\end{verbatim}
303
 
304
\noindent
305
If it is later found that additional fields should be traced, this can be done with a new {\em tracevars} file or for a single field with (the second number being the scaling factor):
306
 
307
\begin{verbatim}
308
> trace wcb.1 wcb.1 -f PV 1.
309
\end{verbatim}
310
 
311
\noindent
312
Note that several fields are allowed for online computation, i.e. with the computation flag set in {\em tracevars}. This is convenient for interactive mode and for few trajectories. However, if tracing is needed for many trajectories, a pre-calculation and saving on the S files is much more efficient! A list of fields for online computation is found in the reference guide for {\em trace}. \\
313
 
314
\noindent
315
It is also possible to trace the surrounding of a trajectory position, i.e. to get for instance not the temperature at the air parcel's position, but at 50 hPa above or below it. This is done with:
316
 
317
\begin{verbatim}
318
> trace wcb.1 wcb.1 -f T:-50HPA 1.
319
> trace wcb.1 wcb.1 -f T:+50HPA 1.
320
\end{verbatim}
321
 
322
\noindent
323
Finally, if it is decided that a field is no longer needed in the trajectory file, or if it has to be corrected, it is possible to remove columns from the trajectory file. This can be achieved with {\em extract} - for instance if only PS, TH and RH should be kept:
324
 
325
\begin{verbatim}
326
> extract wcb.1 wcb.1 -var PS TH RH
327
\end{verbatim}
328
 
329
\noindent
330
Note that {\em extract} can also be used to extract different times or starting positions (see the reference documentation). 
331
 
332
 
333
\section{Final Selection of Trajectories}
334
 
335
In this section the selection of trajectories should be refined, i.e. it is not only based on positional information but also on further meteorological parameters. We look at several questions:\\
336
 
337
\begin{itemize}
338
\item[a)] Is there a trajectory which reaches saturation ($RH>99$\%)? The trajectories should be saved in a new trajectory file.
339
\begin{verbatim}
340
> select wcb.1 sat.1 'GT:RH:99:0 to 72(ANY)'
341
> more sat.1
342
Reference date 19891020_0000 / Time range    4320 min
343
 
344
  time      lon     lat     p        PS         Q        TH        RH        PV   
345
---------------------------------------------------------------------------------
346
   0.00   -72.98   40.45   918  1018.161     9.503   292.020   100.722     0.920    
347
   6.00   -76.45   43.14   879   986.319     8.723   294.933    92.837     1.101    
348
  12.00   -78.53   46.69   808   972.550     7.621   297.875    97.737     0.794    
349
  18.00   -80.08   48.70   770   973.957     6.147   297.914    97.912     1.078    
350
  24.00   -84.49   48.71   563   962.279     2.327   307.548    87.923     1.034    
351
  30.00   -87.89   43.32   377   977.415     0.319   314.210    65.759     0.108    
352
  36.00   -80.69   37.24   396   939.606     0.303   312.705    52.845     0.323    
353
  42.00   -73.05   39.00   477  1013.693     0.298   314.614    16.248     0.309    
354
  48.00   -67.62   47.21   488   970.025     0.442   312.975    23.890     0.463    
355
  54.00   -63.53   54.61   455   950.011     0.386   313.047    30.182     0.479    
356
  60.00   -53.79   58.53   447  1007.039     0.392   311.951    36.578     0.487    
357
  66.00   -38.79   59.08   452  1006.532     0.319   311.316    29.286     0.443    
358
  72.00   -27.72   55.51   493  1009.871     0.279   311.428    15.950     0.513    
359
\end{verbatim}
360
 
361
\item[b)] Get a list of all trajectories which pass through a circle around 20\,W/40\,N and radius 500\,km.
362
\begin{verbatim}
363
> select wcb.1 indlist 'TRUE:INCIRCLE:-20,40,500:ALL(ANY)' -index
364
> more indlist
365
            4
366
            5
367
            6
368
           11
369
           12
370
           13
371
           14
372
           19
373
           20
374
           21
375
           22
376
           47
377
\end{verbatim}
378
Hence, the trajectories 4,5,... pass through the circle. The trajectories themselves can be extracted in a second step with
379
\begin{verbatim}
380
> extract wcb.1 pass.1 -index indlist
381
\end{verbatim}
382
where now the selected trajectories are written to the trajectory file {\em pass.1}.
383
 
384
\item[c)] Select all trajectories which are in the stratosphere after 48 h. The dynamical tropopause is defined as the 2-PVU isosurface?
385
\begin{verbatim}
386
> select wcb.1 out 'GT:PV:2:48 to 72(ALL) & LT:P:500:48 to 72(ALL)' 
387
\end{verbatim}
388
The second criterion guarantees that the air parcel is at a height above 500 hPa; indeed, low-level high-PV regions might mimick a stratosphere, although they are of diabatic origin.  
389
 
390
\item[d)] Select all trajectories which are within 2000 km distance of their starting position after 72 h.
391
\begin{verbatim}
392
> select wcb.1 sel.1 'LT:DIST0:2000:LAST'
393
\end{verbatim}
394
Note that the fields {\em DIST0} needs not to be available on the trajectory file - it is calculated during the selection. {\em DIST0} refers to the spherical distance (in km) from the strting position. If the path length (in km) is needed, {\em DIST} can be used instead.
395
 
396
\item[e)] How many trajectories ascend more than 550 hPa between 12 h and 54 h? We are only interested in the number of selected trajectories.
397
\begin{verbatim}
398
> select wcb.1 count 'GT:P(DIFF):550:12,54' -count
399
> more count
400
3
401
\end{verbatim}
402
 
403
\item[f)] We would like to select all trajectories which reach potential vorticity (PV) greather than 2 PVU at levels above 500 hPa. In a first attempt, this might be accomplished with the criterion
404
\begin{verbatim}
405
> select wcb.1 wcb.1 'GT:PV:2:ALL(ANY) & LT:p:500:ALL(ANY)' 
406
\end{verbatim} 
407
But note that this is not exactly what we want - the first criterion might be fulfilled at a time 48 h, for instance, whereas the second criterion is fulfilled at another time, say 72 h. Hence they are not both fulfilled at the same time! A way around this problem is possible if a TRIGGER column is used to mark the two events. The original trajectory file looks as follows:
408
\begin{verbatim}
409
> more wcb.1
410
Reference date 19891020_0000 / Time range    4320 min
411
 
412
  time      lon     lat     p        PS          RH        PV
413
------------------------------------------------------------
414
 
415
   0.00   -19.56   46.94   905  1005.242     83.514     0.291
416
   6.00   -14.72   48.17   892   999.182     88.325     0.242
417
  12.00   -10.58   50.53   862   993.145     97.718     0.293
418
  18.00    -7.22   53.02   792   972.076     99.216     0.738
419
  24.00    -3.71   55.89   724   956.135     93.218     1.076
420
  30.00    -0.19   58.87   629   971.334     70.088     1.076
421
  36.00     1.46   61.62   452   966.406     66.056     0.558
422
  42.00     0.01   62.49   328   977.209     65.319     1.754
423
  48.00    -1.54   63.41   313   983.930     56.822     2.727
424
  54.00    -3.59   64.77   322   984.627     58.328     1.874
425
  60.00    -9.91   66.07   323   988.185     57.894     2.052
426
  66.00   -20.91   66.02   316   976.560     57.989     2.565
427
  72.00   -28.89   66.19   319  1007.175     54.477     2.693
428
\end{verbatim}
429
Then we mark the two events with a TRIGGER:
430
\begin{verbatim}
431
> select wcb.1 wcb.1 'GT:PV:2:1(TRIGGER) & LT:p:500:2(TRIGGER)' -trigger
432
\end{verbatim}
433
The first criterion (PV) gets the trigger 1 (in binary system 01), the second one (pressure) get the trigger 2 (in binary system 10). If both criteria are fulfilled, the trigger column becomes 3, which corresponds in the binary system to 11 - i.e. each flag corresponds to a bit in the trigger value. With the option '-trigger' the trigger column is written to the output trajectory file:
434
\begin{verbatim}
435
> more wcb.1
436
Reference date 19891020_0000 / Time range    4320 min
437
 
438
  time      lon     lat     p        PS         RH        PV     TRIGGER
439
-------------------------------------------------------------------------
440
 
441
   0.00   -19.56   46.94   905  1005.242    83.514     0.291      0.000
442
   6.00   -14.72   48.17   892   999.182    88.325     0.242      0.000
443
  12.00   -10.58   50.53   862   993.145    97.718     0.293      0.000
444
  18.00    -7.22   53.02   792   972.076    99.216     0.738      0.000
445
  24.00    -3.71   55.89   724   956.135    93.218     1.076      0.000
446
  30.00    -0.19   58.87   629   971.334    70.088     1.076      0.000
447
  36.00     1.46   61.62   452   966.406    66.056     0.558      2.000
448
  42.00     0.01   62.49   328   977.200    65.319     1.754      2.000
449
  48.00    -1.54   63.41   313   983.930    56.822     2.727      3.000
450
  54.00    -3.59   64.77   322   984.627    58.328     1.874      2.000
451
  60.00    -9.91   66.07   323   988.185    57.894     2.052      3.000
452
  66.00   -20.91   66.02   316   976.560    57.989     2.565      3.000
453
  72.00   -28.89   66.19   319  1007.175    54.477     2.693      3.000
454
\end{verbatim}
455
Now the selection can be achieved by refering to the TRIGGER column:
456
\begin{verbatim}
457
> select wcb.1 wcb.1 'ALL:TRIGGER:1,2:ALL(ANY)'
458
\end{verbatim}
459
This selection means that the trigger values 1 and 2 must be set - the operator ALL (first term) guaranteeing that all selected triggers are set. The time specification ALL(ANY) is as before, i.e. the check is performed for all times and he criterion must be fulfilled at any of these times.\\
460
 
461
\item[g)] Select all trajectories which pass at time 60 h over Switzerland! The coordinates of the Swiss boundary are listed in a file {\em borders.dat}:
462
\begin{verbatim}
463
> more borders.dat
464
8.55 47.45
465
7.942863 46.002075
466
7.949024 46.001195
467
7.956945 46.000022
468
7.984226 46.000022
469
7.989800 46.001489
470
8.000068 46.007356
471
8.011508 46.018503
472
...
473
\end{verbatim}
474
The first line is a point (longitude, latitude) within Switzerland (Zurich), the other lines define the boundary of Switzerland (as 1373 points). With this polygon file, the selection command becomes
475
\begin{verbatim}
476
> select wcb.1 out.1 'TRUE:INPOLYGON:borders.dat:60'
477
\end{verbatim} 
478
Note that in a criterion only one polygon can be specified.
479
 
480
\item[g)] New criteria can easily be implemented into the Fortran code; to this aim the file {\em special.f} in directory {\em select} must be edited. The following example shows the implementation of an identification for warm conveyor belts (WCB). The calling sequence for the criterion is {\em SPECIAL:WCB:ascent,first,last}, the air stream must ascend at least {\em ascent} hPa between time {\em first} and time {\em last}. The corresponding Fortran looks as follows:
481
\begin{small}
482
\begin{verbatim}
483
> more select/special.f
484
 
485
      SUBROUTINE special (flag,cmd,tra,ntim,ncol,
486
     >                    vars,times,param,nparam)
487
 
488
c     ***************************************************************************
489
c     *                                                                         *
490
c     * OUTPUT:  flag           -> 1 if trajectory is selected, 0 if not        *
491
c     *                                                                         *
492
c     * INPUT:   cmd            <- command string (e.g. WCB)                    *
493
c     *          tra(ntim,ncol) <- single trajectory: indices time,column       *
494
c     *          ntim           <- number of times                              *
495
c     *          ncol           <- number of columns (including time,lon,lat,p) *
496
c     *          vars(ncol)     <- names of columns                             *
497
c     *          times(ntim)    <- List of times
498
c     *          param(nparam)  <- parameter values                             *
499
c     *          nparam         <- number of parameters                         *
500
c     *                                                                         *
501
c     ***************************************************************************
502
 
503
      implicit none
504
 
505
c     ---------------------------------------------------------------------------
506
c     Declaration of subroutine parameters
507
c     ---------------------------------------------------------------------------
508
 
509
      integer       flag           ! Boolean flag whether trajectory is selected
510
      character*80  cmd            ! Command string
511
      integer       ntim,ncol      ! Dimension of single trajectory
512
      real          tra(ntim,ncol) ! Single trajectory
513
      character*80  vars(ncol)     ! Name of columns
514
      real          times(ntim)    ! List of times
515
      integer       nparam         ! # parameters
516
      real          param(nparam)  ! List of parameters
517
 
518
c     ---------------------------------------------------------------------------
519
c     Declaration of local variables
520
c     ---------------------------------------------------------------------------
521
 
522
      integer       i
523
      integer       ip,i0,i1
524
 
525
c     --------------------------------------------------------------------------  %)
526
c     SPECIAL:WCB:ascent,first,last                                               %)
527
c         : Detect Warm Conveyor Belts (WCB); the air stream must ascend at least %)
528
c         : <ascent=param(1)> hPa between the two times <first=param(2)> and      %)
529
c         : <last=param(3)>. Note, the lowest pressure is allowed to occur at any %)
530
c         : time between <first> and <last>.                                      %)
531
c     --------------------------------------------------------------------------- %)
532
 
533
      if ( cmd.eq.'WCB' ) then
534
 
535
c        Reset the flag for selection
536
         flag = 0
537
 
538
c        Pressure is in the 4th column
539
         ip = 4
540
 
541
c        Get indices for times <first> and <last>
542
         i0 = 0
543
         i1 = 0
544
         do i=1,ntim
545
            if ( param(2).eq.times(i) ) i0 = i
546
            if ( param(3).eq.times(i) ) i1 = i
547
         enddo
548
         if ( (i0.eq.0).or.(i1.eq.0) ) then
549
            print*,' ERROR: invalid times in SPECIAL:WCB... Stop'
550
            stop
551
         endif
552
 
553
c        Check for ascent 
554
         do i=i0+1,i1
555
            if ( ( tra(1,ip)-tra(i,ip) ) .gt. param(1) ) flag = 1
556
         enddo
557
 
558
      endif
559
 
560
c     ---------------------------------------------------------------------------
561
 
562
      end
563
\end{verbatim}
564
\end{small}
565
 
566
\end{itemize}
567
 
568
\section{Trajectory Densities}
569
Single trajectories can be visualised e.g. with Matlab or wit NCL (see template scripts in the Lagranto folder). If many trajectories should be visualised instead, it is much more convenient to show trajectory densities. The easiiest way to get trajectory densities is:
570
\begin{verbatim}
571
> density wcb.1 densisty
572
> ncview density 
573
\end{verbatim}
574
This will project the trajectories in the trajectory file {\em wcb.1} onto a global longitude/latitude grid with 1 degree horizontal resolution. A filter radius of 100 km will be used
575
 
576
\includegraphics[width=0.9\textwidth]{screen1.ps}
577
 
578
\noindent
579
The CF-netCDF file contains several fields: a) the number of trajectory points associated to each grid point (COUNT); b) the residence time of the trajectories (in hours) associated to each grid point - the residence time being the time a trajectory stays at a certain grid cell (RESIDENCE); c) the area (im $km^2$) associated with each grid cell. The area allows to change the unit of the gridded trajectory from counts per grid point to counts per $km^2$.\\
580
 
581
\noindent
582
Often the trajectories do not spread over the whole globe; then the subdomain can be specified with
583
\begin{verbatim}
584
> density traj.1 density -latlon 300 150 -100 10 0.5 0.5 -create
585
\end{verbatim}
586
where the new grid has 300x150 grid points in zonal and meridional direction with south-eastern corner at 100\,W/10\,S and resolution of 0.5 degree in zonal and 0.5 degree in meridional direction. The flag {\em create} forces the netCDF file to be created anew, even if it already exists.
587
 
588
\includegraphics[width=0.9\textwidth]{screen2.ps}
589
 
590
\noindent
591
It is also possible to re-parameterise the trajectories before they are gridded. For instance, the following command interpolates the positions to a 1-h time interval and then performs the grissing. This results in a smoother density plot:
592
\begin{verbatim}
593
> density traj.1 density -latlon 300 150 -100 10 0.5 0.5 -create -interp 1 h
594
\end{verbatim}
595
 
596
\includegraphics[width=0.9\textwidth]{screen3.ps}
597
 
598
\noindent
599
In addition to a gridding of the complete trajectories, the single trajectory times can be gridded. 
600
\begin{verbatim}
601
> density traj.1 density -create -time  0.00 -create
602
> density traj.1 density -create -time  6.00 
603
> density traj.1 density -create -time 12.00 
604
> density traj.1 density -create -time 18.00 
605
> density traj.1 density -create -time 24.00 
606
\end{verbatim}
607
 
608
\noindent
609
In the previous figures, the density of the trajectories was determined - i.e. the number of trajectories associated with a grid point or the residence associated with a grid cell was determined. In addition to this most basic information, it is also possible to perform a gridding of any trajectory field. For instance, the trajectory file contains
610
 
611
\begin{verbatim}
612
Reference date 19891020_0000 / Time range    4320 min
613
 
614
  time      lon     lat     p
615
-----------------------------
616
 
617
   0.00   -79.61   40.45   862
618
   6.00   -80.57   43.23   791
619
  12.00   -82.23   45.89   782
620
  18.00   -84.94   47.07   744
621
\end{verbatim}
622
 
623
\noindent
624
and we would like to know at what height the trajectories typically (in the mean) are at a specific grid point. Then we would grid the pressure ``p'' instead of the position:
625
 
626
\begin{verbatim}
627
> density traj.1 density -create -latlon 300 150 -100 10 0.5 0.5 -field p -time 0.00
628
> density traj.1 density  -field p -time 24.00
629
> density traj.1 density  -field p -time 48.00
630
\end{verbatim}
631
 
632
\noindent
633
The following figures show the gridded pressure at time 0.00 and 24.00 h; note that the trajectories were initialised 100 hPa above ground.
634
 
635
\noindent
636
\begin{center}
637
\includegraphics[width=0.8\textwidth]{screen4.ps}
638
\\
639
\includegraphics[width=0.8\textwidth]{screen5.ps}
640
\end{center}
641
 
642
\section{Interface Script}
643
 
644
\subsection{Start from local directory}
645
 
646
In addition to the programs {\em create\_startf}, {\em caltra}, {\em trace}, {\em select}, Lagranto offers a ``master'' script which combines the call to the individual programs into one single call. For instance,
647
\begin{verbatim}
648
> lagranto local 19891020_00 19891024_18 startf nil -changet
649
\end{verbatim}
650
will start a Lagranto run starting from 00,UTC 20 Octiober 2010 to  18,UTC 24 Octiober 2010, based upon the starting positions in the file {\em startf}. No selection of trajectories is applied, as specified with the flag {\em nil}, and the times on the netCDF P and S files are set relative to the starting date prior  to the Lagranto run. Finally, {\em local} means that all input files are expected in the directory where Lagranto was called.\\
651
 
652
\noindent
653
The output for the above Lagranto call is saved in a newly created directory, which is located in the calling directory:
654
\begin{verbatim}
655
> ls -l ntr_19891020_00_f114_local_startf_nil/
656
-rw-r--r-- 1 michaesp wheel 5328945 2011-03-21 14:03 lsl_19891020_00
657
-rw-r--r-- 1 michaesp wheel   68195 2011-03-21 14:03 runscript.logfile
658
-rwxr--r-- 1 michaesp wheel    1025 2011-03-21 14:02 runscript.sh*
659
\end{verbatim}
660
The three different files are:
661
\begin{itemize}
662
\item[a)]  {\bf lsl\_19891020\_00}: the output trajectory file -the file name starts with {\em lsl} and contains the starting date of the Lagranto run:
663
\begin{verbatim}
664
> more lsl_19891020_00
665
Reference date 19891020_0000 / Time range    6840 min
666
 
667
  time      lon     lat     p        PS         Q        TH        RH
668
---------------------------------------------------------------------
669
 
670
   0.00   -79.61   40.45   862   961.659     6.434   290.838    97.722
671
   6.00   -80.57   43.23   791   980.984     5.334   293.824    98.773
672
   ...
673
\end{verbatim}
674
Note that additional fields have been traced along the trajectories, as specified in the tracing file {\em tracevars}:
675
\begin{verbatim}
676
> more tracevars
677
PS    1.  0 P
678
Q  1000.  0 P
679
TH    1.  0 S
680
RH    1.  1 *
681
\end{verbatim}
682
\item[b)] {\bf runscript.logfile}: a log file with all status and error information of the Lagranto run. If the flag {\em -log} is set in a Lagranto call, the log will be written to screen.
683
\item[c)] {\bf runscript.sh}: the calling script for the programs {\em create\_startf}, {\em caltra}, {\em trace} and {\em select}. The basic idea of {\em lagranto} is to create the output directory, to prepare all netCDF and other files in this output directory and to create a Shell script with name {\em runscript.sh}. If all these preparations were successfull, Lagranto will change into the output directory and launch {\em runscript.sh}. The {\em runscript.sh} for the previous Lagranto call looks as follows:
684
\begin{small}
685
\begin{verbatim}
686
#!/bin/csh
687
#
688
#----- Calling command
689
#
690
# lagranto local 19891020_00 19891024_18 startf nil -changet -log 
691
#
692
#----- Output file    
693
#
694
# lsl_19891020_00                
695
#
696
#------ Abort if no startf is available
697
#
698
if ( ! -f startf ) then
699
  echo " ERROR: no start file available .... Stop"
700
  exit 1
701
endif
702
#
703
#------ Remove existing trajectory files
704
#
705
if ( -f lsl_19891020_00.4 ) then
706
  \rm -f lsl_19891020_00.4
707
endif
708
if ( -f lsl_19891020_00 ) then
709
  \rm -f lsl_19891020_00
710
endif
711
#
712
#------ Run <caltra>
713
#
714
/home/sprenger/lagranto//bin/caltra.sh 19891020_00 19891024_18 startf lsl_19891020_00.4
715
#
716
#------ Abort if caltra was not successful
717
#
718
if ( ! -f lsl_19891020_00.4 ) then
719
  echo " ERROR: caltra failed .... Stop"
720
  exit 1
721
endif
722
#
723
#------ Run <trace>
724
#
725
/home/sprenger/lagranto//bin/trace.sh lsl_19891020_00.4 lsl_19891020_00 -v tracevars
726
#
727
#------ Abort if trace was not successful
728
#
729
if ( ! -f lsl_19891020_00 ) then
730
  echo " ERROR: trace failed .... Stop"
731
  exit 1
732
endif
733
\end{verbatim}
734
\end{small}
735
\end{itemize}
736
 
737
\noindent
738
Note that you are free to change to the output directory and manually launch {\em runscript.sh}, possibly after having modified it to your needs. This way of working is uspported by the optional flag {\em -prep} which will only prepare all files and then changes to the output directory:
739
\begin{verbatim}
740
> lagranto local 19891020_00 19891024_18 startf nil -changet -prep
741
\end{verbatim}
742
At the end of this call you will be asked to change to the output directory, which -after having agreed- will open a new {\em xterm} window. Note that you can always easily change to a output directory by calling
743
\begin{verbatim}
744
> lagranto -open local 
745
\end{verbatim}
746
If several trajectory runs are available in the local directory, you are asked to select one. Often, you would like to see the outcome of a run without changing to the output directory. This is most easily accomplished with the following call:
747
\begin{verbatim}
748
> lagranto -show local 
749
\end{verbatim}
750
 
751
\subsection{Start from case directory}
752
 
753
In this calling sequence, for instance
754
\begin{verbatim}
755
> lagranto tutorial 19891020_00 19891024_18 startf nil -changet
756
\end{verbatim}
757
the input files are not expected in the local directory, but are specified by means of a case identifier. For instance, a case has the identifier {\em tutorial}. Then Lagranto will expect the input netCDF P and S files to be located in
758
\begin{verbatim}
759
> ls -l ${HOME}/cdf/tutorial
760
/home/sprenger/cdf/tutorial/P19891020_00
761
/home/sprenger/cdf/tutorial/P19891020_06
762
/home/sprenger/cdf/tutorial/P19891020_12
763
/home/sprenger/cdf/tutorial/P19891020_18
764
/home/sprenger/cdf/tutorial/P19891021_00
765
/home/sprenger/cdf/tutorial/P19891021_06
766
/home/sprenger/cdf/tutorial/P19891021_12
767
/hom
768
\end{verbatim}
769
and all the other input files (starting positions, tracing file, region file, polygon specification) are expected in 
770
\begin{verbatim}
771
> ls -l ${HOME}/tra/tutorial
772
startf
773
tracevars
774
\end{verbatim}
775
The output of the trajectory calculation will be written to the following output directory, where now the case identifier {\em tutorial} is part of the directory name:
776
\begin{verbatim}
777
> cd /home/michaesp/tra/tutorial/ntr_19891020_00_f114_tutorial_startf_nil
778
> ls -1
779
-rw-r--r-- 1 michaesp wheel 5328945 2011-03-21 14:03 lsl_19891020_00
780
-rw-r--r-- 1 michaesp wheel   68195 2011-03-21 14:03 runscript.logfile
781
-rwxr--r-- 1 michaesp wheel    1025 2011-03-21 14:02 runscript.sh*
782
\end{verbatim}
783
All other aspects are identical to the ones described in the previous section. 
784
 
785
 
786
\section{Installation}
787
 
788
In this section you will find some hints how to install Lagranto on a Linux platform. Everthing is handled with the installation script {\em install.sh} which comes with the Lagranto distribution:
789
\begin{verbatim}
790
> install.csh 
791
install.sh [lib|core|goodies|links|all]
792
\end{verbatim}
793
The installation should proceed in several distinct steps:
794
\begin{itemize}
795
\item[a)] Find the place of the netCDF (http://www.unidata.ucar.edu/software/netcdf/) installation on your system - note that the netCDF comes as a pre-compiled package for many Linux distributions and most often can be installed with the Linux software management. Define an environmental variable {\em NETCDF} which directs to your installation, e.g.
796
\begin{verbatim}
797
> setenv NETCDF /usr/local/netcdf/
798
\end{verbatim}
799
\item[b)] Set the environmental variable {\em LAGRANTO} to the place where you have stored the Lagranto source code and include Lagranto in your search path. In {\em csh} this might look as follows:
800
\begin{verbatim}
801
> setenv LAGRANTO   /home/michaesp/lagranto/
802
> set isLAGRANTO=`echo $PATH | grep $LAGRANTO | wc -l`
803
> if ( $isLAGRANTO == 0 ) then
804
>   setenv PATH $LAGRANTO/bin:${PATH}
805
> endif
806
\end{verbatim}
807
You might include these statements also in your {\em .cshrc} file. If successful, you will then be able to open the lagranto help, e.g. with
808
\begin{verbatim}
809
> lagrantohelp
810
\end{verbatim} 
811
\item[c)] Install the different components of Lagranto and create links - proceed step by step to ensure that each one was successfully completed:
812
\begin{verbatim}
813
> install lib
814
> install core
815
> install goodies
816
> install links
817
\end{verbatim} 
818
\end{itemize}
819
Lagranto should now be ready to run! As a next step you might want to consider the tutorial, which can be invoked with the command:
820
\begin{verbatim}
821
> lagrantohelp tutorial
822
\end{verbatim}
823
If you are familiar with the most basic aspects of Lagranto, please refer to the reference guide which enlists all options of Lagranto:
824
\begin{verbatim}
825
> lagrantohelp refernce
826
\end{verbatim}
827
The contents of the reference guide can also be called from within the Linux shell, e.g. the documentation of {\em caltra} can be seen in man page format with:
828
\begin{verbatim}
829
> lagrantohelp caltra 
830
\end{verbatim}
831
 
832
\end{document}