source: branches/petra/src/readcommand.f90 @ 37

Last change on this file since 37 was 37, checked in by pesei, 9 years ago

Wet dep quick fix and other small changes. Wet depo quick fix not final yet.

File size: 22.9 KB
Line 
1!**********************************************************************
2! Copyright 1998-2015                                                 *
3! Andreas Stohl, Petra Seibert, A. Frank, Gerhard Wotawa,             *
4! Caroline Forster, Sabine Eckhardt, John Burkhart, Harald Sodemann   *
5!                                                                     *
6! This file is part of FLEXPART.                                      *
7!                                                                     *
8! FLEXPART is free software: you can redistribute it and/or modify    *
9! it under the terms of the GNU General Public License as published by*
10! the Free Software Foundation, either version 3 of the License, or   *
11! (at your option) any later version.                                 *
12!                                                                     *
13! FLEXPART is distributed in the hope that it will be useful,         *
14! but WITHOUT ANY WARRANTY; without even the implied warranty of      *
15! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the       *
16! GNU General Public License for more details.                        *
17!                                                                     *
18! You should have received a copy of the GNU General Public License   *
19! along with FLEXPART.  If not, see <http://www.gnu.org/licenses/>.   *
20!**********************************************************************
21
22subroutine readcommand
23
24  !*****************************************************************************
25  !                                                                            *
26  !     This routine reads the user specifications for the current model run.  *
27  !                                                                            *
28  !     Author: A. Stohl                                                       *
29  !                                                                            *
30  !     18 May 1996                                                            *
31  !
32  !     HSO, 1 July 2014: Added optional namelist input
33  !     PS 2/2015: add ldep_incr as optional command input
34  !                make new parameters (last 3) optional for bwd compatibility
35  !                                                                            *
36  !*****************************************************************************
37  !                                                                            *
38  ! Variables:                                                                 *
39  ! bdate                beginning date as Julian date                         *
40  ! ctl                  factor by which time step must be smaller than        *
41  !                      Lagrangian time scale                                 *
42  ! ibdate,ibtime        beginnning date and time (YYYYMMDD, HHMISS)           *
43  ! ideltas [s]          modelling period                                      *
44  ! iedate,ietime        ending date and time (YYYYMMDD, HHMISS)               *
45  ! ifine                reduction factor for vertical wind time step          *
46  ! outputforeachrel     for forward runs it is possible either to create      *
47  !                      one outputfield or several for each releasepoint      *
48  ! iflux                switch to turn on (1)/off (0) flux calculations       *
49  ! iout                 1 for conc. (residence time for backward runs) output,*
50  !                      2 for mixing ratio output, 3 both, 4 for plume        *
51  !                      trajectory output, 5 = options 1 and 4                *
52  ! ipin                 1 continue simulation with dumped particle data, 0 no *
53  ! ipout                0 no particle dump, 1 every output time, 3 only at end*
54  ! itsplit [s]          time constant for particle splitting                  *
55  ! loutaver [s]         concentration output is an average over loutaver      *
56  !                      seconds                                               *
57  ! loutsample [s]       average is computed from samples taken every [s]      *
58  !                      seconds                                               *
59  ! loutstep [s]         time interval of concentration output                 *
60  ! lsynctime [s]        synchronisation time interval for all particles       *
61  ! lagespectra          switch to turn on (1)/off (0) calculation of age      *
62  !                      spectra                                               *
63  ! lconvection          value of either 0 and 1 indicating mixing by          *
64  !                      convection                                            *
65  !                      = 0 .. no convection                                  *
66  !                      + 1 .. parameterisation of mixing by subgrid-scale    *
67  !                              convection = on                               *
68  ! lsubgrid             switch to turn on (1)/off (0) subgrid topography      *
69  !                      parameterization                                      *
70  ! method               method used to compute the particle pseudovelocities  *
71  ! mdomainfill          1 use domain-filling option, 0 not, 2 use strat. O3   *
72  !                                                                            *
73  ! ldep_incr            .true.  for incremental,
74  !                      .false. for accumulated deposition output (DEFAULT)
75  ! Constants:                                                                 *
76  ! unitcommand          unit connected to file COMMAND                        *
77  !                                                                            *
78  !*****************************************************************************
79
80  use par_mod
81  use com_mod
82
83  implicit none
84
85  real (kind=dp) :: juldate
86  character (len=50) :: line
87  logical :: old
88  integer :: ireaderror, icmdstat
89
90  namelist /command/ &
91  ldirect, &
92  ibdate,ibtime, &
93  iedate,ietime, &
94  loutstep, &
95  loutaver, &
96  loutsample, &
97  itsplit, &
98  lsynctime, &
99  ctl, &
100  ifine, &
101  iout, &
102  ipout, &
103  lsubgrid, &
104  lconvection, &
105  lagespectra, &
106  ipin, &
107  ioutputforeachrelease, &
108  iflux, &
109  mdomainfill, &
110  ind_source, &
111  ind_receptor, &
112  mquasilag, &
113  nested_output, &
114  linit_cond, &
115  lnetcdfout, &
116  surf_only, &
117  ldep_incr
118
119  ! Set default values for namelist command
120  ldirect = 0
121  ibdate = 20000101
122  ibtime = 0
123  iedate = 20000102
124  ietime = 0
125  loutstep = 10800
126  loutaver = 10800
127  loutsample = 900
128  itsplit = 999999999
129  lsynctime = 900
130  ctl = -5.0
131  ifine = 4
132  iout = 3
133  ipout = 0
134  lsubgrid = 1
135  lconvection = 1
136  lagespectra = 0
137  ipin = 1
138  ioutputforeachrelease = 1
139  iflux = 1
140  mdomainfill = 0
141  ind_source = 1
142  ind_receptor = 1
143  mquasilag = 0
144  nested_output = 0
145  linit_cond = 0
146  lnetcdfout = 0
147  surf_only = 0
148  ldep_incr = .false.
149
150  ! Open the command file and read user options
151  ! Namelist input first: try to read as namelist file
152  !**************************************************************************
153  open(unitcommand,file=path(1)(1:length(1))//'COMMAND',status='old',form='formatted',err=999)
154
155  ! try namelist input (default)
156  read(unitcommand,command,iostat=ireaderror)
157  close(unitcommand)
158
159  ! distinguish namelist from fixed text input
160  if ((ireaderror .ne. 0) .or. (ldirect .eq. 0)) then
161   
162    ! parse as text file format
163 
164    open(unitcommand,file=path(1)(1:length(1))//'COMMAND',status='old', err=999)
165
166    ! Check the format of the COMMAND file (either in free format,
167    ! or using formatted mask)
168    ! Use of formatted mask is assumed if line 10 contains the word 'DIRECTION'
169    !**************************************************************************
170
171    call skplin(9,unitcommand)
172    read (unitcommand,901) line
173  901   format (a)
174    if (index(line,'LDIRECT') .eq. 0) then
175      old = .false.
176      write(*,*) 'COMMAND in old short format, please update to namelist format'
177    else
178      old = .true.
179      write(*,*) 'COMMAND in old long format, please update to namelist format'
180    endif
181    rewind(unitcommand)
182
183
184    ! Read parameters
185    !****************
186
187    call skplin(7,unitcommand)
188    if (old) call skplin(1,unitcommand)
189    read(unitcommand,*) ldirect
190    if (old) call skplin(3,unitcommand)
191    read(unitcommand,*) ibdate,ibtime
192    if (old) call skplin(3,unitcommand)
193    read(unitcommand,*) iedate,ietime
194    if (old) call skplin(3,unitcommand)
195    read(unitcommand,*) loutstep
196    if (old) call skplin(3,unitcommand)
197    read(unitcommand,*) loutaver
198    if (old) call skplin(3,unitcommand)
199    read(unitcommand,*) loutsample
200    if (old) call skplin(3,unitcommand)
201    read(unitcommand,*) itsplit
202    if (old) call skplin(3,unitcommand)
203    read(unitcommand,*) lsynctime
204    if (old) call skplin(3,unitcommand)
205    read(unitcommand,*) ctl
206    if (old) call skplin(3,unitcommand)
207    read(unitcommand,*) ifine
208    if (old) call skplin(3,unitcommand)
209    read(unitcommand,*) iout
210    if (old) call skplin(3,unitcommand)
211    read(unitcommand,*) ipout
212    if (old) call skplin(3,unitcommand)
213    read(unitcommand,*) lsubgrid
214    if (old) call skplin(3,unitcommand)
215    read(unitcommand,*) lconvection
216    if (old) call skplin(3,unitcommand)
217    read(unitcommand,*) lagespectra
218    if (old) call skplin(3,unitcommand)
219    read(unitcommand,*) ipin
220    if (old) call skplin(3,unitcommand)
221    read(unitcommand,*) ioutputforeachrelease
222    if (old) call skplin(3,unitcommand)
223    read(unitcommand,*) iflux
224    if (old) call skplin(3,unitcommand)
225    read(unitcommand,*) mdomainfill
226    if (old) call skplin(3,unitcommand)
227    read(unitcommand,*) ind_source
228    if (old) call skplin(3,unitcommand)
229    read(unitcommand,*) ind_receptor
230    if (old) call skplin(3,unitcommand)
231    read(unitcommand,*) mquasilag
232    if (old) call skplin(3,unitcommand)
233    read(unitcommand,*) nested_output
234    if (old) call skplin(3,unitcommand)
235    read(unitcommand,*,iostat=icmdstat) linit_cond
236    if (icmdstat .gt. 0) &
237      print*, 'readcommand: linit_cond not present or no read properly,&
238       & - check values: ',icmdstat,linit_cond
239    if (old) call skplin(3,unitcommand)
240    read(unitcommand,*,iostat=icmdstat) icmdstat,surf_only
241    if (icmdstat .gt. 0) &
242      print*, 'readcommand: surf_only not present or not read properly&
243       & - check values: ',icmdstat,surf_only
244    if (old) call skplin(3,unitcommand)
245    read(unitcommand,*,iostat=icmdstat) icmdstat,ldep_incr
246    if (icmdstat .gt. 0) &
247      print*, 'readcommand: ldep_incr not present or not read properly&
248       & - check values: ',icmdstat, ldep_incr
249    close(unitcommand)
250
251  endif ! input format
252
253  ! write command file in namelist format to output directory if requested
254  if (nmlout .eqv. .true.) then
255    open(unitcommand,file=path(2)(1:length(2))//'COMMAND.namelist',err=1000)
256    write(unitcommand,nml=command)
257    close(unitcommand)
258  endif
259
260  ifine=max(ifine,1)
261
262  ! Determine how Markov chain is formulated (for w or for w/sigw)
263  !***************************************************************
264
265  if (ctl .ge. 0.1) then
266    turbswitch=.true.
267  else
268    turbswitch=.false.
269    ifine=1
270  endif
271  fine=1./real(ifine)
272  ctl=1./ctl
273
274  ! Set the switches required for the various options for input/output units
275  !*************************************************************************
276  !AF Set the switches IND_REL and IND_SAMP for the release and sampling
277  !Af switches for the releasefile:
278  !Af IND_REL =  1 : xmass * rho
279  !Af IND_REL =  0 : xmass * 1
280
281  !Af switches for the conccalcfile:
282  !AF IND_SAMP =  0 : xmass * 1
283  !Af IND_SAMP = -1 : xmass / rho
284
285  !AF IND_SOURCE switches between different units for concentrations at the source
286  !Af   NOTE that in backward simulations the release of computational particles
287  !Af   takes place at the "receptor" and the sampling of p[articles at the "source".
288  !Af          1 = mass units
289  !Af          2 = mass mixing ratio units
290  !Af IND_RECEPTOR switches between different units for concentrations at the receptor
291  !Af          1 = mass units
292  !Af          2 = mass mixing ratio units
293
294  if ( ldirect .eq. 1 ) then  ! FWD-Run
295  !Af set release-switch
296     if (ind_source .eq. 1 ) then !mass
297        ind_rel = 0
298     else ! mass mix
299        ind_rel = 1
300     endif
301  !Af set sampling switch
302     if (ind_receptor .eq. 1) then !mass
303        ind_samp = 0
304     else ! mass mix
305        ind_samp = -1
306     endif
307  elseif (ldirect .eq. -1 ) then !BWD-Run
308  !Af set sampling switch
309     if (ind_source .eq. 1 ) then !mass
310        ind_samp = -1
311     else ! mass mix
312        ind_samp = 0
313     endif
314  !Af set release-switch
315     if (ind_receptor .eq. 1) then !mass
316        ind_rel = 1
317     else ! mass mix
318        ind_rel = 0
319     endif
320  endif
321
322  !*************************************************************
323  ! Check whether valid options have been chosen in file COMMAND
324  !*************************************************************
325
326  ! Check options for initial condition output: Switch off for forward runs
327  !************************************************************************
328
329  if (ldirect .eq. 1) linit_cond=0
330  if ((linit_cond .lt. 0) .or. (linit_cond .gt. 2)) then
331    write(*,*) ' #### FLEXPART MODEL ERROR! INVALID OPTION    #### '
332    write(*,*) ' #### FOR LINIT_COND IN FILE "COMMAND".       #### '
333    stop
334  endif
335
336  ! Check input dates
337  !******************
338
339  if (iedate .lt. ibdate) then
340    write(*,*) ' #### FLEXPART MODEL ERROR! BEGINNING DATE    #### '
341    write(*,*) ' #### IS LARGER THAN ENDING DATE. CHANGE      #### '
342    write(*,*) ' #### EITHER POINT 2 OR POINT 3 IN FILE       #### '
343    write(*,*) ' #### "COMMAND".                              #### '
344    stop
345  else if (iedate .eq. ibdate) then
346    if (ietime .lt. ibtime) then
347    write(*,*) ' #### FLEXPART MODEL ERROR! BEGINNING TIME    #### '
348    write(*,*) ' #### IS LARGER THAN ENDING TIME. CHANGE      #### '
349    write(*,*) ' #### EITHER POINT 2 OR POINT 3 IN FILE       #### '
350    write(*,*) ' #### "COMMAND".                              #### '
351    stop
352    endif
353  endif
354
355
356  ! Determine kind of dispersion method
357  !************************************
358
359  if (ctl .gt. 0.) then
360    method=1
361    mintime=minstep
362  else
363    method=0
364    mintime=lsynctime
365  endif
366
367!  check for netcdf output switch (use for non-namelist input only!)
368  if (iout .ge. 8) then
369     lnetcdfout = 1
370     iout = iout - 8
371#ifndef NETCDF_OUTPUT
372     print*,'ERROR: netcdf output not activated during compile time but used in COMMAND file!'
373     print*,'Please recompile with netcdf library or use standard output format.'
374     stop
375#endif
376  endif
377
378  ! Check whether a valid option for gridded model output has been chosen
379  !**********************************************************************
380
381  if ((iout .lt. 1) .or.  (iout .gt. 5)) then
382    write(*,*) ' #### FLEXPART MODEL ERROR! FILE COMMAND:     #### '
383    write(*,*) ' #### IOUT MUST BE 1, 2, 3, 4, OR 5 FOR       #### '
384    write(*,*) ' #### STANDARD FLEXPART OUTPUT OR  9 - 13     #### '
385    write(*,*) ' #### FOR NETCDF OUTPUT                       #### '
386    stop
387  endif
388
389  !AF check consistency between units and volume mixing ratio
390  if ( ((iout .eq. 2) .or.  (iout .eq. 3)) .and.  &
391       (ind_source .gt. 1 .or. ind_receptor .gt. 1) ) then
392    write(*,*) ' #### FLEXPART MODEL ERROR! FILE COMMAND:     #### '
393    write(*,*) ' #### VOLUME MIXING RATIO ONLY SUPPORTED      #### '
394    write(*,*) ' #### FOR MASS UNITS (at the moment)          #### '
395    stop
396  endif
397
398
399  ! For quasilag output for each release is forbidden
400  !*****************************************************************************
401
402  if ((ioutputforeachrelease .eq. 1) .and. (mquasilag .eq. 1)) then
403      write(*,*) '#### FLEXPART MODEL ERROR! FILE COMMAND:     ####'
404      write(*,*) '#### OUTPUTFOREACHRELEASE AND QUASILAGRANGIAN####'
405      write(*,*) '#### MODE IS NOT POSSIBLE   !                ####'
406      stop
407  endif
408
409
410  ! For quasilag backward is forbidden
411  !*****************************************************************************
412
413  if ((ldirect .lt. 0) .and. (mquasilag .eq. 1)) then
414      write(*,*) '#### FLEXPART MODEL ERROR! FILE COMMAND:     ####'
415      write(*,*) '#### FOR BACKWARD RUNS, QUASILAGRANGIAN MODE ####'
416      write(*,*) '#### IS NOT POSSIBLE   !                     ####'
417      stop
418  endif
419
420
421  ! For domainfill runs one releasefield for all releases makes no sense,
422  ! and is "forbidden"
423  !*****************************************************************************
424
425  if ((mdomainfill .eq. 1) .and. (ioutputforeachrelease .eq. 1)) then
426      write(*,*) '#### FLEXPART MODEL ERROR! FILE COMMAND:     ####'
427      write(*,*) '#### FOR DOMAIN FILLING RUNS OUTPUT FOR      ####'
428      write(*,*) '#### EACH RELEASE IS FORBIDDEN !             ####'
429      stop
430  endif
431
432
433  ! For domain-filling trajectories, a plume centroid trajectory makes no sense,
434  ! For backward runs, only residence time output (iout=1) or plume trajectories (iout=4),
435  ! or both (iout=5) makes sense; other output options are "forbidden"
436  !*****************************************************************************
437
438  if (ldirect .lt. 0) then
439    if ((iout .eq. 2) .or. (iout .eq. 3)) then
440      write(*,*) '#### FLEXPART MODEL ERROR! FILE COMMAND:     ####'
441      write(*,*) '#### FOR BACKWARD RUNS, IOUT MUST BE 1,4,OR 5####'
442      stop
443    endif
444  endif
445
446
447  ! For domain-filling trajectories, a plume centroid trajectory makes no sense,
448  ! and is "forbidden"
449  !*****************************************************************************
450
451  if (mdomainfill .ge. 1) then
452    if ((iout .eq. 4) .or. (iout .eq. 5)) then
453      write(*,*) '#### FLEXPART MODEL ERROR! FILE COMMAND:     ####'
454      write(*,*) '#### FOR DOMAIN-FILLING TRAJECTORY OPTION,   ####'
455      write(*,*) '#### IOUT MUST NOT BE SET TO 4 OR 5.         ####'
456      stop
457    endif
458  endif
459
460
461
462  ! Check whether a valid options for particle dump has been chosen
463  !****************************************************************
464
465  if ((ipout .ne. 0) .and. (ipout .ne. 1) .and. (ipout .ne. 2)) then
466    write(*,*) ' #### FLEXPART MODEL ERROR! FILE COMMAND:     #### '
467    write(*,*) ' #### IPOUT MUST BE 1, 2 OR 3!                #### '
468    stop
469  endif
470
471  if(lsubgrid .ne. 1 .and. verbosity .eq. 0) then
472    write(*,*) '             ----------------               '
473    write(*,*) ' INFORMATION: SUBGRIDSCALE TERRAIN EFFECT IS'
474    write(*,*) ' NOT PARAMETERIZED DURING THIS SIMULATION.  '
475    write(*,*) '             ----------------               '
476  endif
477
478
479  ! Check whether convection scheme is either turned on or off
480  !***********************************************************
481
482  if ((lconvection .ne. 0) .and. (lconvection .ne. 1)) then
483    write(*,*) ' #### FLEXPART MODEL ERROR! FILE COMMAND:     #### '
484    write(*,*) ' #### LCONVECTION MUST BE SET TO EITHER 1 OR 0#### '
485    stop
486  endif
487
488
489  ! Check whether synchronisation interval is sufficiently short
490  !*************************************************************
491
492  if (lsynctime .gt. (idiffnorm/2)) then
493    write(*,*) ' #### FLEXPART MODEL ERROR! SYNCHRONISATION   #### '
494    write(*,*) ' #### TIME IS TOO LONG. MAKE IT SHORTER.      #### '
495    write(*,*) ' #### MINIMUM HAS TO BE: ', idiffnorm/2
496    stop
497  endif
498
499
500  ! Check consistency of the intervals, sampling periods, etc., for model output
501  !*****************************************************************************
502
503  if (loutaver .eq. 0) then
504    write(*,*) ' #### FLEXPART MODEL ERROR! TIME AVERAGE OF   #### '
505    write(*,*) ' #### CONCENTRATION FIELD OUTPUT MUST NOT BE  #### '
506    write(*,*) ' #### ZERO.                                   #### '
507    write(*,*) ' #### CHANGE INPUT IN FILE COMMAND.           #### '
508    stop
509  endif
510
511  if (loutaver .gt. loutstep) then
512    write(*,*) ' #### FLEXPART MODEL ERROR! TIME AVERAGE OF   #### '
513    write(*,*) ' #### CONCENTRATION FIELD OUTPUT MUST NOT BE  #### '
514    write(*,*) ' #### GREATER THAN INTERVAL OF OUTPUT.        #### '
515    write(*,*) ' #### CHANGE INPUT IN FILE COMMAND.           #### '
516    stop
517  endif
518
519  if (loutsample .gt. loutaver) then
520    write(*,*) ' #### FLEXPART MODEL ERROR! SAMPLING TIME OF  #### '
521    write(*,*) ' #### CONCENTRATION FIELD OUTPUT MUST NOT BE  #### '
522    write(*,*) ' #### GREATER THAN TIME AVERAGE OF OUTPUT.    #### '
523    write(*,*) ' #### CHANGE INPUT IN FILE COMMAND.           #### '
524    stop
525  endif
526
527  if (mod(loutaver,lsynctime) .ne. 0) then
528    write(*,*) ' #### FLEXPART MODEL ERROR! AVERAGING TIME OF #### '
529    write(*,*) ' #### CONCENTRATION FIELD MUST BE A MULTIPLE  #### '
530    write(*,*) ' #### OF THE SYNCHRONISATION INTERVAL         #### '
531    stop
532  endif
533
534  if ((loutaver/lsynctime) .lt. 2) then
535    write(*,*) ' #### FLEXPART MODEL ERROR! AVERAGING TIME OF #### '
536    write(*,*) ' #### CONCENTRATION FIELD MUST BE AT LEAST    #### '
537    write(*,*) ' #### TWICE THE SYNCHRONISATION INTERVAL      #### '
538    stop
539  endif
540
541  if (mod(loutstep,lsynctime) .ne. 0) then
542    write(*,*) ' #### FLEXPART MODEL ERROR! INTERVAL BETWEEN  #### '
543    write(*,*) ' #### CONCENTRATION FIELDS MUST BE A MULTIPLE #### '
544    write(*,*) ' #### OF THE SYNCHRONISATION INTERVAL         #### '
545    stop
546  endif
547
548  if ((loutstep/lsynctime) .lt. 2) then
549    write(*,*) ' #### FLEXPART MODEL ERROR! INTERVAL BETWEEN  #### '
550    write(*,*) ' #### CONCENTRATION FIELDS MUST BE AT LEAST   #### '
551    write(*,*) ' #### TWICE THE SYNCHRONISATION INTERVAL      #### '
552    stop
553  endif
554
555  if (mod(loutsample,lsynctime) .ne. 0) then
556    write(*,*) ' #### FLEXPART MODEL ERROR! SAMPLING TIME OF  #### '
557    write(*,*) ' #### CONCENTRATION FIELD MUST BE A MULTIPLE  #### '
558    write(*,*) ' #### OF THE SYNCHRONISATION INTERVAL         #### '
559    stop
560  endif
561
562  if (itsplit .lt. loutaver) then
563    write(*,*) ' #### FLEXPART MODEL ERROR! SPLITTING TIME FOR#### '
564    write(*,*) ' #### PARTICLES IS TOO SHORT. PLEASE INCREASE #### '
565    write(*,*) ' #### SPLITTING TIME CONSTANT.                #### '
566    stop
567  endif
568
569  if ((mquasilag .eq. 1) .and. (iout .ge. 4)) then
570    write(*,*) ' #### FLEXPART MODEL ERROR! CONFLICTING       #### '
571    write(*,*) ' #### OPTIONS: IF MQUASILAG=1, PLUME          #### '
572    write(*,*) ' #### TRAJECTORY OUTPUT IS IMPOSSIBLE.        #### '
573    stop
574  endif
575
576  ! Compute modeling time in seconds and beginning date in Julian date
577  !*******************************************************************
578
579  outstep=real(abs(loutstep))
580  if (ldirect .eq. 1) then
581    bdate=juldate(ibdate,ibtime)
582    edate=juldate(iedate,ietime)
583    ideltas=nint((edate-bdate)*86400.)
584  else if (ldirect .eq. -1) then
585    loutaver=-1*loutaver
586    loutstep=-1*loutstep
587    loutsample=-1*loutsample
588    lsynctime=-1*lsynctime
589    bdate=juldate(iedate,ietime)
590    edate=juldate(ibdate,ibtime)
591    ideltas=nint((edate-bdate)*86400.)
592  else
593    write(*,*) ' #### FLEXPART MODEL ERROR! DIRECTION IN      #### '
594    write(*,*) ' #### FILE "COMMAND" MUST BE EITHER -1 OR 1.  #### '
595    stop
596  endif
597
598  return
599
600999   write(*,*) ' #### FLEXPART MODEL ERROR! FILE "COMMAND"    #### '
601  write(*,*) ' #### CANNOT BE OPENED IN THE DIRECTORY       #### '
602  write(*,'(a)') path(1)(1:length(1))
603  stop
604
6051000   write(*,*) ' #### FLEXPART MODEL ERROR! FILE "COMMAND"    #### '
606  write(*,*) ' #### CANNOT BE OPENED IN THE DIRECTORY       #### '
607  write(*,'(a)') path(2)(1:length(1))
608  stop
609end subroutine readcommand
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG