Changeset 3c1da52 in flexpart.git for src/readspecies.f90
- Timestamp:
- Jul 15, 2015, 9:45:31 AM (9 years ago)
- Branches:
- svn-petra
- Children:
- c07b09e
- Parents:
- 31674b5
- git-author:
- Ignacio Pisso <Ignacio.Pisso@…> (07/15/15 09:23:56)
- git-committer:
- Ignacio Pisso <Ignacio.Pisso@…> (07/15/15 09:45:31)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
src/readspecies.f90
r31674b5 r3c1da52 20 20 !********************************************************************** 21 21 22 subroutine readspecies(id_spec, pos_spec)22 subroutine readspecies(id_spec,id_pos) 23 23 24 24 !***************************************************************************** 25 25 ! * 26 26 ! This routine reads names and physical constants of chemical species/ * 27 ! radionuclides given in the parameter pos_spec*27 ! radionuclides given in the parameter id_pos * 28 28 ! * 29 29 ! Author: A. Stohl * … … 32 32 ! 33 33 ! Changes: * 34 ! * 34 35 ! N. Kristiansen, 31.01.2013: Including parameters for in-cloud scavenging * 35 36 ! * 36 37 ! HSO, 13 August 2013: added optional namelist input 37 38 ! PS, 2/2015: access= -> position= 39 ! PS, 4/2015: quick fix for bug associated with end=22 38 40 ! * 39 41 !***************************************************************************** … … 60 62 implicit none 61 63 62 integer :: i, pos_spec,j64 integer :: i, id_pos,j 63 65 integer :: idow,ihour,id_spec 64 66 character(len=3) :: aspecnumb … … 69 71 real :: pdsigma, pdryvel, pweightmolar, pohreact, pspec_ass, pkao 70 72 real :: pweta_in, pwetb_in, pwetc_in, pwetd_in 71 integer :: readerror73 integer :: ios 72 74 73 75 ! declare namelist … … 100 102 ! Open the SPECIES file and read species names and properties 101 103 !************************************************************ 102 specnum( pos_spec)=id_spec103 write(aspecnumb,'(i3.3)') specnum( pos_spec)104 specnum(id_pos)=id_spec 105 write(aspecnumb,'(i3.3)') specnum(id_pos) 104 106 open(unitspecies,file=path(1)(1:length(1))//'SPECIES/SPECIES_'//aspecnumb,status='old',form='formatted',err=998) 105 !write(*,*) 'reading SPECIES',specnum( pos_spec)107 !write(*,*) 'reading SPECIES',specnum(id_pos) 106 108 107 109 ASSSPEC=.FALSE. 108 110 109 111 ! try namelist input 110 read(unitspecies,species_params,iostat= readerror)112 read(unitspecies,species_params,iostat=ios) 111 113 close(unitspecies) 112 114 113 if ((pweightmolar.eq.-789.0).or.(readerror.ne.0)) then ! no namelist found 114 115 readerror=1 116 117 open(unitspecies,file=path(1)(1:length(1))//'SPECIES/SPECIES_'//aspecnumb,status='old',err=998) 115 if ((pweightmolar.eq.-789.0).or.(ios.ne.0)) then ! no namelist found 116 117 open(unitspecies,file=path(1) (1:length(1))//'SPECIES/SPECIES_'//aspecnumb, status='old', err=998) 118 118 119 119 do i=1,6 120 120 read(unitspecies,*) 121 121 end do 122 123 read(unitspecies,'(a10)',end=22) species(pos_spec) 124 ! write(*,*) species(pos_spec) 125 read(unitspecies,'(f18.1)',end=22) decay(pos_spec) 126 ! write(*,*) decay(pos_spec) 127 read(unitspecies,'(e18.1)',end=22) weta(pos_spec) 128 ! write(*,*) weta(pos_spec) 129 read(unitspecies,'(f18.2)',end=22) wetb(pos_spec) 130 ! write(*,*) wetb(pos_spec) 122 read(unitspecies,'(a10)',end=20) species(id_pos) 123 ! write(*,*) species(id_pos) 124 read(unitspecies,'(f18.1)',end=20) decay(id_pos) 125 ! write(*,*) decay(id_pos) 126 read(unitspecies,'(e18.1)',end=20) weta(id_pos) 127 ! write(*,*) weta(id_pos) 128 read(unitspecies,'(f18.2)',end=20) wetb(id_pos) 129 ! write(*,*) wetb(id_pos) 131 130 132 131 !*** NIK 31.01.2013: including in-cloud scavening parameters 133 read(unitspecies,'(e18.1)',end=2 2) weta_in(pos_spec)134 ! write(*,*) weta_in( pos_spec)135 read(unitspecies,'(f18.2)',end=2 2) wetb_in(pos_spec)136 ! write(*,*) wetb_in( pos_spec)137 read(unitspecies,'(f18.2)',end=2 2) wetc_in(pos_spec)138 ! write(*,*) wetc_in( pos_spec)139 read(unitspecies,'(f18.2)',end=2 2) wetd_in(pos_spec)140 ! write(*,*) wetd_in( pos_spec)141 142 read(unitspecies,'(f18.1)',end=2 2) reldiff(pos_spec)143 ! write(*,*) reldiff( pos_spec)144 read(unitspecies,'(e18.1)',end=2 2) henry(pos_spec)145 ! write(*,*) henry( pos_spec)146 read(unitspecies,'(f18.1)',end=2 2) f0(pos_spec)147 ! write(*,*) f0( pos_spec)148 read(unitspecies,'(e18.1)',end=2 2) density(pos_spec)149 ! write(*,*) density( pos_spec)150 read(unitspecies,'(e18.1)',end=2 2) dquer(pos_spec)151 ! write(*,*) dquer( pos_spec)152 read(unitspecies,'(e18.1)',end=2 2) dsigma(pos_spec)153 ! write(*,*) dsigma( pos_spec)154 read(unitspecies,'(f18.2)',end=2 2) dryvel(pos_spec)155 ! write(*,*) dryvel( pos_spec)156 read(unitspecies,'(f18.2)',end=2 2) weightmolar(pos_spec)157 ! write(*,*) weightmolar( pos_spec)158 read(unitspecies,'(e18.1)',end=2 2) ohreact(pos_spec)159 ! write(*,*) ohreact( pos_spec)160 read(unitspecies,'(i18)',end=2 2) spec_ass(pos_spec)161 ! write(*,*) spec_ass( pos_spec)162 read(unitspecies,'(f18.2)',end=2 2) kao(pos_spec)163 ! write(*,*) kao( pos_spec)164 165 pspecies=species( pos_spec)166 pdecay=decay( pos_spec)167 pweta=weta( pos_spec)168 pwetb=wetb( pos_spec)169 pweta_in=weta_in( pos_spec)170 pwetb_in=wetb_in( pos_spec)171 pwetc_in=wetc_in( pos_spec)172 pwetd_in=wetd_in( pos_spec)173 preldiff=reldiff( pos_spec)174 phenry=henry( pos_spec)175 pf0=f0( pos_spec)176 pdensity=density( pos_spec)177 pdquer=dquer( pos_spec)178 pdsigma=dsigma( pos_spec)179 pdryvel=dryvel( pos_spec)180 pweightmolar=weightmolar( pos_spec)181 pohreact=ohreact( pos_spec)182 pspec_ass=spec_ass( pos_spec)183 pkao=kao( pos_spec)132 read(unitspecies,'(e18.1)',end=20) weta_in(id_pos) 133 ! write(*,*) weta_in(id_pos) 134 read(unitspecies,'(f18.2)',end=20) wetb_in(id_pos) 135 ! write(*,*) wetb_in(id_pos) 136 read(unitspecies,'(f18.2)',end=20) wetc_in(id_pos) 137 ! write(*,*) wetc_in(id_pos) 138 read(unitspecies,'(f18.2)',end=20) wetd_in(id_pos) 139 ! write(*,*) wetd_in(id_pos) 140 141 read(unitspecies,'(f18.1)',end=20) reldiff(id_pos) 142 ! write(*,*) reldiff(id_pos) 143 read(unitspecies,'(e18.1)',end=20) henry(id_pos) 144 ! write(*,*) henry(id_pos) 145 read(unitspecies,'(f18.1)',end=20) f0(id_pos) 146 ! write(*,*) f0(id_pos) 147 read(unitspecies,'(e18.1)',end=20) density(id_pos) 148 ! write(*,*) density(id_pos) 149 read(unitspecies,'(e18.1)',end=20) dquer(id_pos) 150 ! write(*,*) dquer(id_pos) 151 read(unitspecies,'(e18.1)',end=20) dsigma(id_pos) 152 ! write(*,*) dsigma(id_pos) 153 read(unitspecies,'(f18.2)',end=20) dryvel(id_pos) 154 ! write(*,*) dryvel(id_pos) 155 read(unitspecies,'(f18.2)',end=20) weightmolar(id_pos) 156 ! write(*,*) weightmolar(id_pos) 157 read(unitspecies,'(e18.1)',end=20) ohreact(id_pos) 158 ! write(*,*) ohreact(id_pos) 159 read(unitspecies,'(i18)',end=20) spec_ass(id_pos) 160 ! write(*,*) spec_ass(id_pos) 161 read(unitspecies,'(f18.2)',end=20) kao(id_pos) 162 ! write(*,*) kao(id_pos) 163 164 pspecies=species(id_pos) 165 pdecay=decay(id_pos) 166 pweta=weta(id_pos) 167 pwetb=wetb(id_pos) 168 pweta_in=weta_in(id_pos) 169 pwetb_in=wetb_in(id_pos) 170 pwetc_in=wetc_in(id_pos) 171 pwetd_in=wetd_in(id_pos) 172 preldiff=reldiff(id_pos) 173 phenry=henry(id_pos) 174 pf0=f0(id_pos) 175 pdensity=density(id_pos) 176 pdquer=dquer(id_pos) 177 pdsigma=dsigma(id_pos) 178 pdryvel=dryvel(id_pos) 179 pweightmolar=weightmolar(id_pos) 180 pohreact=ohreact(id_pos) 181 pspec_ass=spec_ass(id_pos) 182 pkao=kao(id_pos) 184 183 185 184 else 186 185 187 species(pos_spec)=pspecies 188 decay(pos_spec)=pdecay 189 weta(pos_spec)=pweta 190 wetb(pos_spec)=pwetb 191 weta_in(pos_spec)=pweta_in 192 wetb_in(pos_spec)=pwetb_in 193 wetc_in(pos_spec)=pwetc_in 194 wetd_in(pos_spec)=pwetd_in 195 reldiff(pos_spec)=preldiff 196 henry(pos_spec)=phenry 197 f0(pos_spec)=pf0 198 density(pos_spec)=pdensity 199 dquer(pos_spec)=pdquer 200 dsigma(pos_spec)=pdsigma 201 dryvel(pos_spec)=pdryvel 202 weightmolar(pos_spec)=pweightmolar 203 ohreact(pos_spec)=pohreact 204 spec_ass(pos_spec)=pspec_ass 205 kao(pos_spec)=pkao 206 207 endif 208 209 i=pos_spec 210 211 if ((weta(pos_spec).gt.0).and.(henry(pos_spec).le.0)) then 212 if (dquer(pos_spec).le.0) goto 996 ! no particle, no henry set 213 endif 214 215 if (spec_ass(pos_spec).gt.0) then 186 species(id_pos)=pspecies 187 decay(id_pos)=pdecay 188 weta(id_pos)=pweta 189 wetb(id_pos)=pwetb 190 weta_in(id_pos)=pweta_in 191 wetb_in(id_pos)=pwetb_in 192 wetc_in(id_pos)=pwetc_in 193 wetd_in(id_pos)=pwetd_in 194 reldiff(id_pos)=preldiff 195 henry(id_pos)=phenry 196 f0(id_pos)=pf0 197 density(id_pos)=pdensity 198 dquer(id_pos)=pdquer 199 dsigma(id_pos)=pdsigma 200 dryvel(id_pos)=pdryvel 201 weightmolar(id_pos)=pweightmolar 202 ohreact(id_pos)=pohreact 203 spec_ass(id_pos)=pspec_ass 204 kao(id_pos)=pkao 205 206 endif 207 208 if ((weta(id_pos).gt.0).and.(henry(id_pos).le.0)) then 209 if (dquer(id_pos).le.0) goto 996 ! no particle, no henry set 210 endif 211 212 if (spec_ass(id_pos).gt.0) then 216 213 spec_found=.FALSE. 217 do j=1, pos_spec-1218 if (spec_ass( pos_spec).eq.specnum(j)) then219 spec_ass( pos_spec)=j214 do j=1,id_pos-1 215 if (spec_ass(id_pos).eq.specnum(j)) then 216 spec_ass(id_pos)=j 220 217 spec_found=.TRUE. 221 218 ASSSPEC=.TRUE. … … 227 224 endif 228 225 229 if (dsigma(i ).eq.1.) dsigma(i)=1.0001 ! avoid floating exception230 if (dsigma(i ).eq.0.) dsigma(i)=1.0001 ! avoid floating exception231 232 if ((reldiff(i ).gt.0.).and.(density(i).gt.0.)) then226 if (dsigma(id_pos).eq.1.) dsigma(id_pos)=1.0001 ! avoid floating exception 227 if (dsigma(id_pos).eq.0.) dsigma(id_pos)=1.0001 ! avoid floating exception 228 229 if ((reldiff(id_pos).gt.0.).and.(density(id_pos).gt.0.)) then 233 230 write(*,*) '#### FLEXPART MODEL ERROR! FILE "SPECIES" ####' 234 231 write(*,*) '#### IS CORRUPT. SPECIES CANNOT BE BOTH ####' … … 246 243 247 244 do j=1,24 ! initialize everything to no variation 248 area_hour(i ,j)=1.249 point_hour(i ,j)=1.245 area_hour(id_pos,j)=1. 246 point_hour(id_pos,j)=1. 250 247 end do 251 248 do j=1,7 252 area_dow(i ,j)=1.253 point_dow(i ,j)=1.249 area_dow(id_pos,j)=1. 250 point_dow(id_pos,j)=1. 254 251 end do 255 252 256 if (readerror.ne.0) then ! text format input 257 258 read(unitspecies,*,end=22) 253 if (ios.ne.0) then ! text format input 254 255 read(unitspecies,*,iostat=ios) 256 if (ios .ne. 0) goto 22 ! ifort does not like err= 259 257 do j=1,24 ! 24 hours, starting with 0-1 local time 260 read(unitspecies,*) ihour,area_hour(i ,j),point_hour(i,j)258 read(unitspecies,*) ihour,area_hour(id_pos,j),point_hour(id_pos,j) 261 259 end do 262 260 read(unitspecies,*) 263 261 do j=1,7 ! 7 days of the week, starting with Monday 264 read(unitspecies,*) idow,area_dow(i ,j),point_dow(i,j)262 read(unitspecies,*) idow,area_dow(id_pos,j),point_dow(id_pos,j) 265 263 end do 266 264
Note: See TracChangeset
for help on using the changeset viewer.