Changeset ff050cd in flexpart.git for src_parallel/readreleases.f90
- Timestamp:
- Aug 15, 2013, 3:23:48 PM (11 years ago)
- Branches:
- flexpart91_hasod
- Children:
- 31113de
- Parents:
- 7c1fd44
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
src_parallel/readreleases.f90
r3eed2e6 rff050cd 33 33 ! Update: 29 January 2001 * 34 34 ! Release altitude can be either in magl or masl * 35 ! HSO, 12 August 2013 36 ! Added optional namelist input 35 37 ! * 36 38 !***************************************************************************** … … 67 69 implicit none 68 70 69 integer :: numpartmax,i,j,id1,it1,id2,it2,specnum_rel,idum,stat 71 integer :: numpartmax,i,j,id1,it1,id2,it2,idum,stat 72 integer :: specnum_rel(maxspec) 70 73 real :: vsh(ni),fracth(ni),schmih(ni),releaserate,xdum,cun 71 74 real(kind=dp) :: jul1,jul2,juldate … … 73 76 logical :: old 74 77 78 ! help variables for namelist reading 79 integer :: numpoints, parts, readerror 80 integer*2 :: zkind 81 integer :: idate1, itime1, idate2, itime2 82 real :: lon1,lon2,lat1,lat2,z1,z2,mass(nspec) 83 character*40 :: comment 84 85 ! declare namelists 86 namelist /releases_ctrl/ & 87 nspec, & 88 specnum_rel 89 90 namelist /release/ & 91 idate1, itime1, & 92 idate2, itime2, & 93 lon1, lon2, & 94 lat1, lat2, & 95 z1, z2, & 96 zkind, & 97 mass, & 98 parts, & 99 comment 100 101 numpoint=0 102 103 ! presetting namelist releases_ctrl 104 nspec = -1 ! use negative value to determine failed namelist input 105 specnum_rel(1) = 0 106 75 107 !sec, read release to find how many releasepoints should be allocated 76 77 108 open(unitreleases,file=path(1)(1:length(1))//'RELEASES',status='old', & 78 109 err=999) 110 111 ! check if namelist input provided 112 read(unitreleases,releases_ctrl,iostat=readerror) 113 114 if ((readerror.ne.0).or.(nspec.lt.0)) then 115 116 ! no namelist format, close file and allow reopening in old format 117 rewind(unitreleases) 118 119 readerror=1 ! indicates old format 79 120 80 121 ! Check the format of the RELEASES file (either in free format, … … 93 134 rewind(unitreleases) 94 135 95 96 136 ! Skip first 11 lines (file header) 97 137 !********************************** 98 138 99 139 call skplin(11,unitreleases) 100 101 102 140 read(unitreleases,*,err=998) nspec 103 141 if (old) call skplin(2,unitreleases) 104 142 do i=1,nspec 105 read(unitreleases,*,err=998) specnum_rel 143 read(unitreleases,*,err=998) specnum_rel(1) 106 144 if (old) call skplin(2,unitreleases) 107 145 end do 108 146 109 numpoint=0110 147 100 numpoint=numpoint+1 111 148 read(unitreleases,*,end=25) … … 142 179 25 numpoint=numpoint-1 143 180 181 else 182 readerror=0 183 do while (readerror.eq.0) 184 idate1=-1 185 read(unitreleases,release,iostat=readerror) 186 if ((idate1.lt.0).or.(readerror.ne.0)) then 187 readerror=1 188 else 189 numpoint=numpoint+1 190 endif 191 end do 192 readerror=0 193 endif ! if namelist input 194 195 rewind(unitreleases) 196 144 197 !allocate memory for numpoint releaspoint 145 allocate(ireleasestart(numpoint) & 146 ,stat=stat) 147 if (stat.ne.0) write(*,*)'ERROR: could not allocate RELEASPOINT' 148 allocate(ireleaseend(numpoint) & 149 ,stat=stat) 150 if (stat.ne.0) write(*,*)'ERROR: could not allocate RELEASPOINT' 151 allocate(xpoint1(numpoint) & 152 ,stat=stat) 153 if (stat.ne.0) write(*,*)'ERROR: could not allocate RELEASPOINT' 154 allocate(xpoint2(numpoint) & 155 ,stat=stat) 156 if (stat.ne.0) write(*,*)'ERROR: could not allocate RELEASPOINT' 157 allocate(ypoint1(numpoint) & 158 ,stat=stat) 159 if (stat.ne.0) write(*,*)'ERROR: could not allocate RELEASPOINT' 160 allocate(ypoint2(numpoint) & 161 ,stat=stat) 162 if (stat.ne.0) write(*,*)'ERROR: could not allocate RELEASPOINT' 163 allocate(zpoint1(numpoint) & 164 ,stat=stat) 165 if (stat.ne.0) write(*,*)'ERROR: could not allocate RELEASPOINT' 166 allocate(zpoint2(numpoint) & 167 ,stat=stat) 168 if (stat.ne.0) write(*,*)'ERROR: could not allocate RELEASPOINT' 169 allocate(kindz(numpoint) & 170 ,stat=stat) 171 if (stat.ne.0) write(*,*)'ERROR: could not allocate RELEASPOINT' 172 allocate(xmass(numpoint,maxspec) & 173 ,stat=stat) 174 if (stat.ne.0) write(*,*)'ERROR: could not allocate RELEASPOINT' 175 allocate(rho_rel(numpoint) & 176 ,stat=stat) 177 if (stat.ne.0) write(*,*)'ERROR: could not allocate RELEASPOINT' 178 allocate(npart(numpoint) & 179 ,stat=stat) 180 if (stat.ne.0) write(*,*)'ERROR: could not allocate RELEASPOINT' 181 allocate(xmasssave(numpoint) & 182 ,stat=stat) 183 if (stat.ne.0) write(*,*)'ERROR: could not allocate RELEASPOINT' 184 185 write (*,*) ' Releasepoints allocated: ', numpoint 186 187 do i=1,numpoint 188 xmasssave(i)=0. 189 end do 198 allocate(ireleasestart(numpoint),stat=stat) 199 200 if (stat.ne.0) write(*,*)'ERROR: could not allocate releasepoint' 201 allocate(ireleaseend(numpoint),stat=stat) 202 if (stat.ne.0) write(*,*)'ERROR: could not allocate releasepoint' 203 allocate(xpoint1(numpoint),stat=stat) 204 if (stat.ne.0) write(*,*)'ERROR: could not allocate releasepoint' 205 allocate(xpoint2(numpoint),stat=stat) 206 if (stat.ne.0) write(*,*)'ERROR: could not allocate releasepoint' 207 allocate(ypoint1(numpoint),stat=stat) 208 if (stat.ne.0) write(*,*)'ERROR: could not allocate releasepoint' 209 allocate(ypoint2(numpoint),stat=stat) 210 if (stat.ne.0) write(*,*)'ERROR: could not allocate releasepoint' 211 allocate(zpoint1(numpoint),stat=stat) 212 if (stat.ne.0) write(*,*)'ERROR: could not allocate releasepoint' 213 allocate(zpoint2(numpoint),stat=stat) 214 if (stat.ne.0) write(*,*)'ERROR: could not allocate releasepoint' 215 allocate(kindz(numpoint),stat=stat) 216 if (stat.ne.0) write(*,*)'ERROR: could not allocate releasepoint' 217 allocate(xmass(numpoint,maxspec),stat=stat) 218 if (stat.ne.0) write(*,*)'ERROR: could not allocate releasepoint' 219 allocate(rho_rel(numpoint),stat=stat) 220 if (stat.ne.0) write(*,*)'ERROR: could not allocate releasepoint' 221 allocate(npart(numpoint),stat=stat) 222 if (stat.ne.0) write(*,*)'ERROR: could not allocate releasepoint' 223 allocate(xmasssave(numpoint),stat=stat) 224 if (stat.ne.0) write(*,*)'ERROR: could not allocate releasepoint' 225 226 write (*,*) numpoint,' releasepoints allocated.' 227 228 do i=1,numpoint 229 xmasssave(i)=0. 230 end do 190 231 191 232 !now save the information … … 198 239 end do 199 240 200 rewind(unitreleases) 201 202 241 if (readerror.ne.0) then 203 242 ! Skip first 11 lines (file header) 204 243 !********************************** 205 244 206 call skplin(11,unitreleases) 207 245 call skplin(11,unitreleases) 208 246 209 247 ! Assign species-specific parameters needed for physical processes 210 248 !************************************************************************* 211 249 212 read(unitreleases,*,err=998) nspec 213 if (nspec.gt.maxspec) goto 994 214 if (old) call skplin(2,unitreleases) 250 read(unitreleases,*,err=998) nspec 251 if (nspec.gt.maxspec) goto 994 252 if (old) call skplin(2,unitreleases) 253 endif 215 254 do i=1,nspec 216 read(unitreleases,*,err=998) specnum_rel 217 if (old) call skplin(2,unitreleases) 218 call readspecies(specnum_rel,i) 255 if (readerror.ne.0) then 256 read(unitreleases,*,err=998) specnum_rel(1) 257 if (old) call skplin(2,unitreleases) 258 call readspecies(specnum_rel(1),i) 259 else 260 call readspecies(specnum_rel(i),i) 261 endif 219 262 220 263 ! For backward runs, only 1 species is allowed … … 301 344 ! Read specifications for each release point 302 345 !******************************************* 303 346 numpoints=numpoint 304 347 numpoint=0 305 348 numpartmax=0 306 349 releaserate=0. 307 350 1000 numpoint=numpoint+1 351 352 if (readerror.eq.0) then ! reading namelist format 353 354 if (numpoint.gt.numpoints) goto 250 355 zkind = 1 356 mass = 0 357 parts = 0 358 comment = ' ' 359 read(unitreleases,release,iostat=readerror) 360 id1=idate1 361 it1=itime1 362 id2=idate2 363 it2=itime2 364 xpoint1(numpoint)=lon1 365 xpoint2(numpoint)=lon2 366 ypoint1(numpoint)=lat1 367 ypoint2(numpoint)=lat2 368 zpoint1(numpoint)=z1 369 zpoint2(numpoint)=z2 370 kindz(numpoint)=zkind 371 do i=1,nspec 372 xmass(numpoint,i)=mass(i) 373 end do 374 npart(numpoint)=parts 375 compoint(min(1001,numpoint))=comment 376 377 else 378 308 379 read(unitreleases,*,end=250) 309 380 read(unitreleases,*,err=998,end=250) id1,it1 … … 347 418 endif 348 419 420 endif ! if namelist format 349 421 350 422 ! If a release point contains no particles, stop and issue error message
Note: See TracChangeset
for help on using the changeset viewer.