source: flexpart.git/src/readpartpositions.f90 @ e200b7a

10.4.1_peseiFPv9.3.1FPv9.3.1b_testingFPv9.3.2GFS_025NetCDFbugfixes+enhancementsdepositiondevflexpart-noresmfp9.3.1-20161214-nc4grib2nc4_repairinputlistlaptoprelease-10release-10.4.1scaling-bugsvn-petrasvn-trunkunivie
Last change on this file since e200b7a was e200b7a, checked in by Matthias Langer <matthias.langer@…>, 11 years ago

git-svn-id: http://flexpart.flexpart.eu:8088/svn/FlexPart90/trunk@4 ef8cc7e1-21b7-489e-abab-c1baa636049d

  • Property mode set to 100644
File size: 6.3 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 readpartpositions
23
24  !*****************************************************************************
25  !                                                                            *
26  !   This routine opens the particle dump file and reads all the particle     *
27  !   positions from a previous run to initialize the current run.             *
28  !                                                                            *
29  !                                                                            *
30  !     Author: A. Stohl                                                       *
31  !                                                                            *
32  !     24 March 2000                                                          *
33  !                                                                            *
34  !*****************************************************************************
35  !                                                                            *
36  ! Variables:                                                                 *
37  !                                                                            *
38  !*****************************************************************************
39
40  use par_mod
41  use com_mod
42
43  implicit none
44
45  integer :: ibdatein,ibtimein,nspecin,itimein,numpointin,i,j,ix
46  integer :: id1,id2,it1,it2
47  real :: xlonin,ylatin,ran1,topo,hmixi,pvi,qvi,rhoi,tri,tti
48  character :: specin*7
49  real(kind=dp) :: julin,julpartin,juldate
50
51  integer :: idummy = -8
52
53  numparticlecount=0
54
55  ! Open header file of dumped particle data
56  !*****************************************
57
58  open(unitpartin,file=path(2)(1:length(2))//'header', &
59       form='unformatted',err=998)
60
61  read(unitpartin) ibdatein,ibtimein
62  read(unitpartin)
63  read(unitpartin)
64
65  read(unitpartin)
66  read(unitpartin)
67  read(unitpartin) nspecin
68  nspecin=nspecin/3
69  if ((ldirect.eq.1).and.(nspec.ne.nspecin)) goto 997
70
71  do i=1,nspecin
72    read(unitpartin)
73    read(unitpartin)
74    read(unitpartin) j,specin
75    if ((ldirect.eq.1).and.(species(i)(1:7).ne.specin)) goto 996
76  end do
77
78  read(unitpartin) numpointin
79  if (numpointin.ne.numpoint) then
80     write(*,*) ' ####              WARNING                     #### '
81     write(*,*) ' #### NUMBER OF RELEASE LOCATIONS DOES NOT     #### '
82     write(*,*) ' #### AGREE WITH CURRENT SETTINGS!             #### '
83     write(*,*) ' #### THIS IS JUST A CHECK TO SEE IF YOU READ  #### '
84     write(*,*) ' #### THE CORRECT PARTICLE DUMP FILE.          #### '
85  endif
86  do i=1,numpointin
87    read(unitpartin)
88    read(unitpartin)
89    read(unitpartin)
90    read(unitpartin)
91    do j=1,nspec
92      read(unitpartin)
93      read(unitpartin)
94      read(unitpartin)
95    end do
96  end do
97  read(unitpartin)
98  read(unitpartin)
99
100  do ix=0,numxgrid-1
101    read(unitpartin)
102  end do
103
104
105  ! Open data file of dumped particle data
106  !***************************************
107
108  close(unitpartin)
109  open(unitpartin,file=path(2)(1:length(2))//'partposit_end', &
110       form='unformatted',err=998)
111
112
113100   read(unitpartin,end=99) itimein
114    i=0
115200   i=i+1
116    read(unitpartin) npoint(i),xlonin,ylatin,ztra1(i),itramem(i), &
117         topo,pvi,qvi,rhoi,hmixi,tri,tti,(xmass1(i,j),j=1,nspec)
118
119    if (xlonin.eq.-9999.9) goto 100
120    xtra1(i)=(xlonin-xlon0)/dx
121    ytra1(i)=(ylatin-ylat0)/dy
122    numparticlecount=max(numparticlecount,npoint(i))
123    goto 200
124
12599   numpart=i-1
126
127  close(unitpartin)
128
129  julin=juldate(ibdatein,ibtimein)+real(itimein,kind=dp)/86400._dp
130  if (abs(julin-bdate).gt.1.e-5) goto 994
131  do i=1,numpart
132    julpartin=juldate(ibdatein,ibtimein)+ &
133         real(itramem(i),kind=dp)/86400._dp
134    nclass(i)=min(int(ran1(idummy)*real(nclassunc))+1, &
135         nclassunc)
136    idt(i)=mintime
137    itra1(i)=0
138    itramem(i)=nint((julpartin-bdate)*86400.)
139    itrasplit(i)=ldirect*itsplit
140  end do
141
142  return
143
144
145994   write(*,*) ' #### FLEXPART MODEL ERROR IN READPARTPOSITIONS#### '
146  write(*,*) ' #### ENDING TIME OF PREVIOUS MODEL RUN DOES   #### '
147  write(*,*) ' #### NOT AGREE WITH STARTING TIME OF THIS RUN.#### '
148  call caldate(julin,id1,it1)
149  call caldate(bdate,id2,it2)
150  write(*,*) 'julin: ',julin,id1,it1
151  write(*,*) 'bdate: ',bdate,id2,it2
152  stop
153
154
155996   write(*,*) ' #### FLEXPART MODEL ERROR IN READPARTPOSITIONS#### '
156  write(*,*) ' #### SPECIES NAMES TO BE READ IN DO NOT       #### '
157  write(*,*) ' #### AGREE WITH CURRENT SETTINGS!             #### '
158  stop
159
160997   write(*,*) ' #### FLEXPART MODEL ERROR IN READPARTPOSITIONS#### '
161  write(*,*) ' #### THE NUMBER OF SPECIES TO BE READ IN DOES #### '
162  write(*,*) ' #### NOT AGREE WITH CURRENT SETTINGS!         #### '
163  stop
164
165998   write(*,*) ' #### FLEXPART MODEL ERROR!   THE FILE         #### '
166  write(*,*) ' #### '//path(2)(1:length(2))//'grid'//' #### '
167  write(*,*) ' #### CANNOT BE OPENED. IF A FILE WITH THIS    #### '
168  write(*,*) ' #### NAME ALREADY EXISTS, DELETE IT AND START #### '
169  write(*,*) ' #### THE PROGRAM AGAIN.                       #### '
170  stop
171
172end subroutine readpartpositions
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG