#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" ! !---------------------------------------------------------------------- ! READ_CAR reads the 7 coefficients for the ozone photochemistry ! linear parametrisation for each month and put it on a ! grid (regular along the latitudes) ! ! Parametrisation of the ozone photochemistry ! ! d r ! --- = < P - L > + ! dt ! d (P - L) ! < --------- > ( r - < r > ) + ! d r ! ! d (P - L) ! < --------- > ( T - < T > ) + ! d T ! ! d (P - L) ! < --------- > ( S - < S > ) ! d S ! ! with r : volumic ozone mixing ratio ! T : temperature ! S : ozone column above a point ! P : photochemical production of ozone ! L : photochemical destruction of ozone ! ! Terms between < > denote two-dimensional variables at photochemical ! equilibrium (from the 2D photochemical model). Other terms are ! three-dimensional variables (from the General Circulation Model) ! ! ! coefficient 1 : < P - L > units : s-1 ! ! d (P - L) ! coefficient 2 : < --------- > units : kg kg-1 s-1 ! d r ! ! coefficient 3 : < r > units : kg kg-1 ! ! d (P - L) ! coefficient 4 : < --------- > units : kg kg-1 s-1 K-1 ! d T ! ! coefficient 5 : < T > units : K ! ! d (P - L) ! coefficient 6 : < --------- > units : kg cm2 kg-1 s-1 molecule-1 ! d S ! ! coefficient 7 : < S > units : molecule cm-2 ! ! output : PHI (64) array containing latitudes of the grid ! PRES (64x26/45x12) array containing pressures of the grid ! COEFF(64x26/45x7x12) array containing the coefficients on the ! grid !---------------------------------------------------------------------- module file_cariolle implicit none ! --- in/out --------------------------------- private public :: Cariolle_Info public :: Cariolle_Read ! --- const ----------------------------- character(len=*), parameter :: mname = 'file_cariolle' ! dimensions integer, parameter :: nlat = 64 integer, parameter :: nlev = 45 ! 45 = v2.1 ! 26 = v1.0 integer, parameter :: nc = 7 integer, parameter :: nt = 12 contains ! ===================================================================== subroutine Cariolle_Info( status, nlats, nlevs, ncs ) ! --- in/out ------------------------- integer, intent(inout) :: status integer, intent(out), optional :: nlats integer, intent(out), optional :: nlevs integer, intent(out), optional :: ncs ! --- const --------------------------------- character(len=*), parameter :: rname = mname//'/Cariolle_Info' ! --- begin --------------------------- ! trap previous errors: if (status>0) return ! return info: if ( present(nlats) ) nlats = nlat if ( present(nlevs) ) nlevs = nlev if ( present(ncs ) ) ncs = nc ! ok status = 0 end subroutine Cariolle_Info ! *** subroutine Cariolle_Read( fname, imonth, phi, pres, coeff, status ) use GO , only : gol, goErr use Binas, only : Dobs ! --- in/out --------------------------------- character(len=*), intent(in) :: fname integer, intent(in) :: imonth real, intent(out) :: phi(nlat) ! deg real, intent(out) :: pres(nlat,nlev) ! Pa , top->surf real, intent(out) :: coeff(nlat,nlev,nc) integer, intent(inout) :: status ! --- const --------------------------------- character(len=*), parameter :: rname = mname//'/Cariolle_Read' ! --- local ----------------------------------- integer :: fu logical :: opened integer :: k integer :: ic, it ! --- begin ----------------------------------- ! trap previous errors: if (status>0) return ! select free file unit: fu = 123 do inquire( unit=fu, opened=opened ) if ( .not. opened ) exit fu = fu + 1 end do ! ! open file: ! open( fu, file=trim(fname), form='formatted', & status='old', iostat=status ) if ( status /= 0 ) then write (gol,'("from open Cariolle coeff file:")'); call goErr write (gol,'(" ",a)') trim( fname); call goErr write (gol,'("in ",a)') rname; call goErr; status=1; return end if ! ! read the header ! ! 98 ==> CARIOLLE V1.0 ! 122 ==> CARIOLLE V2.1 ! do k = 1, 122 read (fu,*,iostat=status) if ( status /= 0 ) then write (gol,'("reading from Cariolle coeff file:")'); call goErr write (gol,'(" file : ",a)') trim( fname); call goErr write (gol,'(" line : ",i4)') k; call goErr write (gol,'("in ",a)') rname; call goErr; status=1; return end if end do ! ! read the latitudes ! read (fu,*) read (fu,'(10f7.2)',iostat=status) phi if ( status /= 0 ) then write (gol,'("reading latitudes Cariolle coeff file:")'); call goErr write (gol,'(" file : ",a)') trim( fname); call goErr write (gol,'("in ",a)') rname; call goErr; status=1; return end if ! ! read the data (pressure + 7 coefficients) ! ! loop for each month of the year do it = 1, imonth ! read (fu,*) read (fu,*) read (fu,'(10f9.2)',iostat=status) pres ! hPa if ( status /= 0 ) then write (gol,'("reading pressures from Cariolle coeff file:")'); call goErr write (gol,'(" file : ", a)') trim( fname); call goErr write (gol,'(" month : ",i2)') it ; call goErr write (gol,'("in ",a)') rname; call goErr; status=1; return end if ! do ic = 1, nc ! ! loop on the 7 coefficients read (fu,*) read (fu,'(10e13.5)',iostat=status) coeff(:,:,ic) if ( status /= 0 ) then write (gol,'("reading Cariolle coeffs:")'); call goErr write (gol,'(" file : ",a )') trim( fname); call goErr write (gol,'(" month : ",i2)') it ; call goErr write (gol,'(" coef : ",i2)') ic ; call goErr write (gol,'("in ",a)') rname; call goErr; status=1; return end if end do ! end do ! ! close file ! close( fu, iostat=status ) if ( status /= 0 ) then write (gol,'("closing Cariolle coeff file:")'); call goErr write (gol,'(" file : ",a)') trim( fname); call goErr write (gol,'("in ",a)') rname; call goErr; status=1; return end if ! ! unit conversion ! ! pressures from hPa to Pa : pres = pres * 1.0e2 ! ! done ! ! ok status = 0 end subroutine Cariolle_Read end module file_cariolle