Ignore:
Timestamp:
Jun 3, 2011, 7:28:17 PM (13 years ago)
Author:
musat
Message:

Ajouts CFMIP2/CMIP5

  • 6eme fichier de sortie "stations" histstn.nc qui necessite 2 fichiers: PARAM/npCFMIP_param.data contenant le nombre de points (120 pour simulations AMIP, 73 pour aqua) PARAM/pointlocations.txt contenat le numero, les coordonnees (lon,lat) et le nom de chaque station
  • flag LOGICAL dans tous les appels histwrite_phy pour pouvoir sortir le fichier histstn.nc

NB: 1) les flags de type phys_ que l'on met dans le physiq.def_L* pour ajouter plus de sorties

necessitent dorenavant 6 valeurs, la 6eme correspondant au fichier histstn.nc

2) par defaut le fichier histstn.nc ne sort pas; pour le sortir ajouter les lignes suivantes

dans physiq.def_L*

### Type de fichier : global (n) ou stations (y)
phys_out_filestations = n n n n n y

  • introduction de 2 jeux de flags pour les taux des GES; taux actuels avec suffixes _act, taux futurs avec "_per" avec 2 appels au rayonnement si taux "_per" different des taux "_act" (utiles pour diags. CFMIP 4CO2)
  • flags "betaCRF" pour calculs CRF pour experiences sensibilite proprietes optiques eau liquide nuageuse avec initialisations par defaut; sinon besoin de fichier beta_crf.data

IM

File:
1 edited

Legend:

Unmodified
Added
Removed

  • LMDZ4/branches/LMDZ4_AR5/libf/phylmd/iophy.F90

    r1331 r1534  
    99  REAL,allocatable,dimension(:),save :: io_lon
    1010  INTEGER, save :: phys_domain_id
     11  INTEGER, save :: npstn
     12  INTEGER, allocatable, dimension(:), save :: nptabij
    1113 
    1214  INTERFACE histwrite_phy
    1315    MODULE PROCEDURE histwrite2d_phy,histwrite3d_phy
     16  END INTERFACE
     17
     18  INTERFACE histbeg_phy_all
     19    MODULE PROCEDURE histbeg_phy,histbeg_phy_points
    1420  END INTERFACE
    1521
    … …  
    144150 
    145151  end subroutine histbeg_phy
    146  
    147   subroutine histwrite2d_phy(nid,name,itau,field)
     152
     153  subroutine histbeg_phy_points(rlon,rlat,pim,tabij,ipt,jpt, &
     154             plon,plat,plon_bounds,plat_bounds, &
     155             nname,itau0,zjulian,dtime,nnhori,nnid_day)
     156  USE dimphy
     157  USE mod_phys_lmdz_para
     158  USE mod_grid_phy_lmdz
     159  use ioipsl
     160  use write_field
     161  implicit none
     162  include 'dimensions.h'
     163
     164    real,dimension(klon),intent(in) :: rlon
     165    real,dimension(klon),intent(in) :: rlat
     166    integer, intent(in) :: itau0
     167    real,intent(in) :: zjulian
     168    real,intent(in) :: dtime
     169    integer, intent(in) :: pim
     170    integer, intent(out) :: nnhori
     171    character(len=20), intent(in) :: nname
     172    INTEGER, intent(out) :: nnid_day
     173    integer :: i
     174    REAL,dimension(klon_glo)        :: rlat_glo
     175    REAL,dimension(klon_glo)        :: rlon_glo
     176    INTEGER, DIMENSION(pim), intent(in)  :: tabij
     177    REAL,dimension(pim), intent(in) :: plat, plon
     178    INTEGER,dimension(pim), intent(in) :: ipt, jpt
     179    REAL,dimension(pim,2), intent(out) :: plat_bounds, plon_bounds
     180
     181    INTEGER, SAVE :: tabprocbeg, tabprocend
     182!$OMP THREADPRIVATE(tabprocbeg, tabprocend)
     183    INTEGER :: ip
     184    INTEGER, PARAMETER :: nip=1
     185    INTEGER :: npproc
     186    REAL, allocatable, dimension(:) :: npplat, npplon
     187    REAL, allocatable, dimension(:,:) :: npplat_bounds, npplon_bounds
     188    INTEGER, PARAMETER :: jjmp1=jjm+1-1/jjm
     189    REAL, dimension(iim,jjmp1) :: zx_lon, zx_lat
     190
     191    CALL gather(rlat,rlat_glo)
     192    CALL bcast(rlat_glo)
     193    CALL gather(rlon,rlon_glo)
     194    CALL bcast(rlon_glo)
     195
     196!$OMP MASTER
     197    DO i=1,pim
     198
     199!    print*,'CFMIP_iophy i tabij lon lat',i,tabij(i),plon(i),plat(i)
     200
     201     plon_bounds(i,1)=rlon_glo(tabij(i)-1)
     202     plon_bounds(i,2)=rlon_glo(tabij(i)+1)
     203     if(plon_bounds(i,2).LE.0..AND.plon_bounds(i,1).GE.0.) THEN
     204      if(rlon_glo(tabij(i)).GE.0.) THEN
     205       plon_bounds(i,2)=-1*plon_bounds(i,2)
     206      endif
     207     endif
     208     if(plon_bounds(i,2).GE.0..AND.plon_bounds(i,1).LE.0.) THEN
     209      if(rlon_glo(tabij(i)).LE.0.) THEN
     210       plon_bounds(i,2)=-1*plon_bounds(i,2)
     211      endif
     212     endif
     213!
     214     IF ( tabij(i).LE.iim) THEN
     215      plat_bounds(i,1)=rlat_glo(tabij(i))
     216     ELSE
     217      plat_bounds(i,1)=rlat_glo(tabij(i)-iim)
     218     ENDIF
     219     plat_bounds(i,2)=rlat_glo(tabij(i)+iim)
     220!
     221!    print*,'CFMIP_iophy point i lon lon_bds',i,plon_bounds(i,1),rlon_glo(tabij(i)),plon_bounds(i,2)
     222!    print*,'CFMIP_iophy point i lat lat_bds',i,plat_bounds(i,1),rlat_glo(tabij(i)),plat_bounds(i,2)
     223!
     224    ENDDO
     225    if (is_sequential) then
     226
     227     npstn=pim
     228     IF(.NOT. ALLOCATED(nptabij)) THEN
     229      ALLOCATE(nptabij(pim))
     230     ENDIF
     231     DO i=1,pim
     232      nptabij(i)=tabij(i)
     233     ENDDO
     234
     235       CALL gr_fi_ecrit(1,klon,iim,jjmp1,rlon_glo,zx_lon)
     236       if (iim.gt.1) then
     237       DO i = 1, iim
     238         zx_lon(i,1) = rlon_glo(i+1)
     239         zx_lon(i,jjmp1) = rlon_glo(i+1)
     240       ENDDO
     241       endif
     242       CALL gr_fi_ecrit(1,klon,iim,jjmp1,rlat_glo,zx_lat)
     243
     244    DO i=1,pim
     245!    print*,'CFMIP_iophy i tabij lon lat',i,tabij(i),plon(i),plat(i)
     246
     247     plon_bounds(i,1)=zx_lon(ipt(i)-1,jpt(i))
     248     plon_bounds(i,2)=zx_lon(ipt(i)+1,jpt(i))
     249
     250     if (ipt(i).EQ.1) then
     251      plon_bounds(i,1)=zx_lon(iim,jpt(i))
     252      plon_bounds(i,2)=360.+zx_lon(ipt(i)+1,jpt(i))
     253     endif
     254 
     255     if (ipt(i).EQ.iim) then
     256      plon_bounds(i,2)=360.+zx_lon(1,jpt(i))
     257     endif
     258
     259     plat_bounds(i,1)=zx_lat(ipt(i),jpt(i)-1)
     260     plat_bounds(i,2)=zx_lat(ipt(i),jpt(i)+1)
     261
     262     if (jpt(i).EQ.1) then
     263      plat_bounds(i,1)=zx_lat(ipt(i),1)+0.001
     264      plat_bounds(i,2)=zx_lat(ipt(i),1)-0.001
     265     endif
     266 
     267     if (jpt(i).EQ.jjmp1) then
     268      plat_bounds(i,1)=zx_lat(ipt(i),jjmp1)+0.001
     269      plat_bounds(i,2)=zx_lat(ipt(i),jjmp1)-0.001
     270     endif
     271!
     272!    print*,'CFMIP_iophy point i lon lon_bds',i,plon_bounds(i,1),rlon(tabij(i)),plon_bounds(i,2)
     273!    print*,'CFMIP_iophy point i lat lat_bds',i,plat_bounds(i,1),rlat(tabij(i)),plat_bounds(i,2)
     274!
     275    ENDDO
     276!    print*,'iophy is_sequential nname, nnhori, nnid_day=',trim(nname),nnhori,nnid_day
     277     call histbeg(nname,pim,plon,plon_bounds, &
     278                           plat,plat_bounds, &
     279                           itau0, zjulian, dtime, nnhori, nnid_day)
     280    else
     281     npproc=0
     282     DO ip=1, pim
     283      tabprocbeg=klon_mpi_begin
     284      tabprocend=klon_mpi_end
     285      IF(tabij(ip).GE.tabprocbeg.AND.tabij(ip).LE.tabprocend) THEN
     286       npproc=npproc+1
     287       npstn=npproc
     288      ENDIF
     289     ENDDO
     290!    print*,'CFMIP_iophy mpi_rank npstn',mpi_rank,npstn
     291     IF(.NOT. ALLOCATED(nptabij)) THEN
     292      ALLOCATE(nptabij(npstn))
     293      ALLOCATE(npplon(npstn), npplat(npstn))
     294      ALLOCATE(npplon_bounds(npstn,2), npplat_bounds(npstn,2))
     295     ENDIF
     296     npproc=0
     297     DO ip=1, pim
     298      IF(tabij(ip).GE.tabprocbeg.AND.tabij(ip).LE.tabprocend) THEN
     299       npproc=npproc+1
     300       nptabij(npproc)=tabij(ip)
     301!      print*,'mpi_rank npproc ip plon plat tabij=',mpi_rank,npproc,ip, &
     302!      plon(ip),plat(ip),tabij(ip)
     303       npplon(npproc)=plon(ip)
     304       npplat(npproc)=plat(ip)
     305       npplon_bounds(npproc,1)=plon_bounds(ip,1)
     306       npplon_bounds(npproc,2)=plon_bounds(ip,2)
     307       npplat_bounds(npproc,1)=plat_bounds(ip,1)
     308       npplat_bounds(npproc,2)=plat_bounds(ip,2)
     309!!!
     310!!! print qui sert a reordonner les points stations selon l'ordre CFMIP
     311!!! ne pas enlever
     312        print*,'iophy_mpi rank ip lon lat',mpi_rank,ip,plon(ip),plat(ip)
     313!!!
     314      ENDIF
     315     ENDDO
     316     call histbeg(nname,npstn,npplon,npplon_bounds, &
     317                            npplat,npplat_bounds, &
     318                            itau0,zjulian,dtime,nnhori,nnid_day,phys_domain_id)
     319    endif
     320!$OMP END MASTER
     321
     322  end subroutine histbeg_phy_points
     323 
     324  subroutine histwrite2d_phy(nid,lpoint,name,itau,field)
    148325  USE dimphy
    149326  USE mod_phys_lmdz_para
    … …  
    153330   
    154331    integer,intent(in) :: nid
     332    logical,intent(in) :: lpoint
    155333    character*(*), intent(IN) :: name
    156334    integer, intent(in) :: itau
    157335    real,dimension(:),intent(in) :: field
    158336    REAL,dimension(klon_mpi) :: buffer_omp
    159     INTEGER :: index2d(iim*jj_nb)
     337    INTEGER, allocatable, dimension(:) :: index2d
    160338    REAL :: Field2d(iim,jj_nb)
     339
     340    integer :: ip
     341    real,allocatable,dimension(:) :: fieldok
    161342
    162343    IF (size(field)/=klon) CALL abort_gcm('iophy::histwrite2d','Field first dimension not equal to klon',1)
    … …  
    165346!$OMP MASTER
    166347    CALL grid1Dto2D_mpi(buffer_omp,Field2d)
    167     CALL histwrite(nid,name,itau,Field2d,iim*jj_nb,index2d)
     348    if(.NOT.lpoint) THEN
     349     ALLOCATE(index2d(iim*jj_nb))
     350     ALLOCATE(fieldok(iim*jj_nb))
     351     CALL histwrite(nid,name,itau,Field2d,iim*jj_nb,index2d)
     352    else
     353     ALLOCATE(fieldok(npstn))
     354     ALLOCATE(index2d(npstn))
     355
     356     if(is_sequential) then
     357!     klon_mpi_begin=1
     358!     klon_mpi_end=klon
     359      DO ip=1, npstn
     360       fieldok(ip)=buffer_omp(nptabij(ip))
     361      ENDDO
     362     else
     363      DO ip=1, npstn
     364!     print*,'histwrite2d is_sequential npstn ip name nptabij',npstn,ip,name,nptabij(ip)
     365       IF(nptabij(ip).GE.klon_mpi_begin.AND. &
     366          nptabij(ip).LE.klon_mpi_end) THEN
     367         fieldok(ip)=buffer_omp(nptabij(ip)-klon_mpi_begin+1)
     368       ENDIF
     369      ENDDO
     370     endif
     371     CALL histwrite(nid,name,itau,fieldok,npstn,index2d)
     372!
     373    endif
     374    deallocate(index2d)
     375    deallocate(fieldok)
    168376!$OMP END MASTER   
    169377  end subroutine histwrite2d_phy
    170378
    171 
    172  
    173   subroutine histwrite3d_phy(nid,name,itau,field)
     379  subroutine histwrite3d_phy(nid,lpoint,name,itau,field)
    174380  USE dimphy
    175381  USE mod_phys_lmdz_para
    … …  
    180386   
    181387    integer,intent(in) :: nid
     388    logical,intent(in) :: lpoint
    182389    character*(*), intent(IN) :: name
    183390    integer, intent(in) :: itau
    184391    real,dimension(:,:),intent(in) :: field  ! --> field(klon,:)
    185392    REAL,dimension(klon_mpi,size(field,2)) :: buffer_omp
    186     INTEGER :: nlev
    187     INTEGER :: index3d(iim*jj_nb*size(field,2))
    188393    REAL :: Field3d(iim,jj_nb,size(field,2))
    189    
     394    INTEGER :: ip, n, nlev
     395    INTEGER, ALLOCATABLE, dimension(:) :: index3d
     396    real,allocatable, dimension(:,:) :: fieldok
     397
    190398    IF (size(field,1)/=klon) CALL abort_gcm('iophy::histwrite3d','Field first dimension not equal to klon',1)
    191399    nlev=size(field,2)
    192    
     400
     401!   print*,'hist3d_phy mpi_rank npstn=',mpi_rank,npstn
     402
     403!   DO ip=1, npstn
     404!    print*,'hist3d_phy mpi_rank nptabij',mpi_rank,nptabij(ip)
     405!   ENDDO
     406
    193407    CALL Gather_omp(field,buffer_omp)
    194408!$OMP MASTER
    195409    CALL grid1Dto2D_mpi(buffer_omp,field3d)
    196     CALL histwrite(nid,name,itau,Field3d,iim*jj_nb*nlev,index3d)
     410    if(.NOT.lpoint) THEN
     411     ALLOCATE(index3d(iim*jj_nb*nlev))
     412     ALLOCATE(fieldok(iim*jj_nb,nlev))
     413     CALL histwrite(nid,name,itau,Field3d,iim*jj_nb*nlev,index3d)
     414    else
     415      nlev=size(field,2)
     416      ALLOCATE(index3d(npstn*nlev))
     417      ALLOCATE(fieldok(npstn,nlev))
     418
     419      if(is_sequential) then
     420!      klon_mpi_begin=1
     421!      klon_mpi_end=klon
     422       DO n=1, nlev
     423       DO ip=1, npstn
     424        fieldok(ip,n)=buffer_omp(nptabij(ip),n)
     425       ENDDO
     426       ENDDO
     427      else
     428       DO n=1, nlev
     429       DO ip=1, npstn
     430        IF(nptabij(ip).GE.klon_mpi_begin.AND. &
     431         nptabij(ip).LE.klon_mpi_end) THEN
     432         fieldok(ip,n)=buffer_omp(nptabij(ip)-klon_mpi_begin+1,n)
     433        ENDIF
     434       ENDDO
     435       ENDDO
     436      endif
     437      CALL histwrite(nid,name,itau,fieldok,npstn*nlev,index3d)
     438    endif
     439  deallocate(index3d)
     440  deallocate(fieldok)
    197441!$OMP END MASTER   
    198442  end subroutine histwrite3d_phy
    199443 
    200  
    201 
    202444end module iophy
Note: See TracChangeset for help on using the changeset viewer.