! First include the set of model-wide compiler flags #include "tm5.inc" Module MTObsFcInfo ! ! This module defines the structure to store the output data of the NRT retrieval ! Henk Eskes, KNMI, Feb 2015, Aug 2015 ! implicit none private public :: TObsFcInfo, ObsFcInfoAllocate, ObsFcInfoDeallocate ! ---------------------------------------------------- ! structure containing observation, forecast and the observation operator ! ---------------------------------------------------- ! ! count : number of NO2 observations (ground pixels) ! copy of the number of observations provided as input ! flag : Quality flag ! ! --> retrieval results ! note: vertical, slant columns and errors have units: 10^15 molecules cm^-2 ! no2vcd,no2vcdsig : vertical column and error ! Note: this is not equal to the sum of trop and strat ! no2vcdsum : = no2vcdtrop + no2vcdstrat ! no2vcdsumsig : vertical column error (for no2vcdtrop + no2vcdstrat) ! no2vcdsigak : vertical column error, when kernels are used ! no2vcdtrop : tropospheric vertical column ! no2vcdtropsig : tropospheric vertical column error (including smoothing) ! no2vcdtropsigak : tropospheric vertical column error, ! when kernels are used (excluding smoothing) ! no2vcdstrat : stratospheric vertical column and error ! no2vcdstratsig : stratospheric vertical column and error ! no2slcstrat : (model) stratospheric slant column ! amf,amftrop,amfstrat : air-mass factors: total, tropospere, stratosphere (unit 1) ! avkernel : averaging kernel vector (unit 1) ! dimension ("lm","count") ! ghostcol : ghost column (model predicted part of column as obscured by clouds) ! ! --> model derived values ! icell,jcell : pixel longitude, latitude grid index ! psurf : model surface pressure at pixel, unit "hPa" ! no2vcdfc,no2vcdfctrop : model forecast vertical column (10^15 mol/cm^2) ! no2slcfc,no2slcfctrop : model slant column (10^15 mol/cm^2) ! no2slcstrat : model stratospheric slant column (10^15 mol/cm^2) ! levtropopause : index of highest model level in the troposphere ! modelTerrainHeight : orography at TM5 model resolution, unit "m" ! ! -- extra fields -- ! cloudradfraction : (optional) computed when not available as input, [0-1] ! cloudradfraction_modified : .true. when CRF is computed in the observation operator module ! amfgeo : geometrical airmass factor (used internally) ! amfclear : clear-sky AMF (cloud fraction neglected) ! ! The one-dimensional arrays all have dimension "count" ! The unit of the vertical and slant columns is "10^15 molecules / cm^2" ! The air mass factors have unit "1" ! type TObsFcInfo integer :: count ! flagging integer,dimension(:),allocatable :: flag ! vcd, amf real,dimension(:),allocatable :: no2vcd real,dimension(:),allocatable :: no2vcdsig, no2vcdsigak real,dimension(:),allocatable :: no2vcdsum, no2vcdsumsig real,dimension(:),allocatable :: no2vcdtrop, no2vcdtropsig real,dimension(:),allocatable :: no2vcdtropsigak real,dimension(:),allocatable :: no2vcdstrat, no2vcdstratsig real,dimension(:),allocatable :: amf, amftrop, amfstrat real,dimension(:,:),allocatable :: avkernel real,dimension(:),allocatable :: ghostcol ! model derived quantities integer,dimension(:),allocatable :: icell, jcell real,dimension(:),allocatable :: psurf real,dimension(:),allocatable :: no2vcdfc, no2vcdfctrop real,dimension(:),allocatable :: no2slcfc, no2slcfctrop real,dimension(:),allocatable :: no2slcstrat integer,dimension(:),allocatable :: levtropopause real,dimension(:),allocatable :: modelTerrainHeight ! extra data real,dimension(:),allocatable :: amfgeo real,dimension(:),allocatable :: amfclear real,dimension(:),allocatable :: cloudRadFraction logical :: cloudRadFraction_computed end type TObsFcInfo contains subroutine ObsFcInfoAllocate ( obsFcInfo, npix, nlev ) ! ! Allocate arrays in structure ! use pqf_module, only : PQF_E_GENERIC_EXCEPTION implicit none ! type(TObsFcInfo),intent(inout) :: obsFcInfo integer, intent(in) :: npix, nlev ! real,parameter :: nf_fill_double = 9.9692099683868690d+36 integer, parameter :: nf_fill_int = -2147483647 ! print *, 'Allocating ObsFcInfo, dimensions ',npix,nlev allocate ( obsFcInfo%icell(npix) ) obsFcInfo%icell(:) = nf_fill_int allocate ( obsFcInfo%jcell(npix) ) obsFcInfo%jcell(:) = nf_fill_int allocate ( obsFcInfo%no2vcd(npix) ) obsFcInfo%no2vcd(:) = nf_fill_double allocate ( obsFcInfo%no2vcdsum(npix) ) obsFcInfo%no2vcdsum(:) = nf_fill_double allocate ( obsFcInfo%no2vcdsumsig(npix) ) obsFcInfo%no2vcdsumsig(:) = nf_fill_double allocate ( obsFcInfo%no2vcdsig(npix) ) obsFcInfo%no2vcdsig(:) = nf_fill_double allocate ( obsFcInfo%no2vcdsigak(npix) ) obsFcInfo%no2vcdsigak(:) = nf_fill_double allocate ( obsFcInfo%no2vcdtrop(npix) ) obsFcInfo%no2vcdtrop(:) = nf_fill_double allocate ( obsFcInfo%no2vcdtropsig(npix) ) obsFcInfo%no2vcdtropsig(:) = nf_fill_double allocate ( obsFcInfo%no2vcdtropsigak(npix) ) obsFcInfo%no2vcdtropsigak(:) = nf_fill_double allocate ( obsFcInfo%no2vcdstrat(npix) ) obsFcInfo%no2vcdstrat(:) = nf_fill_double allocate ( obsFcInfo%no2vcdstratsig(npix) ) obsFcInfo%no2vcdstratsig(:) = nf_fill_double allocate ( obsFcInfo%no2slcstrat(npix) ) obsFcInfo%no2vcdstrat(:) = nf_fill_double allocate ( obsFcInfo%amf(npix) ) obsFcInfo%amf(:) = nf_fill_double allocate ( obsFcInfo%amftrop(npix) ) obsFcInfo%amftrop(:) = nf_fill_double allocate ( obsFcInfo%amfstrat(npix) ) obsFcInfo%amfstrat(:) = nf_fill_double allocate ( obsFcInfo%psurf(npix) ) obsFcInfo%psurf(:) = nf_fill_double allocate ( obsFcInfo%no2slcfc(npix) ) obsFcInfo%no2slcfc(:) = nf_fill_double allocate ( obsFcInfo%no2slcfctrop(npix) ) obsFcInfo%no2slcfctrop(:) = nf_fill_double allocate ( obsFcInfo%no2vcdfc(npix) ) obsFcInfo%no2vcdfc(:) = nf_fill_double allocate ( obsFcInfo%no2vcdfctrop(npix) ) obsFcInfo%no2vcdfctrop(:) = nf_fill_double allocate ( obsFcInfo%flag(npix) ) obsFcInfo%flag(:) = PQF_E_GENERIC_EXCEPTION allocate ( obsFcInfo%levtropopause(npix) ) obsFcInfo%levtropopause(:) = nf_fill_int allocate ( obsFcInfo%avkernel(nlev,npix) ) obsFcInfo%avkernel(:,:) = nf_fill_double allocate ( obsFcInfo%ghostcol(npix) ) obsFcInfo%ghostcol(:) = nf_fill_double allocate ( obsFcInfo%modelterrainheight(npix) ) obsFcInfo%modelterrainheight(:) = nf_fill_double ! extra allocate ( obsFcInfo%cloudradfraction(npix) ) obsFcInfo%cloudradfraction(:) = nf_fill_double allocate ( obsFcInfo%amfgeo(npix) ) obsFcInfo%amfgeo(:) = nf_fill_double allocate ( obsFcInfo%amfclear(npix) ) obsFcInfo%amfclear(:) = nf_fill_double ! obsFcInfo%count = npix ! end subroutine ObsFcInfoAllocate subroutine ObsFcInfoDeallocate ( obsFcInfo ) ! ! De-allocate arrays in structure ! implicit none ! type(TObsFcInfo),intent(inout) :: obsFcInfo ! if ( allocated( obsFcInfo%icell ) ) deallocate ( obsFcInfo%icell ) if ( allocated( obsFcInfo%jcell ) ) deallocate ( obsFcInfo%jcell ) if ( allocated( obsFcInfo%no2vcd ) ) deallocate ( obsFcInfo%no2vcd ) if ( allocated( obsFcInfo%no2vcdsum ) ) deallocate ( obsFcInfo%no2vcdsum ) if ( allocated( obsFcInfo%no2vcdsumsig ) ) deallocate ( obsFcInfo%no2vcdsumsig ) if ( allocated( obsFcInfo%no2vcdsig ) ) deallocate ( obsFcInfo%no2vcdsig ) if ( allocated( obsFcInfo%no2vcdsigak ) ) deallocate ( obsFcInfo%no2vcdsigak ) if ( allocated( obsFcInfo%no2vcdtrop ) ) deallocate ( obsFcInfo%no2vcdtrop ) if ( allocated( obsFcInfo%no2vcdtropsig ) ) deallocate ( obsFcInfo%no2vcdtropsig ) if ( allocated( obsFcInfo%no2vcdtropsigak ) ) deallocate ( obsFcInfo%no2vcdtropsigak ) if ( allocated( obsFcInfo%no2vcdstrat ) ) deallocate ( obsFcInfo%no2vcdstrat ) if ( allocated( obsFcInfo%no2vcdstratsig ) ) deallocate ( obsFcInfo%no2vcdstratsig ) if ( allocated( obsFcInfo%no2slcstrat ) ) deallocate ( obsFcInfo%no2slcstrat ) if ( allocated( obsFcInfo%amf ) ) deallocate ( obsFcInfo%amf ) if ( allocated( obsFcInfo%amftrop ) ) deallocate ( obsFcInfo%amftrop ) if ( allocated( obsFcInfo%amfstrat ) ) deallocate ( obsFcInfo%amfstrat ) if ( allocated( obsFcInfo%psurf ) ) deallocate ( obsFcInfo%psurf ) if ( allocated( obsFcInfo%no2slcfc ) ) deallocate ( obsFcInfo%no2slcfc ) if ( allocated( obsFcInfo%no2slcfctrop ) ) deallocate ( obsFcInfo%no2slcfctrop ) if ( allocated( obsFcInfo%no2vcdfc ) ) deallocate ( obsFcInfo%no2vcdfc ) if ( allocated( obsFcInfo%no2vcdfctrop ) ) deallocate ( obsFcInfo%no2vcdfctrop ) if ( allocated( obsFcInfo%flag ) ) deallocate ( obsFcInfo%flag ) if ( allocated( obsFcInfo%levtropopause ) ) deallocate ( obsFcInfo%levtropopause ) if ( allocated( obsFcInfo%avkernel ) ) deallocate ( obsFcInfo%avkernel ) if ( allocated( obsFcInfo%ghostcol ) ) deallocate ( obsFcInfo%ghostcol ) if ( allocated( obsFcInfo%modelterrainheight ) ) deallocate ( obsFcInfo%modelterrainheight ) ! extra if ( allocated( obsFcInfo%cloudradfraction ) ) deallocate ( obsFcInfo%cloudradfraction ) if ( allocated( obsFcInfo%amfgeo ) ) deallocate ( obsFcInfo%amfgeo ) if ( allocated( obsFcInfo%amfclear ) ) deallocate ( obsFcInfo%amfclear ) ! obsFcInfo%count = -1 ! end subroutine ObsFcInfoDeallocate end Module MTObsFcInfo