!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !! !! !! GNU General Public License !! !! !! !! This file is part of the Flexible Modeling System (FMS). !! !! !! !! FMS is free software; you can redistribute it and/or modify !! !! it and are expected to follow the terms of the GNU General Public !! !! License as published by the Free Software Foundation. !! !! !! !! FMS is distributed in the hope that it will be useful, !! !! but WITHOUT ANY WARRANTY; without even the implied warranty of !! !! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the !! !! GNU General Public License for more details. !! !! !! !! You should have received a copy of the GNU General Public License !! !! along with FMS; if not, write to: !! !! Free Software Foundation, Inc. !! !! 59 Temple Place, Suite 330 !! !! Boston, MA 02111-1307 USA !! !! or see: !! !! http://www.gnu.org/licenses/gpl.txt !! !! !! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !--------------------------------------------------------------------- !------------ FMS version number and tagname for this file ----------- ! $Id: cosp_radar.f90,v 1.1.2.1.4.1 2009/08/10 10:47:56 rsh Exp $ ! $Name: mom4p1_pubrel_dec2009_nnz $ ! (c) British Crown Copyright 2008, the Met Office. ! All rights reserved. ! ! Redistribution and use in source and binary forms, with or without modification, are permitted ! provided that the following conditions are met: ! ! * Redistributions of source code must retain the above copyright notice, this list ! of conditions and the following disclaimer. ! * Redistributions in binary form must reproduce the above copyright notice, this list ! of conditions and the following disclaimer in the documentation and/or other materials ! provided with the distribution. ! * Neither the name of the Met Office nor the names of its contributors may be used ! to endorse or promote products derived from this software without specific prior written ! permission. ! ! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR ! IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND ! FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR ! CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL ! DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, ! DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER ! IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT ! OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. MODULE MOD_COSP_RADAR USE MOD_COSP_CONSTANTS USE MOD_COSP_TYPES use radar_simulator_types use array_lib use atmos_lib use format_input IMPLICIT NONE INTERFACE subroutine radar_simulator(me,freq,k2,do_ray,use_gas_abs,use_mie_table,mt, & nhclass,hp,nprof,ngate,nsizes,D,hgt_matrix,hm_matrix,re_matrix,p_matrix,t_matrix, & rh_matrix,Ze_non,Ze_ray,h_atten_to_vol,g_atten_to_vol,dBZe, & g_to_vol_in,g_to_vol_out) use m_mrgrnk use array_lib use math_lib use optics_lib use radar_simulator_types implicit none ! ----- INPUTS ----- integer, intent(in) :: me type(mie), intent(in) :: mt type(class_param) :: hp real*8, intent(in) :: freq,k2 integer, intent(in) :: do_ray,use_gas_abs,use_mie_table, & nhclass,nprof,ngate,nsizes real*8, dimension(nsizes), intent(in) :: D real*8, dimension(nprof,ngate), intent(in) :: hgt_matrix, p_matrix, & t_matrix,rh_matrix real*8, dimension(nhclass,nprof,ngate), intent(in) :: hm_matrix real*8, dimension(nhclass,nprof,ngate), intent(inout) :: re_matrix ! ----- OUTPUTS ----- real*8, dimension(nprof,ngate), intent(out) :: Ze_non,Ze_ray, & g_atten_to_vol,dBZe,h_atten_to_vol ! ----- OPTIONAL ----- real*8, optional, dimension(ngate,nprof) :: & g_to_vol_in,g_to_vol_out end subroutine radar_simulator END INTERFACE CONTAINS !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% !------------------- SUBROUTINE COSP_RADAR ------------------------ !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% SUBROUTINE COSP_RADAR(me, gbx,sgx,sghydro,z) IMPLICIT NONE ! Arguments integer, intent(in) :: me type(cosp_gridbox),intent(in) :: gbx ! Gridbox info type(cosp_subgrid),intent(in) :: sgx ! Subgrid info type(cosp_sghydro),intent(in) :: sghydro ! Subgrid info for hydrometeors type(cosp_sgradar),intent(inout) :: z ! Output from simulator, subgrid ! Local variables integer :: & nsizes ! num of discrete drop sizes real*8 :: & freq, & ! radar frequency (GHz) k2 ! |K|^2, -1=use frequency dependent default real*8, dimension(:,:), allocatable :: & g_to_vol ! integrated atten due to gases, r>v (dB) real*8, dimension(:,:), allocatable :: & Ze_non, & ! radar reflectivity withOUT attenuation (dBZ) Ze_ray, & ! Rayleigh reflectivity (dBZ) h_atten_to_vol, & ! attenuation by hydromets, radar to vol (dB) g_atten_to_vol, & ! gaseous atteunation, radar to vol (dB) dBZe, & ! effective radar reflectivity factor (dBZ) hgt_matrix, & ! height of hydrometeors (km) t_matrix, & !temperature (k) p_matrix, & !pressure (hPa) rh_matrix !relative humidity (%) real*8, dimension(:,:,:), allocatable :: & hm_matrix, & ! hydrometeor mixing ratio (g kg^-1) re_matrix integer, parameter :: one = 1 logical :: hgt_reversed integer :: pr,i,j,k,unt ! ----- main program settings ------ freq = gbx%radar_freq k2 = gbx%k2 ! ! note: intitialization section that was here has been relocated to SUBROUTINE CONSTRUCT_COSP_GRIDBOX by roj, Feb 2008 ! mt_ttl=gbx%mt_ttl ! these variables really should be moved into the mt structure rather than kept as global arrays. mt_tti=gbx%mt_tti ! Inputs to Quickbeam allocate(hgt_matrix(gbx%Npoints,gbx%Nlevels),p_matrix(gbx%Npoints,gbx%Nlevels), & t_matrix(gbx%Npoints,gbx%Nlevels),rh_matrix(gbx%Npoints,gbx%Nlevels)) allocate(hm_matrix(gbx%Nhydro,gbx%Npoints,gbx%Nlevels)) allocate(re_matrix(gbx%Nhydro,gbx%Npoints,gbx%Nlevels)) ! Outputs from Quickbeam allocate(Ze_non(gbx%Npoints,gbx%Nlevels)) allocate(Ze_ray(gbx%Npoints,gbx%Nlevels)) allocate(h_atten_to_vol(gbx%Npoints,gbx%Nlevels)) allocate(g_atten_to_vol(gbx%Npoints,gbx%Nlevels)) allocate(dBZe(gbx%Npoints,gbx%Nlevels)) ! Optional argument. It is computed and returned in the first call to ! radar_simulator, and passed as input in the rest allocate(g_to_vol(gbx%Nlevels,gbx%Npoints)) p_matrix = gbx%p/100.0 ! From Pa to hPa hgt_matrix = gbx%zlev/1000.0 ! From m to km t_matrix = gbx%T-273.15 ! From K to C rh_matrix = gbx%q re_matrix = 0.0 ! Quickbeam assumes the first row is closest to the radar call order_data(hgt_matrix,hm_matrix,p_matrix,t_matrix, & rh_matrix,gbx%surface_radar,hgt_reversed) ! ----- loop over subcolumns ----- do pr=1,sgx%Ncolumns ! atmospheric profiles are the same within the same gridbox ! only hydrometeor profiles will be different if (hgt_reversed) then do i=1,gbx%Nhydro hm_matrix(i,:,:) = sghydro%mr_hydro(:,pr,gbx%Nlevels:1:-1,i)*1000.0 ! Units from kg/kg to g/kg if (gbx%use_reff) then re_matrix(i,:,:) = sghydro%Reff(:,pr,gbx%Nlevels:1:-1,i)*1.e6 ! Units from m to micron endif enddo else do i=1,gbx%Nhydro hm_matrix(i,:,:) = sghydro%mr_hydro(:,pr,:,i)*1000.0 ! Units from kg/kg to g/kg if (gbx%use_reff) then re_matrix(i,:,:) = sghydro%Reff(:,pr,:,i)*1.e6 ! Units from m to micron endif enddo endif ! ----- call radar simulator ----- if (pr == 1) then ! Compute gaseous attenuation for all profiles j=0 if (gbx%Npoints == 53) then unt=10 j=1 endif if (gbx%Npoints == 153) then unt=11 j=101 endif call radar_simulator(me,freq,k2,gbx%do_ray,gbx%use_gas_abs,gbx%use_mie_tables,gbx%mt, & ! v0.2: mt changed to gbx%mt, roj gbx%Nhydro,gbx%hp,gbx%Npoints,gbx%Nlevels,gbx%nsizes,gbx%D, & ! v0.2: hp->gbx%hp, D->gbx%d, nsizes->gbx%nsizes, roj hgt_matrix,hm_matrix,re_matrix,p_matrix,t_matrix,rh_matrix, & Ze_non,Ze_ray,h_atten_to_vol,g_atten_to_vol,dBZe,g_to_vol_out=g_to_vol) else ! Use gaseous atteunuation for pr = 1 call radar_simulator(me, freq,k2,gbx%do_ray,gbx%use_gas_abs,gbx%use_mie_tables,gbx%mt, & gbx%Nhydro,gbx%hp,gbx%Npoints,gbx%Nlevels,gbx%nsizes,gbx%D, & hgt_matrix,hm_matrix,re_matrix,p_matrix,t_matrix,rh_matrix, & Ze_non,Ze_ray,h_atten_to_vol,g_atten_to_vol,dBZe,g_to_vol_in=g_to_vol) endif ! ----- BEGIN output section ----- ! spaceborne radar : from TOA to SURFACE if (gbx%surface_radar == 1) then z%Ze_tot(:,pr,:)=dBZe(:,:) else if (gbx%surface_radar == 0) then ! Spaceborne z%Ze_tot(:,pr,:)=dBZe(:,gbx%Nlevels:1:-1) endif enddo !pr ! Change undefined value to one defined in COSP where (z%Ze_tot == -999.0) z%Ze_tot = R_UNDEF deallocate(hgt_matrix,p_matrix,t_matrix,rh_matrix) deallocate(hm_matrix,re_matrix, & Ze_non,Ze_ray,h_atten_to_vol,g_atten_to_vol,dBZe) deallocate(g_to_vol) ! deallocate(mt_ttl,mt_tti) !v0.2: roj feb 2008 can not be done here, !these variables now part of gbx structure and dealocated later END SUBROUTINE COSP_RADAR END MODULE MOD_COSP_RADAR