source: flexpart.git/src/readcommand.f90 @ 462f74b

10.4.1_peseiGFS_025bugfixes+enhancementsdevrelease-10release-10.4.1scaling-bugunivie
Last change on this file since 462f74b was 462f74b, checked in by Sabine Eckhardt <sabine@…>, 8 years ago

first version of the backward scavenging

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