Changeset 2625


Ignore:
Timestamp:
Feb 22, 2022, 3:21:24 PM (3 years ago)
Author:
gmilcareck
Message:

GENERIC GCM
Minor bug fix for aerave_new.F when input wavelenght in data file are in descending order.
GM

Location:
trunk/LMDZ.GENERIC
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • trunk/LMDZ.GENERIC/README

    r2603 r2625  
    17111711Minor bug fix for the (rather specific) case when there is no surface and
    17121712inertiedat() is not defined.
     1713
     1714== 22/02/2022 == GM
     1715Minor bug fix for aerave_new.F when input wavelenght in data file are in
     1716descending order.

  • trunk/LMDZ.GENERIC/libf/phystd/aerave_new.F

    r135 r2625  
    5555c.......................................................................
    5656c
    57       INTEGER iir,nirmx
    58       PARAMETER (nirmx=100)
     57      INTEGER iir
     58      INTEGER,PARAMETER :: nirmx=100
    5959      INTEGER idata,ndata
    6060c
    … …  
    6666     &    ,omegdata(ndata),gdata(ndata)
    6767      REAL qextcorrdata(ndata)
    68       INTEGER ibande,nbande
    69       PARAMETER (nbande=1000)
     68      INTEGER ibande
     69      INTEGER,PARAMETER :: nbande=1000
    7070      REAL long,deltalong
    7171      INTEGER ilong
    … …  
    8383      long=longref
    8484
    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
     85c 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
     99c
     100     
     101       
    110102
    111103
    112104c********************************************************
    113105c 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))
     106c 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
     129c 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
    128152      ENDIF
     153
    129154c********************************************************
    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.......................................................................
     155c.......................................................................
     156c wavelengths (longdata) from data file in ascending order
     157c.......................................................................
     158      IF (longdata(1).LT.longdata(ndata)) THEN
     159        DO iir=1,nir
     160c
     161c.......................................................................
     162c
     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
     168c
     169c.......................................................................
     170c
     171          DO ibande=1,nbande
     172c
     173c.......................................................................
     174c
     175            long=longir(iir) + (ibande-0.5E+0) * deltalong
     176            CALL blackl(DBLE(long),DBLE(temp),tmp1)
     177            emit=REAL(tmp1)
     178c
     179c.......................................................................
     180c
     181c 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)
     200c
     201c.......................................................................
     202c
     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
     207c
     208c.......................................................................
     209c
     210          ENDDO
     211c
     212c.......................................................................
     213c
     214          gir(iir)=gir(iir)/omegir(iir)
     215          omegir(iir)=omegir(iir)/epir(iir)
     216          epir(iir)=epir(iir)/totalemit(iir)
     217c
     218c.......................................................................
     219c
     220        ENDDO
     221c.......................................................................
     222c wavelengths (longdata) from data file in descending order
     223c.......................................................................
     224      ELSEIF (longdata(1).GT.longdata(ndata)) THEN
     225        DO iir=1,nir
     226c
     227c.......................................................................
     228c
     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
     234c
     235c.......................................................................
     236c
     237          DO ibande=1,nbande
     238c
     239c.......................................................................
     240c
     241            long=longir(iir) + (ibande-0.5E+0) * deltalong
     242            CALL blackl(DBLE(long),DBLE(temp),tmp1)
     243            emit=REAL(tmp1)
     244c
     245c.......................................................................
     246c
     247c 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)
     266c
     267c.......................................................................
     268c
     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
     273c
     274c.......................................................................
     275c
     276          ENDDO
     277c
     278c.......................................................................
     279c
     280          gir(iir)=gir(iir)/omegir(iir)
     281          omegir(iir)=omegir(iir)/epir(iir)
     282          epir(iir)=epir(iir)/totalemit(iir)
     283c
     284c.......................................................................
     285c
     286        ENDDO
     287      ENDIF
    161288c
    162289c********************************************************
    163 c interpolation
    164       ilong=1
    165       DO idata=2,ndata
    166         IF (long.gt.longdata(idata)) ilong=idata
    167       ENDDO
    168       i1=ilong
    169       i2=ilong+1
    170       IF (i2.gt.ndata) i2=ndata
    171       IF (long.lt.longdata(1)) i2=1
    172       IF (i1.eq.i2) THEN
    173         c1=1.E+0
    174         c2=0.E+0
    175       ELSE
    176         c1=(longdata(i2)-long) / (longdata(i2)-longdata(i1))
    177         c2=(longdata(i1)-long) / (longdata(i1)-longdata(i2))
    178       ENDIF
    179 c********************************************************
    180 c
    181           qextcorr=c1*qextcorrdata(i1)+c2*qextcorrdata(i2)
    182           omeg=c1*omegdata(i1)+c2*omegdata(i2)
    183           g=c1*gdata(i1)+c2*gdata(i2)
    184 c
    185 c.......................................................................
    186 c
    187           totalemit(iir)=totalemit(iir)+deltalong*emit
    188           epir(iir)=epir(iir)+deltalong*emit*qextcorr
    189           omegir(iir)=omegir(iir)+deltalong*emit*omeg*qextcorr
    190           gir(iir)=gir(iir)+deltalong*emit*omeg*qextcorr*g
    191 c
    192 c.......................................................................
    193 c
    194         ENDDO
    195 c
    196 c.......................................................................
    197 c
    198         gir(iir)=gir(iir)/omegir(iir)
    199         omegir(iir)=omegir(iir)/epir(iir)
    200         epir(iir)=epir(iir)/totalemit(iir)
    201 c
    202 c.......................................................................
    203 c
    204       ENDDO
    205290c
    206291c......................................................................
    … …  
    220305c......................................................................
    221306c
    222       RETURN
    223307      END
Note: See TracChangeset for help on using the changeset viewer.