source: LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/mod_llnl_stats.F90 @ 5205

Last change on this file since 5205 was 5185, checked in by abarral, 2 months ago

Replace REPROBUS CPP KEY by logical using handmade wonky wrapper

  • Property copyright set to
    Name of program: LMDZ
    Creation date: 1984
    Version: LMDZ5
    License: CeCILL version 2
    Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539
    See the license file in the root directory
File size: 5.9 KB
Line 
1! (c) 2008, Lawrence Livermore National Security Limited Liability Corporation.
2! All rights reserved.
3! $Revision: 88 $, $Date: 2013-11-13 15:08:38 +0100 (mer. 13 nov. 2013) $
4! $URL: http://cfmip-obs-sim.googlecode.com/svn/stable/v1.4.0/llnl/llnl_stats.F90 $
5
6! Redistribution and use in source and binary forms, with or without modification, are permitted
7! provided that the following conditions are met:
8
9!     * Redistributions of source code must retain the above copyright notice, this list
10!       of conditions and the following disclaimer.
11!     * Redistributions in binary form must reproduce the above copyright notice, this list
12!       of conditions and the following disclaimer in the documentation and/or other materials
13!       provided with the distribution.
14!     * Neither the name of the Lawrence Livermore National Security Limited Liability Corporation
15!       nor the names of its contributors may be used to endorse or promote products derived from
16!       this software without specific prior written permission.
17
18! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR
19! IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
20! FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR
21! CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
22! DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
23! DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER
24! IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT
25! OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
26
27! History
28
29! Jan 2013 - G. Cesana        - Added betaperp_tot and temp_tot arguments
30
31MODULE MOD_LLNL_STATS
32  USE MOD_COSP_CONSTANTS
33  IMPLICIT NONE
34
35CONTAINS
36
37!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
38!-------------------- FUNCTION COSP_CFAD ------------------------
39!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
40FUNCTION COSP_CFAD(Npoints,Ncolumns,Nlevels,Nbins,x,xmin,xmax,bmin,bwidth)
41   ! Input arguments
42   integer,intent(in) :: Npoints,Ncolumns,Nlevels,Nbins
43   real,dimension(Npoints,Ncolumns,Nlevels),intent(in) :: x
44   real,intent(in) :: xmin,xmax
45   real,intent(in) :: bmin,bwidth
46   
47   real,dimension(Npoints,Nbins,Nlevels) :: cosp_cfad
48   ! Local variables
49   integer :: i, j, k
50   integer :: ibin
51   
52   !--- Input arguments
53   ! Npoints: Number of horizontal points
54   ! Ncolumns: Number of subcolumns
55   ! Nlevels: Number of levels
56   ! Nbins: Number of x axis bins
57   ! x: variable to process (Npoints,Ncolumns,Nlevels)
58   ! xmin: minimum value allowed for x
59   ! xmax: minimum value allowed for x
60   ! bmin: mimumum value of first bin
61   ! bwidth: bin width
62
63   ! Output: 2D histogram on each horizontal point (Npoints,Nbins,Nlevels)
64   
65   cosp_cfad = 0.0
66   ! bwidth intervals in the range [bmin,bmax=bmin+Nbins*hwidth]
67   ! Valid x values smaller than bmin and larger than bmax are set
68   ! into the smallest bin and largest bin, respectively.
69   DO j = 1, Nlevels, 1
70      DO k = 1, Ncolumns, 1
71         DO i = 1, Npoints, 1
72            if (x(i,k,j) == R_GROUND) then
73               cosp_cfad(i,:,j) = R_UNDEF
74            elseif ((x(i,k,j) >= xmin) .AND. (x(i,k,j) <= xmax)) then
75               ibin = ceiling((x(i,k,j) - bmin)/bwidth)
76               if (ibin > Nbins) ibin = Nbins
77               if (ibin < 1)     ibin = 1
78               cosp_cfad(i,ibin,j) = cosp_cfad(i,ibin,j) + 1.0
79            end if
80         enddo  !i
81      enddo  !k
82   enddo  !j
83   where ((cosp_cfad /= R_UNDEF).AND.(cosp_cfad /= 0.0)) cosp_cfad = cosp_cfad / Ncolumns
84END FUNCTION COSP_CFAD
85
86!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
87!------------- SUBROUTINE COSP_LIDAR_ONLY_CLOUD -----------------
88!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
89SUBROUTINE COSP_LIDAR_ONLY_CLOUD(Npoints,Ncolumns,Nlevels,temp_tot,beta_tot, &
90                   betaperp_tot,beta_mol,Ze_tot,lidar_only_freq_cloud,tcc)
91   ! Input arguments
92   integer,intent(in) :: Npoints,Ncolumns,Nlevels
93   real,dimension(Npoints,Nlevels),intent(in) :: beta_mol   ! Molecular backscatter
94   real,dimension(Npoints,Ncolumns,Nlevels),intent(in) :: beta_tot   ! Total backscattered signal
95   real,dimension(Npoints,Ncolumns,Nlevels),intent(in) :: temp_tot   ! Total backscattered signal
96   real,dimension(Npoints,Ncolumns,Nlevels),intent(in) :: betaperp_tot   ! perpendicular Total backscattered signal
97   real,dimension(Npoints,Ncolumns,Nlevels),intent(in) :: Ze_tot     ! Radar reflectivity
98   ! Output arguments
99   real,dimension(Npoints,Nlevels),intent(out) :: lidar_only_freq_cloud
100   real,dimension(Npoints),intent(out) :: tcc
101
102   ! local variables
103   REAL :: sc_ratio
104   REAL :: s_cld, s_att
105   parameter (S_cld = 5.0)
106   parameter (s_att = 0.01)
107   integer :: flag_sat !first saturated level encountered from top
108   integer :: flag_cld !cloudy column
109   integer :: pr,i,j
110
111   lidar_only_freq_cloud = 0.0
112   tcc = 0.0
113   DO pr=1,Npoints
114     DO i=1,Ncolumns
115       flag_sat = 0
116       flag_cld = 0
117       DO j=Nlevels,1,-1 !top->surf
118        sc_ratio = beta_tot(pr,i,j)/beta_mol(pr,j)
119        if ((sc_ratio .le. s_att) .AND. (flag_sat .eq. 0)) flag_sat = j
120        if (Ze_tot(pr,i,j) .lt. -30.) then  !radar can't detect cloud
121         if ( (sc_ratio .gt. s_cld) .or. (flag_sat .eq. j) ) then  !lidar sense cloud
122            lidar_only_freq_cloud(pr,j)=lidar_only_freq_cloud(pr,j)+1. !top->surf
123            flag_cld=1
124         endif
125        else  !radar sense cloud (z%Ze_tot(pr,i,j) .ge. -30.)
126           flag_cld=1
127        endif
128       enddo !levels
129       if (flag_cld .eq. 1) tcc(pr)=tcc(pr)+1.
130     enddo !columns
131   enddo !points
132   lidar_only_freq_cloud=lidar_only_freq_cloud/Ncolumns
133   tcc=tcc/Ncolumns
134
135END SUBROUTINE COSP_LIDAR_ONLY_CLOUD
136END MODULE MOD_LLNL_STATS
Note: See TracBrowser for help on using the repository browser.