Changeset 2625 for trunk/LMDZ.GENERIC/libf/phystd/aerave_new.F
- Timestamp:
- Feb 22, 2022, 3:21:24 PM (3 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/LMDZ.GENERIC/libf/phystd/aerave_new.F
r135 r2625 55 55 c....................................................................... 56 56 c 57 INTEGER iir ,nirmx58 PARAMETER (nirmx=100)57 INTEGER iir 58 INTEGER,PARAMETER :: nirmx=100 59 59 INTEGER idata,ndata 60 60 c … … 66 66 & ,omegdata(ndata),gdata(ndata) 67 67 REAL qextcorrdata(ndata) 68 INTEGER ibande ,nbande69 PARAMETER (nbande=1000)68 INTEGER ibande 69 INTEGER,PARAMETER :: nbande=1000 70 70 REAL long,deltalong 71 71 INTEGER ilong … … 83 83 long=longref 84 84 85 86 !if(nir.eq.27)then 87 !print*,'long',long 88 !print*,'longdata',longdata 89 !print*,'epdata',epdata 90 !print*,'omegdata',omegdata 91 !print*,'gdata',gdata 92 !print*,'data looks aok!' 93 94 !print*,'ndata=',ndata 95 !print*,'longdata',shape(longdata) 96 !print*,'epdata',shape(epdata) 97 !print*,'omegdata',shape(omegdata) 98 !print*,'gdata',shape(gdata) 99 ! print*,'longref',longref 100 !print*,'epref',epref 101 !print*,'temp',temp 102 !print*,'nir',nir 103 !print*,'longir',longir 104 !print*,'epir',epir 105 !print*,'omegir',gir 106 !print*,'qref',qref 107 !print*,'longir=',longir 108 !stop 109 !endif 85 c check ordering of longdata 86 DO idata=1,ndata-1 87 IF (longdata(1).LT.longdata(ndata)) THEN 88 IF (.not.(longdata(idata).LT.longdata(idata+1))) THEN 89 call abort_physic("aerave_new", 90 & "Non descending order in longdata",1) 91 ENDIF 92 ELSEIF (longdata(1).GT.longdata(ndata)) THEN 93 IF (.not.(longdata(idata).GT.longdata(idata+1))) THEN 94 call abort_physic("aerave_new", 95 & "Non ascending order in longdata",1) 96 ENDIF 97 ENDIF 98 ENDDO 99 c 100 101 110 102 111 103 112 104 c******************************************************** 113 105 c interpolation 114 ilong=1 115 DO idata=2,ndata 116 IF (long.gt.longdata(idata)) ilong=idata 117 ENDDO 118 i1=ilong 119 i2=ilong+1 120 IF (i2.gt.ndata) i2=ndata 121 IF (long.lt.longdata(1)) i2=1 122 IF (i1.eq.i2) THEN 123 c1=1.E+0 124 c2=0.E+0 125 ELSE 126 c1=(longdata(i2)-long) / (longdata(i2)-longdata(i1)) 127 c2=(longdata(i1)-long) / (longdata(i1)-longdata(i2)) 106 c wavelengths (longdata) from data file in ascending order 107 IF (longdata(1).LT.longdata(ndata)) THEN 108 ilong=1 109 DO idata=2,ndata 110 IF (long.gt.longdata(idata)) ilong=idata 111 ENDDO 112 i1=ilong 113 i2=ilong+1 114 IF (i2.gt.ndata) i2=ndata 115 IF (long.lt.longdata(1)) i2=1 116 IF (i1.eq.i2) THEN 117 c1=1.E+0 118 c2=0.E+0 119 ELSE 120 c1=(longdata(i2)-long) / (longdata(i2)-longdata(i1)) 121 c2=(longdata(i1)-long) / (longdata(i1)-longdata(i2)) 122 ENDIF 123 qref=c1*epdata(i1)+c2*epdata(i2) 124 omegaref=c1*omegdata(i1)+c2*omegdata(i2) 125 factep=qref/epref 126 DO idata=1,ndata 127 qextcorrdata(idata)=epdata(idata)/factep 128 ENDDO 129 c wavelengths (longdata) from data file in descending order 130 ELSEIF (longdata(1).GT.longdata(ndata)) THEN 131 ilong=1 132 DO idata=2,ndata 133 IF (long.lt.longdata(idata)) ilong=idata 134 ENDDO 135 i1=ilong+1 136 i2=ilong 137 IF (i1.gt.ndata) i1=ndata 138 IF (long.gt.longdata(1)) i1=1 139 IF (i1.eq.i2) THEN 140 c1=1.E+0 141 c2=0.E+0 142 ELSE 143 c1=(longdata(i2)-long) / (longdata(i2)-longdata(i1)) 144 c2=(longdata(i1)-long) / (longdata(i1)-longdata(i2)) 145 ENDIF 146 qref=c1*epdata(i1)+c2*epdata(i2) 147 omegaref=c1*omegdata(i1)+c2*omegdata(i2) 148 factep=qref/epref 149 DO idata=1,ndata 150 qextcorrdata(idata)=epdata(idata)/factep 151 ENDDO 128 152 ENDIF 153 129 154 c******************************************************** 130 c 131 qref=c1*epdata(i1)+c2*epdata(i2) 132 omegaref=c1*omegdata(i1)+c2*omegdata(i2) 133 factep=qref/epref 134 DO idata=1,ndata 135 qextcorrdata(idata)=epdata(idata)/factep 136 ENDDO 137 c 138 c....................................................................... 139 c 140 DO iir=1,nir 141 c 142 c....................................................................... 143 c 144 deltalong=(longir(iir+1)-longir(iir)) / nbande 145 totalemit(iir)=0.E+0 146 epir(iir)=0.E+0 147 omegir(iir)=0.E+0 148 gir(iir)=0.E+0 149 c 150 c....................................................................... 151 c 152 DO ibande=1,nbande 153 c 154 c....................................................................... 155 c 156 long=longir(iir) + (ibande-0.5E+0) * deltalong 157 CALL blackl(DBLE(long),DBLE(temp),tmp1) 158 emit=REAL(tmp1) 159 c 160 c....................................................................... 155 c....................................................................... 156 c wavelengths (longdata) from data file in ascending order 157 c....................................................................... 158 IF (longdata(1).LT.longdata(ndata)) THEN 159 DO iir=1,nir 160 c 161 c....................................................................... 162 c 163 deltalong=(longir(iir+1)-longir(iir)) / nbande 164 totalemit(iir)=0.E+0 165 epir(iir)=0.E+0 166 omegir(iir)=0.E+0 167 gir(iir)=0.E+0 168 c 169 c....................................................................... 170 c 171 DO ibande=1,nbande 172 c 173 c....................................................................... 174 c 175 long=longir(iir) + (ibande-0.5E+0) * deltalong 176 CALL blackl(DBLE(long),DBLE(temp),tmp1) 177 emit=REAL(tmp1) 178 c 179 c....................................................................... 180 c 181 c interpolation 182 ilong=1 183 DO idata=2,ndata 184 IF (long.gt.longdata(idata)) ilong=idata 185 ENDDO 186 i1=ilong 187 i2=ilong+1 188 IF (i2.gt.ndata) i2=ndata 189 IF (long.lt.longdata(1)) i2=1 190 IF (i1.eq.i2) THEN 191 c1=1.E+0 192 c2=0.E+0 193 ELSE 194 c1=(longdata(i2)-long) / (longdata(i2)-longdata(i1)) 195 c2=(longdata(i1)-long) / (longdata(i1)-longdata(i2)) 196 ENDIF 197 qextcorr=c1*qextcorrdata(i1)+c2*qextcorrdata(i2) 198 omeg=c1*omegdata(i1)+c2*omegdata(i2) 199 g=c1*gdata(i1)+c2*gdata(i2) 200 c 201 c....................................................................... 202 c 203 totalemit(iir)=totalemit(iir)+deltalong*emit 204 epir(iir)=epir(iir)+deltalong*emit*qextcorr 205 omegir(iir)=omegir(iir)+deltalong*emit*omeg*qextcorr 206 gir(iir)=gir(iir)+deltalong*emit*omeg*qextcorr*g 207 c 208 c....................................................................... 209 c 210 ENDDO 211 c 212 c....................................................................... 213 c 214 gir(iir)=gir(iir)/omegir(iir) 215 omegir(iir)=omegir(iir)/epir(iir) 216 epir(iir)=epir(iir)/totalemit(iir) 217 c 218 c....................................................................... 219 c 220 ENDDO 221 c....................................................................... 222 c wavelengths (longdata) from data file in descending order 223 c....................................................................... 224 ELSEIF (longdata(1).GT.longdata(ndata)) THEN 225 DO iir=1,nir 226 c 227 c....................................................................... 228 c 229 deltalong=(longir(iir+1)-longir(iir)) / nbande 230 totalemit(iir)=0.E+0 231 epir(iir)=0.E+0 232 omegir(iir)=0.E+0 233 gir(iir)=0.E+0 234 c 235 c....................................................................... 236 c 237 DO ibande=1,nbande 238 c 239 c....................................................................... 240 c 241 long=longir(iir) + (ibande-0.5E+0) * deltalong 242 CALL blackl(DBLE(long),DBLE(temp),tmp1) 243 emit=REAL(tmp1) 244 c 245 c....................................................................... 246 c 247 c interpolation 248 ilong=1 249 DO idata=2,ndata 250 IF (long.lt.longdata(idata)) ilong=idata 251 ENDDO 252 i1=ilong+1 253 i2=ilong 254 IF (i1.gt.ndata) i1=ndata 255 IF (long.gt.longdata(1)) i1=1 256 IF (i1.eq.i2) THEN 257 c1=1.E+0 258 c2=0.E+0 259 ELSE 260 c1=(longdata(i2)-long) / (longdata(i2)-longdata(i1)) 261 c2=(longdata(i1)-long) / (longdata(i1)-longdata(i2)) 262 ENDIF 263 qextcorr=c1*qextcorrdata(i1)+c2*qextcorrdata(i2) 264 omeg=c1*omegdata(i1)+c2*omegdata(i2) 265 g=c1*gdata(i1)+c2*gdata(i2) 266 c 267 c....................................................................... 268 c 269 totalemit(iir)=totalemit(iir)+deltalong*emit 270 epir(iir)=epir(iir)+deltalong*emit*qextcorr 271 omegir(iir)=omegir(iir)+deltalong*emit*omeg*qextcorr 272 gir(iir)=gir(iir)+deltalong*emit*omeg*qextcorr*g 273 c 274 c....................................................................... 275 c 276 ENDDO 277 c 278 c....................................................................... 279 c 280 gir(iir)=gir(iir)/omegir(iir) 281 omegir(iir)=omegir(iir)/epir(iir) 282 epir(iir)=epir(iir)/totalemit(iir) 283 c 284 c....................................................................... 285 c 286 ENDDO 287 ENDIF 161 288 c 162 289 c******************************************************** 163 c interpolation164 ilong=1165 DO idata=2,ndata166 IF (long.gt.longdata(idata)) ilong=idata167 ENDDO168 i1=ilong169 i2=ilong+1170 IF (i2.gt.ndata) i2=ndata171 IF (long.lt.longdata(1)) i2=1172 IF (i1.eq.i2) THEN173 c1=1.E+0174 c2=0.E+0175 ELSE176 c1=(longdata(i2)-long) / (longdata(i2)-longdata(i1))177 c2=(longdata(i1)-long) / (longdata(i1)-longdata(i2))178 ENDIF179 c********************************************************180 c181 qextcorr=c1*qextcorrdata(i1)+c2*qextcorrdata(i2)182 omeg=c1*omegdata(i1)+c2*omegdata(i2)183 g=c1*gdata(i1)+c2*gdata(i2)184 c185 c.......................................................................186 c187 totalemit(iir)=totalemit(iir)+deltalong*emit188 epir(iir)=epir(iir)+deltalong*emit*qextcorr189 omegir(iir)=omegir(iir)+deltalong*emit*omeg*qextcorr190 gir(iir)=gir(iir)+deltalong*emit*omeg*qextcorr*g191 c192 c.......................................................................193 c194 ENDDO195 c196 c.......................................................................197 c198 gir(iir)=gir(iir)/omegir(iir)199 omegir(iir)=omegir(iir)/epir(iir)200 epir(iir)=epir(iir)/totalemit(iir)201 c202 c.......................................................................203 c204 ENDDO205 290 c 206 291 c...................................................................... … … 220 305 c...................................................................... 221 306 c 222 RETURN223 307 END
Note: See TracChangeset
for help on using the changeset viewer.