source: LMDZ6/trunk/libf/phylmdiso/cosp/radar_simulator_init.F90 @ 3927

Last change on this file since 3927 was 3927, checked in by Laurent Fairhead, 3 years ago

Initial import of the physics wih isotopes from Camille Risi
CR

File size: 4.8 KB
Line 
1  subroutine radar_simulator_init(freq,k2,use_gas_abs,do_ray,undef, &
2                  nhclass, &
3                  hclass_type,hclass_phase, &
4                      hclass_dmin,hclass_dmax, &
5                          hclass_apm,hclass_bpm,hclass_rho, &
6                          hclass_p1,hclass_p2,hclass_p3, &
7                          load_scale_LUTs_flag,update_scale_LUTs_flag,LUT_file_name, &
8                  hp &      ! output
9                  )
10  use radar_simulator_types
11  implicit none
12 
13! Purpose:
14!
15!   Initialize variables used by the radar simulator.
16!   Part of QuickBeam v3.0 by John Haynes and Roj Marchand
17!   
18!
19! Inputs: 
20!   []   from data in hydrometeor class input
21!
22!   [freq]            radar frequency (GHz)
23!
24!   [k2]              |K|^2, the dielectric constant, set to -1 to use the
25!                     frequency dependent default
26!
27!   [use_gas_abs]     1=do gaseous abs calcs, 0=no gasesous absorbtion calculated,
28!                     2=calculate (and use) absorption for first profile on all profiles
29!
30!   [undef]           mising data value
31!   [nhclass]         number of hydrometeor types
32!
33!   For each hydrometero type:
34!       hclass_type     Type of distribution (see quickbeam documentation)
35!       hclass_phase            1==ice, 0=liquid
36!
37!       hclass_dmin         minimum diameter allowed is drop size distribution N(D<Dmin)=0
38!       hclass_dmax         maximum diameter allowed is drop size distribution N(D>Dmax)=0
39!
40!   hclass_apm,hclass_bpm   Density of partical apm*D^bpm or constant = rho
41!   hclass_rho,
42!   
43!       hclass_p1,hclass_p2,    Default values of DSD parameters (see quickbeam documentation)
44!   hclass_p3,
45!
46!   load_scale_LUTs_flag    Flag, load scale factors Look Up Table from file at start of run
47!   update_scale_LUTs_flag  Flag, save new scale factors calculated during this run to LUT
48!   LUT_file_name       Name of file containing LUT
49!
50! Outputs:
51!   [hp]            structure that define hydrometeor types
52!
53! Modified:
54!   08/23/2006  placed into subroutine form (Roger Marchand)
55!   June 2010   New interface to support "radar_simulator_params" structure
56   
57! ----- INPUT -----
58
59   real, intent(in)    :: freq,k2
60   integer, intent(in) :: nhclass       ! number of hydrometeor classes in use
61   integer, intent(in) :: use_gas_abs,do_ray
62   real :: undef
63   real,dimension(nhclass),intent(in)     ::    hclass_dmin,hclass_dmax, &
64                                hclass_apm,hclass_bpm,hclass_rho, &
65                                hclass_p1,hclass_p2,hclass_p3
66   integer,dimension(nhclass),intent(in)  ::    hclass_type,hclass_phase
67 
68   logical, intent(in)       :: load_scale_LUTs_flag,update_scale_LUTs_flag
69   character*240, intent(in) :: LUT_file_name
70
71! ----- OUTPUTS ----- 
72  type(class_param), intent(out) :: hp
73
74! ----- INTERNAL ----- 
75  integer :: i,j
76  real*8  :: delt, deltp
77       
78    !
79    ! set radar simulation properites
80    !
81    hp%freq=freq
82    hp%k2=k2
83    hp%use_gas_abs=use_gas_abs
84    hp%do_ray=do_ray
85    hp%nhclass=nhclass
86   
87    hp%load_scale_LUTs=load_scale_LUTs_flag
88    hp%update_scale_LUTs=update_scale_LUTs_flag
89    hp%scale_LUT_file_name=LUT_file_name
90   
91    !
92        ! Store settings for hydrometeor types in hp (class_parameter) structure.   
93        !
94        do i = 1,nhclass
95        hp%dtype(i) = hclass_type(i)
96        hp%phase(i) = hclass_phase(i)
97        hp%dmin(i)  = hclass_dmin(i)
98        hp%dmax(i)  = hclass_dmax(i)
99        hp%apm(i)   = hclass_apm(i)
100        hp%bpm(i)   = hclass_bpm(i)
101        hp%rho(i)   = hclass_rho(i)
102        hp%p1(i)    = hclass_p1(i)
103        hp%p2(i)    = hclass_p2(i)
104        hp%p3(i)    = hclass_p3(i)
105        enddo
106   
107        !
108        ! initialize scaling array
109        !
110        hp%N_scale_flag = .false.
111        hp%fc = undef
112        hp%rho_eff = undef
113   
114        hp%Z_scale_flag = .false.
115        hp%Ze_scaled = 0.0
116        hp%Zr_scaled = 0.0
117        hp%kr_scaled = 0.0
118 
119        !
120    ! set up Re bin "structure" for z_scaling
121        !
122    hp%base_list(1)=0;
123    do j=1,Re_MAX_BIN
124        hp%step_list(j)=0.1+0.1*((j-1)**1.5);
125        if(hp%step_list(j)>Re_BIN_LENGTH) then
126            hp%step_list(j)=Re_BIN_LENGTH;
127        endif
128        if(j>1) then
129            hp%base_list(j)=hp%base_list(j-1)+floor(Re_BIN_LENGTH/hp%step_list(j-1));
130        endif
131    enddo
132
133    !
134    ! set up Temperature bin structure used for z scaling
135    !
136    do i=1,cnt_ice
137        hp%mt_tti(i)=(i-1)*5-90 + 273.15
138    enddo
139   
140    do i=1,cnt_liq
141        hp%mt_ttl(i)=(i-1)*5-60 + 273.15
142    enddo
143 
144    !
145        ! define array discrete diameters used in mie calculations
146        !
147        delt = (log(dmax)-log(dmin))/(nd-1)
148        deltp = exp(delt)
149
150        hp%D(1) = dmin
151        do i=2,nd
152          hp%D(i) = hp%D(i-1)*deltp
153        enddo   
154 
155 
156  end subroutine radar_simulator_init
Note: See TracBrowser for help on using the repository browser.