#define TRACEBACK write (gol,'("in ",a," (",a,i6,")")') rname, __FILE__, __LINE__ ; call goErr #define IF_NOTOK_RETURN(action) if (status/=0) then; TRACEBACK; action; return; end if #define IF_ERROR_RETURN(action) if (status> 0) then; TRACEBACK; action; return; end if ! #include "tm5.inc" ! ! ! Tools for climatologies ! module Climat use GO, only : TDate implicit none ! --- in/out ------------------------- private public :: TClimat public :: Init, Done public :: Set, Setup ! --- const ----------------------------- character(len=*), parameter :: mname = 'Climat' ! --- types -------------------------- type TClimat character(len=16) :: name, unit ! actual climatology: real, pointer :: data(:,:,:) ! current data can be used to interpolate to this interval: type(TDate) :: tr(2) ! fields used for temporal interpolation: character(len=10) :: tinterp integer :: nt type(TDate), pointer :: t_in(:) real, pointer :: data_in(:,:,:,:) end type TClimat ! --- interfaces ------------------------ interface Init module procedure climat_Init end interface interface Done module procedure climat_Done end interface interface Setup module procedure climat_Setup end interface interface Set module procedure climat_Set end interface contains ! ========================================================= subroutine climat_Init( climat, name, unit, tinterp, im, jm, lm, status ) use GO, only : gol, goErr use GO, only : AnyDate ! --- in/out -------------------------------- type(TClimat), intent(out) :: climat character(len=*), intent(in) :: name, unit character(len=*), intent(in) :: tinterp integer, intent(in) :: im, jm, lm integer, intent(inout) :: status ! --- const ------------------------------ character(len=*), parameter :: rname = mname//'/climat_Init' ! --- local ------------------------------- integer :: it ! --- begin -------------------------------- ! trap previous errors: if (status>0) return ! store info: climat%name = name climat%unit = unit climat%tinterp = tinterp ! setup target data: allocate( climat%data (im,jm,lm) ) ! data not valid yet: climat%tr(1) = AnyDate() climat%tr(2) = AnyDate() ! setup interpolation data if necessary: select case ( climat%tinterp ) case ( 'linear' ) climat%nt = 2 allocate( climat%data_in(im,jm,lm,climat%nt) ) allocate( climat%t_in(climat%nt) ) do it = 1, climat%nt climat%t_in(it) = AnyDate() end do case ( 'constant' ) climat%nt = 0 nullify( climat%data_in ) case default write (gol,'("unsupported time interpolation : ",a)') trim(climat%tinterp); call goErr write (gol,'("in ",a)') rname; call goErr; status=1; return end select ! ok status = 0 end subroutine climat_Init ! *** subroutine climat_Done( climat, status ) use GO, only : gol, goErr ! --- in/out -------------------------------- type(TClimat), intent(inout) :: climat integer, intent(inout) :: status ! --- const ------------------------------ character(len=*), parameter :: rname = mname//'/climat_Done' ! --- begin -------------------------------- ! trap previous errors: if (status>0) return ! clear arrays: deallocate( climat%data ) select case ( climat%tinterp ) case ( 'linear' ) deallocate( climat%data_in ) deallocate( climat%t_in ) case ( 'constant' ) ! nothing to clear case default write (gol,'("unsupported time interpolation : ",a)') trim(climat%tinterp); call goErr write (gol,'("in ",a)') rname; call goErr; status=1; return end select ! ok status = 0 end subroutine climat_Done ! *** subroutine climat_Set( climat, status, data, tr, t_in, data_in, it_in ) use GO, only : gol, goErr ! --- in/out -------------------------------- type(TClimat), intent(inout) :: climat integer, intent(inout) :: status real, intent(in), optional :: data(:,:,:) type(TDate), intent(in), optional :: tr(2) type(TDate), intent(in), optional :: t_in real, intent(in), optional :: data_in(:,:,:) integer, optional :: it_in ! --- const ------------------------------ character(len=*), parameter :: rname = mname//'/climat_Set' ! --- begin -------------------------------- ! trap previous errors: if (status>0) return ! store data ? if ( present(data) ) climat%data = data ! store time range within which data is valid: if ( present(tr) ) climat%tr = tr ! input fields for time interpolation ? if ( present(t_in) .or. present(data_in) .or. present(it_in) ) then ! all should be present ... if ( .not. ( present(t_in) .and. present(data_in) .and. present(it_in) ) ) then write (gol,'("all or none input arguments should be present:")'); call goErr write (gol,'(" t_in : ",l1)') present(t_in) write (gol,'(" data_in : ",l1)') present(data_in) write (gol,'(" it_in : ",l1)') present(it_in) write (gol,'("in ",a)') rname; call goErr; status=1; return end if ! check it_inex ... if ( (it_in < 1) .or. (it_in > climat%nt) ) then write (gol,'("it_in not valid:")'); call goErr write (gol,'(" it_in : ",i2)') it_in write (gol,'(" nt : ",i2)') climat%nt write (gol,'(" tinterp : ",a)') trim(climat%tinterp) write (gol,'("in ",a)') rname; call goErr; status=1; return end if ! store: climat%t_in(it_in) = t_in climat%data_in(:,:,:,it_in) = data_in end if ! ok status = 0 end subroutine climat_Set ! *** subroutine climat_Setup( climat, t, status ) use GO, only : gol, goErr use GO, only : TDate, IsAnyDate, operator(<), InterpolFractions ! --- in/out -------------------------------- type(TClimat), intent(inout) :: climat type(TDate), intent(in) :: t integer, intent(inout) :: status ! --- const ------------------------------ character(len=*), parameter :: rname = mname//'/climat_Setup' ! --- local ------------------------------- real :: alfa(2) integer :: it ! --- begin -------------------------------- ! trap previous errors: if (status>0) return ! not filled yet ? if ( IsAnyDate(climat%tr(1)) .or. IsAnyDate(climat%tr(2)) ) then ! return as a warning; calling program should fill data ... status=-1; return end if ! t outside time range for which current data is valid ? if ( (t < climat%tr(1)) .or. (climat%tr(2) < t) ) then ! return as a warning; calling program should fill data ... status=-1; return end if ! apply time interpolation: select case ( climat%tinterp ) case ( 'linear' ) ! determine fractions to be applied to data1 and data2: call InterpolFractions( t, climat%t_in(1), climat%t_in(2), alfa(1), alfa(2), status ) if (status/=0) then; write (gol,'("in ",a)') rname; call goErr; status=1; return; end if ! interpolate and store : climat%data = alfa(1) * climat%data_in(:,:,:,1) + alfa(2) * climat%data_in(:,:,:,2) case ( 'constant' ) ! nothing to be done case default write (gol,'("unsupported time interpolation : ",a)') trim(climat%tinterp); call goErr write (gol,'("in ",a)') rname; call goErr; status=1; return end select ! ok status = 0 end subroutine climat_Setup end module Climat