diff --git a/src/MNH/gps_zenith_grid.f90 b/src/MNH/gps_zenith_grid.f90 new file mode 100644 index 0000000000000000000000000000000000000000..346ae30a7579f45d6320e3f766071231b5e86b01 --- /dev/null +++ b/src/MNH/gps_zenith_grid.f90 @@ -0,0 +1,250 @@ +!MNH_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier +!MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence +!MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt +!MNH_LIC for details. version 1. +!----------------------------------------------------------------- +!--------------- special set of characters for RCS information +!----------------------------------------------------------------- +! $Source$ $Revision$ $Date$ +!----------------------------------------------------------------- +!########################################## +MODULE MODI_GPS_ZENITH_GRID +!########################################## +INTERFACE + SUBROUTINE GPS_ZENITH_GRID(PRT,PTEMP,PPABST,PZTD,PZHD,PZWD) +! +!* 0.1 Declarations of dummy arguments : +! +REAL, DIMENSION(:,:,:), INTENT(IN) :: PRT ! Water vapor mixing ratio +REAL, DIMENSION(:,:,:), INTENT(IN) :: PTEMP ! Air temperature +REAL, DIMENSION(:,:,:), INTENT(IN) :: PPABST ! Absolute pressure +! +REAL, DIMENSION(:,:), INTENT(OUT) :: PZTD ! Zenithal Total Delay +REAL, DIMENSION(:,:), INTENT(OUT) :: PZHD ! Zenithal Hydrostatic Delay +REAL, DIMENSION(:,:), INTENT(OUT) :: PZWD ! Zenithal Wet Delay +! +END SUBROUTINE GPS_ZENITH_GRID +END INTERFACE +END MODULE MODI_GPS_ZENITH_GRID +! +! ################################################################### + SUBROUTINE GPS_ZENITH_GRID(PRT,PTEMP,PPABST,PZTD,PZHD,PZWD) +! ################################################################### +! +!!**** *GPS_ZENITH_GRID * - computes GPS zenithal delays parameters +!! +!! PURPOSE +!! ------- +!! The purpose of this routine is to compute the synthetic ZTD, ZHD, ZWD +!! +!!** METHOD +!! ------ +!! The delays are computed using the refractivity formula and its different +!! approximation. The delay can be separated in two parts. +!! The ZTD is the sum of the ZHD, including all component of the atmosphere, and +!! the ZWD, including the contribution of the non-dipolar moment of the molecule +!! of water vapor. +!! +!! EXTERNAL +!! -------- +!! NONE +!! +!! IMPLICIT ARGUMENTS +!! ------------------ +!! +!! REFERENCE +!! --------- +!! +!! Brenot et al, 2006, JGR +!! +!! AUTHORS +!! ------- +!! H. Brenot * Laboratoire de Geophysique Interne et Tectonophysique * +!! / * Meteo France * +!! & V. Ducrocq * Meteo France * +!! & A. Walpersdorf * Laboratoire de Geophysique Interne et Tectonophysique * +!! +!! MODIFICATIONS +!! ------------- +!! Original 18/11/04 +!! Modified 4/12/2007 +!! July, 2015 (O.Nuissier/F.Duffourg) Add microphysics diagnostic for +!! aircraft, ballon and profiler +!! +!------------------------------------------------------------------------------- +! +!* 0. DECLARATIONS +! ------------ +! +! +USE MODD_PARAMETERS +USE MODD_LUNIT +USE MODD_CST +USE MODE_FM +USE MODD_GR_FIELD_n +USE MODD_DIAG_FLAG +USE MODD_GRID, ONLY: XLONORI,XLATORI +USE MODD_GRID_n +USE MODE_GRIDPROJ +USE MODE_ll +! +USE MODI_IO_ll +! +IMPLICIT NONE +! +! +!* 0.1 Declarations of dummy arguments : +! +REAL, DIMENSION(:,:,:), INTENT(IN) :: PRT ! Water vapour mixing ratio +REAL, DIMENSION(:,:,:), INTENT(IN) :: PTEMP ! Air temperature +REAL, DIMENSION(:,:,:), INTENT(IN) :: PPABST ! Absolute pressure +! +REAL, DIMENSION(:,:), INTENT(OUT) :: PZTD ! Zenithal Total Delay +REAL, DIMENSION(:,:), INTENT(OUT) :: PZHD ! Zenithal Hydrostatic Delay +REAL, DIMENSION(:,:), INTENT(OUT) :: PZWD ! Zenithal Wet Delay +! +!* 0.2 Declarations of local variables : +! +!-------- Calculation parameters -------- +REAL :: ZK1,ZK2,ZK3 ! k1, k2 and K3 atmospheric refractivity constants +REAL :: ZRDSRV ! XRD/XRV +!-------- Loop parameters --------------- +INTEGER :: IIB,IIE ! Loop limits for coordinate X +INTEGER :: IJB,IJE ! Loop limits for coordinate Y +INTEGER :: IKB,IKE ! Loop limits for coordinate Z +INTEGER :: JK ! Loop variables of control +INTEGER :: IIU,IJU,IKU ! Loop variables of model +INTEGER :: ILUOUT0, IRESP ! file unit and return code for output +REAL, DIMENSION(:),ALLOCATABLE :: ZXHATM,ZYHATM ! mass-point positions +REAL, DIMENSION(:,:,:),ALLOCATABLE :: ZZHATM ! mass level altitude +!-------- Physical parameters for the integration ---------------------------- +REAL, DIMENSION(:,:,:),ALLOCATABLE :: ZE ! Partial pressure of water vapor +REAL, DIMENSION(:,:,:),ALLOCATABLE :: ZTV ! Virtual temperature +!-------- External model contributions --------------------------------------- +REAL, DIMENSION(:,:),ALLOCATABLE :: ZZHDX ! Zenith hydrostatic delay external (of the model contribution) +!-------- Top model parameters ----------------------------------------------- +REAL, DIMENSION(:,:), ALLOCATABLE :: ZPTOP ! pressure of top of model +REAL, DIMENSION(:,:), ALLOCATABLE :: ZTVTOP ! Tv of top of model +REAL, DIMENSION(:,:), ALLOCATABLE :: ZETOP ! water vapor partial pressure of top of model +REAL, DIMENSION(:,:), ALLOCATABLE :: ZTEMPTOP ! absolute temperature of top of model +REAL, DIMENSION(:,:),ALLOCATABLE :: ZG0 ! gravity of ground (Vedel approx.) +REAL, DIMENSION(:,:),ALLOCATABLE :: ZR0 ! radius of ground (Vedel approx.) +REAL, DIMENSION(:,:),ALLOCATABLE :: ZGTOP ! gravity of top of model using Vedel approx. +! +!------------------------------------------------------------------------------- +! +!* 1. INTIALIZE DIMENSIONS AND ALLOCATE ARRAYS +! ---------------------------------------- +! +CALL FMLOOK_ll(CLUOUT0,CLUOUT0,ILUOUT0,IRESP) +CALL GET_INDICE_ll(IIB,IJB,IIE,IJE) +IIU = SIZE (PTEMP,1) +IJU = SIZE (PTEMP,2) +IKU = SIZE (PTEMP,3) +IKB = JPVEXT + 1 +IKE = IKU - JPVEXT +! +ALLOCATE(ZXHATM(IIU)) +ALLOCATE(ZYHATM(IJU)) +ALLOCATE(ZZHATM(IIU,IJU,IKU)) +ALLOCATE(ZE(IIU,IJU,IKU)) +ALLOCATE(ZTV(IIU,IJU,IKU)) +ALLOCATE(ZGTOP(IIU,IJU)) +ALLOCATE(ZPTOP(IIU,IJU)) +ALLOCATE(ZTVTOP(IIU,IJU)) +ALLOCATE(ZTEMPTOP(IIU,IJU)) +! +ALLOCATE(ZZHDX(IIU,IJU)) +! +ALLOCATE(ZG0(IIU,IJU)) +ALLOCATE(ZR0(IIU,IJU)) +! +PZTD(:,:) = 0. ! Zenithal Total Delay +PZHD(:,:) = 0. ! Zenithal Hydrostatic Delay +PZWD(:,:) = 0. ! Zenithal Wet Delay +! +ZZHDX(:,:) = 0. +!------------------------------------------------------------------------------- +! +!* 3. REFRACTIVITY COEFFICIENTS AND OTHER CONSTANTS +! --------------------------------------------- +! +! Refractivity coeficients +! Bevis et al. (1994) +ZK1 = 0.776 ! K/Pa +ZK2 = 0.704 ! K/Pa +ZK3 = 3739. ! K2/Pa +ZRDSRV=XRD/XRV +! +!-------------------------------------------------------------------------------! +!* 4. AUXILLARY VARIABLES +! ------------------- +! +! +ZXHATM(1:IIU-1) = 0.5*(XXHAT(1:IIU-1)+XXHAT(2:IIU)) +ZXHATM(IIU) = 2.*XXHAT(IIU)-ZXHATM(IIU-1) +ZYHATM(1:IJU-1) = 0.5*(XYHAT(1:IJU-1)+XYHAT(2:IJU)) +ZYHATM(IJU) = 2.*XXHAT(IJU)-ZXHATM(IJU-1) +ZZHATM(:,:,1:IKU-1)=0.5*(XZZ(:,:,1:IKU-1)+XZZ(:,:,2:IKU)) +ZZHATM(:,:,IKU) = 2.*XZZ(:,:,IKU) -ZZHATM(:,:,IKU-1) +! +! Vapor pressure in Pa : ZE +ZE(:,:,:) = PPABST(:,:,:) * PRT(:,:,:) / ( ZRDSRV + PRT(:,:,:) ) +! Virtual temperature +ZTV(:,:,:) = PTEMP(:,:,:) * ( 1. + PRT(:,:,:) / ZRDSRV ) & + / ( 1. + PRT(:,:,:) ) +ZTVTOP(:,:)= ( ZTV(:,:,IKE) - ZTV(:,:,IKE-1) ) * ( XZZ(:,:,IKE+1) - ZZHATM(:,:,IKE) ) & + /(ZZHATM(:,:,IKE) - ZZHATM(:,:,IKE-1)) + ZTV(:,:,IKE) +! for extrapolation above model top : gtop, Ptop et Ttop +ZG0(:,:) = 9.780356 * ( 1. + 5.2885E-3 * (SIN(XLAT(:,:)*XPI/180.))**2 - & + 5.9E-6 * ( SIN(2.*XLAT(:,:)*XPI/180.))**2 ) +ZR0(:,:) = 6378.1E3 / SQRT ( ( SIN(XLAT(:,:)*XPI/180.)*6378.1/6356.6 )**2 + & + ( COS(XLAT(:,:)*XPI/180.))**2 ) +ZGTOP(:,:) = ZG0(:,:) * ( ZR0(:,:) / ( ZR0(:,:) + XZZ(:,:,IKE+1) ) )**2 +ZPTOP(:,:) = PPABST(:,:,IKE)*EXP(XG*0.5*(XZZ(:,:,IKE)-XZZ(:,:,IKE+1))/(XRD*ZTVTOP(:,:))) +ZTEMPTOP(:,:) = (( PTEMP(:,:,IKE) - PTEMP(:,:,IKE-1) ) * 0.5 * ( XZZ(:,:,IKE+1) - XZZ(:,:,IKE) ) & + / ( XZZ(:,:,IKE) - XZZ(:,:,IKE-1) )) + PTEMP(:,:,IKE) +! +!-------------------------------------------------------------------------------! +!* 5. ZENITH DELAYS +! ------------- +! +! 5.1 Integrated Model Contribution +! +! 5.1.1 ZHD +! +! level contribution for integrated formulation +DO JK=IKB,IKE + PZHD(:,:) = PZHD(:,:) + ( 1.E-6 * ZK1 * PPABST(:,:,JK) * & + ( XZZ(:,:,JK+1) - XZZ(:,:,JK) ) / ZTV(:,:,JK) ) +END DO +! +! 5.1.2 ZWD +! +DO JK=IKB,IKE + PZWD(:,:) = PZWD(:,:) + (1.E-6 * ( (ZK2-ZRDSRV*ZK1) + ( ZK3/PTEMP(:,:,JK) ) ) * & + ZE(:,:,JK)* ( XZZ(:,:,JK+1) - XZZ(:,:,JK) ) / PTEMP(:,:,JK)) +END DO +! +! 5.2 External Model Contribution (hydrostatic delay only, wet delay negligeable) +! +ZZHDX(:,:) = 1.E-6 * ZK1 * ZPTOP(:,:) * XRD * ( 1. + 2. * XRD * ZTEMPTOP(:,:) & + / ( ( XRADIUS + XZZ(:,:,IKE+1) ) * ZGTOP(:,:) ) + 2. * ( XRD * ZTEMPTOP(:,:) & + / ( (XRADIUS + XZZ(:,:,IKE+1)) * ZGTOP(:,:) ))**2 ) / ZGTOP(:,:) + +PZHD(:,:) = PZHD(:,:) + ZZHDX(:,:) +! +! 5.3 Total Zenith delay +! +PZTD(:,:) = PZHD(:,:) + PZWD(:,:) +! +! +DEALLOCATE(ZE) +DEALLOCATE(ZTV) +DEALLOCATE(ZGTOP) +DEALLOCATE(ZPTOP) +DEALLOCATE(ZTVTOP) +DEALLOCATE(ZTEMPTOP) +DEALLOCATE(ZZHDX) +! +END SUBROUTINE GPS_ZENITH_GRID