Ignore:
Timestamp:
Jul 22, 2024, 9:29:09 PM (4 months ago)
Author:
abarral
Message:

Replace most uses of CPP_DUST by the corresponding logical defined in lmdz_cppkeys_wrapper.F90
Convert several files from .F to .f90 to allow Dust to compile w/o rrtm/ecrad
Create lmdz_yoerad.f90
(lint) Remove "!" on otherwise empty line

Location:
LMDZ6/branches/Amaury_dev/libf/phylmd/cosp
Files:
36 edited

Legend:

Unmodified
Added
Removed
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/MISR_simulator.F

    r5095 r5099  
    1 !
     1
    22! Copyright (c) 2009,  Roger Marchand, version 1.2
    33! All rights reserved.
    44! $Revision: 23 $, $Date: 2011-03-31 15:41:37 +0200 (jeu. 31 mars 2011) $
    55! $URL: http://cfmip-obs-sim.googlecode.com/svn/stable/v1.4.0/MISR_simulator/MISR_simulator.f $
    6 !
     6
    77! Redistribution and use in source and binary forms, with or without modification, are permitted
    88! provided that the following conditions are met:
    9 !
     9
    1010!     * Redistributions of source code must retain the above copyright notice, this list of
    1111!       conditions and the following disclaimer.
     
    1515!     * Neither the name of the University of Washington nor the names of its contributors may be used
    1616!       to endorse or promote products derived from this software without specific prior written permission.
    17 !
     17
    1818! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING,
    1919! BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT
     
    2222! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
    2323! NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
    24 !
    2524
    2625      SUBROUTINE MISR_simulator(
     
    119118   
    120119      tauchk = -1.*log(0.9999999)
    121        
    122       !
     120
    123121      ! For each GCM cell or horizontal model grid point ...
    124       !
     122
    125123      do j=1,npoints   
    126124
    127          !
    128125         !  estimate distribution of Model layer tops
    129          ! 
     126
    130127         dist_model_layertops(j,:)=0
    131128
     
    154151     &          dist_model_layertops(j,iMISR_ztop)+1
    155152       enddo
    156    
    157    
    158          !
     153
    159154         ! compute total cloud optical depth for each column
    160          !       
     155
    161156       do ibox=1,ncol     
    162157       
     
    280275       enddo  ! loop of subcolumns
    281276       enddo    ! loop of gridpoints
    282        
    283 
    284         !     
     277
    285278        !   Modify MISR CTH for satellite spatial / pattern matcher effects
    286     !
     279
    287280    !   Code in this region added by roj 5/2006 to account
    288281    !   for spatial effect of the MISR pattern matcher.
     
    291284    !   a lower CTH, THEN misr will tend to but place the
    292285    !   odd column at the same height as it neighbors.
    293     !
     286
    294287    !   This setup assumes the columns represent a about a 1 to 4 km scale
    295288    !   it will need to be modified significantly, otherwise
     
    337330        endif
    338331
    339         !     
    340332    !     DETERMINE CLOUD TYPE FREQUENCIES
    341     !
     333
    342334    !     Now that ztop and tau have been determined,
    343335    !     determine amount of each cloud type
     
    407399
    408400          else
    409            
    410             !
     401
    411402            ! determine index for MISR bin set
    412             !
    413403
    414404            iMISR_ztop=2
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/array_lib.F90

    r5095 r5099  
    44! Compiled/Modified:
    55!   07/01/06  John Haynes (haynes@atmos.colostate.edu)
    6 !
     6
    77! infind (function)
    88! lin_interpolate (function)
     
    1919  use m_mrgrnk
    2020  implicit none
    21 !
     21
    2222! Purpose:
    2323!   Finds the index of an array that is closest to a value, plus the
    2424!   difference between the value found and the value specified
    25 !
     25
    2626! Inputs:
    2727!   [list]   an array of sequential values
     
    2929! Optional input:
    3030!   [sort]   set to 1 if [list] is in unknown/non-sequential order
    31 !
     31
    3232! Returns:
    3333!   index of [list] that is closest to [val]
    34 !
     34
    3535! Optional output:
    3636!   [dist]   set to variable containing [list([result])] - [val]
    37 !
     37
    3838! Requires:
    3939!   mrgrnk library
    40 !
     40
    4141! Created:
    4242!   10/16/03  John Haynes (haynes@atmos.colostate.edu)
     
    101101  use m_mrgrnk
    102102  implicit none
    103 !
     103
    104104! Purpose:
    105105!   linearly interpolate a set of y2 values given a set of y1,x1,x2
    106 !
     106
    107107! Inputs:
    108108!   [yarr]    an array of y1 values
     
    110110!   [xxarr]   an array of x2 values
    111111!   [tol]     maximum distance for a match
    112 !
     112
    113113! Output:
    114114!   [yyarr]   interpolated array of y2 values
    115 !
     115
    116116! Requires:
    117117!   mrgrnk library
    118 !
     118
    119119! Created:
    120120!   06/07/06  John Haynes (haynes@atmos.colostate.edu)
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/atmos_lib.F90

    r5095 r5099  
    44! Compiled/Modified:
    55!   07/01/06  John Haynes (haynes@atmos.colostate.edu)
    6 !
     6
    77! mcclatchey (subroutine)
    88 
     
    1717  subroutine mcclatchey(stype,hgt,prs,tk,rh)
    1818  implicit none
    19 !
     19
    2020! Purpose:
    2121!   returns a standard atmospheric profile
    22 !
     22
    2323! Input:
    2424!   [stype]   type of profile to return
     
    2626!             2 = mid-latitude winter
    2727!             3 = tropical
    28 !
     28
    2929! Outputs:
    3030!   [hgt]     height (m)
     
    3232!   [tk]      temperature (K)
    3333!   [rh]      relative humidity (%)
    34 !
     34
    3535! Created:
    3636!   06/01/2006  John Haynes (haynes@atmos.colostate.edu)
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/calc_Re.F90

    r5095 r5099  
    77! Purpose:
    88!   Calculates Effective Radius (1/2 distribution 3rd moment / 2nd moment).
    9 !
     9
    1010!   For some distribution types, the total number concentration (per kg), Np
    1111!   may be optionally specified.   Should be set to zero, otherwise.
    12 !
     12
    1313!   Roj Marchand July 2010
    1414
    1515
    1616! Inputs:
    17 !
     17
    1818!   [Q]        hydrometeor mixing ratio (g/kg)  ! not needed for some distribution types
    1919!   [Np]       Optional Total number concentration (per kg).  0 = use defaults (p1, p2, p3)
    20 !
     20
    2121!   [rho_a]    ambient air density (kg m^-3)   
    22 !
     22
    2323!   Distribution parameters as per quickbeam documentation.
    2424!   [dtype]    distribution type
     
    2828!   [bmp]      b params for mass
    2929!   [p1],[p2],[p3]  distribution parameters
    30 !
    31 !
     30
     31
    3232! Outputs:
    3333!   [Re]       Effective radius, 1/2 the 3rd moment/2nd moment (um)
    34 !
     34
    3535! Created:
    3636!   July 2010  Roj Marchand
    37 !
    38 
    39  
     37
    4038! ----- INPUTS ----- 
    4139 
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/congvec.h

    r2428 r5099  
    55! $Revision: 23 $, $Date: 2011-03-31 15:41:37 +0200 (jeu. 31 mars 2011) $
    66! $URL: http://cfmip-obs-sim.googlecode.com/svn/stable/v1.4.0/icarus-scops-4.1-bsd/congvec.f $
    7 !
     7
    88! Redistribution and use in source and binary forms, with or without
    99! modification, are permitted provided that the
    1010! following conditions are met:
    11 !
     11
    1212!     * Redistributions of source code must retain the above
    1313!       copyright  notice, this list of conditions and the following
     
    2121!       derived from this software without specific prior written
    2222!       permission.
    23 !
     23
    2424! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
    2525! "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
     
    3333! (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
    3434! OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 
    35 !
     35
    3636! *****************************COPYRIGHT*******************************
    3737! *****************************COPYRIGHT*******************************
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/cosp_defs.h

    r2431 r5099  
    33! $Revision: 88 $, $Date: 2013-11-13 15:08:38 +0100 (mer. 13 nov. 2013) $
    44! $URL: http://cfmip-obs-sim.googlecode.com/svn/stable/v1.4.0/cosp_defs.h $
    5 !
     5
    66! Redistribution and use in source and binary forms, with or without modification, are permitted
    77! provided that the following conditions are met:
    8 !
     8
    99!     * Redistributions of source code must retain the above copyright notice, this list
    1010!       of conditions and the following disclaimer.
     
    1515!       to endorse or promote products derived from this software without specific prior written
    1616!       permission.
    17 !
     17
    1818! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR
    1919! IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/dsd.F90

    r5095 r5099  
    1 !
     1
    22! $Id$
    3 !
     3
    44  subroutine dsd(Q,Re_,Np,D,N,nsizes,dtype,rho_a,tk, &
    55             dmin,dmax,apm,bpm,rho_c,p1,p2,p3)
     
    1010! Purpose:
    1111!   Create a discrete drop size distribution
    12 !
     12
    1313!   Starting with Quickbeam V3, this routine now allows input of
    1414!   both effective radius (Re) and total number concentration (Nt)
    1515!   Roj Marchand July 2010
    16 !
     16
    1717!   The version in Quickbeam v.104 was modified to allow Re but not Nt
    1818!   This is a significantly modified form for the version     
    19 !
     19
    2020!   Originally Part of QuickBeam v1.03 by John Haynes
    2121!   http://reef.atmos.colostate.edu/haynes/radarsim
    22 !
     22
    2323! Inputs:
    24 !
     24
    2525!   [Q]        hydrometeor mixing ratio (g/kg)
    2626!   [Re]       Optional Effective Radius (microns).  0 = use defaults (p1, p2, p3)
    27 !
     27
    2828!   [D]        array of discrete drop sizes (um) where we desire to know the number concentraiton n(D).
    2929!   [nsizes]   number of elements of [D]
    30 !
     30
    3131!   [dtype]    distribution type
    3232!   [rho_a]    ambient air density (kg m^-3)
     
    3636!   [rho_c]    alternate constant density (kg m^-3)
    3737!   [p1],[p2],[p3]  distribution parameters
    38 !
     38
    3939! Input/Output:
    4040!   [apm]      a parameter for mass (kg m^[-bpm])
    4141!   [bmp]      b params for mass
    42 !
     42
    4343! Outputs:
    4444!   [N]        discrete concentrations (cm^-3 um^-1)
    4545!              or, for monodisperse, a constant (1/cm^3)
    46 !
     46
    4747! Requires:
    4848!   function infind
    49 !
     49
    5050! Created:
    5151!   11/28/05  John Haynes (haynes@atmos.colostate.edu)
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/format_input.F90

    r5095 r5099  
    44! Compiled/Modified:
    55!   08/28/2006  John Haynes (haynes@atmos.colostate.edu)
    6 !
     6
    77! irreg_to_grid (subroutine)
    88! order_data (subroutine)
     
    2323!   Linearly interpolate sounding-level data to the hydrometeor-level
    2424!   resolution
    25 !
     25
    2626! Inputs:
    2727!   [hgt_matrix]       hydrometeor-level heights
     
    3030!   [env_p_matrix]     sounding-level pressures
    3131!   [env_rh_matrix]    sounding-level relative humidities
    32 !
     32
    3333! Outputs:
    3434!   [t_matrix]         hydrometeor-level temperatures
    3535!   [p_matrix]         hydrometeor-level pressures
    3636!   [rh_matrix]        hydrometeor-level relative humidities
    37 !
     37
    3838! Created:
    3939!   08/28/2006  John Haynes (haynes@atmos.colostate.edu)
     
    7373!   Ensure that input data is in top-down order/bottom-up order,
    7474!   for space-based/surface based radars, respectively
    75 !
     75
    7676! Inputs:
    7777!   [hgt_matrix]   heights
     
    8181!   [rh_matrix]    relative humidities
    8282!   [sfc_radar]    1=surface radar, 0=spaceborne
    83 !
     83
    8484! Outputs:
    8585!   [hgt_matrix],[hm_matrix],[p_matrix,[t_matrix],[rh_matrix] in proper
    8686!   order for input to the radar simulator routine
    8787!   [hgt_reversed]   T=heights were reordered,F=heights were not reordered
    88 !
     88
    8989! Note:
    9090!   The order for all profiles is assumed to the same as the first profile.
    91 !
     91
    9292! Created:
    9393!   08/28/2006  John Haynes (haynes@atmos.colostate.edu)
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/gases.F90

    r5095 r5099  
    66! Purpose:
    77!   Compute 2-way gaseous attenuation through a volume in microwave
    8 !
     8
    99! Inputs:
    1010!   [PRES_mb]   pressure (mb) (hPa)
     
    1212!   [RH]        relative humidity (%)
    1313!   [f]         frequency (GHz), < 300 GHz
    14 !
     14
    1515! Returns:
    1616!   2-way gaseous attenuation (dB/km)
    17 !
     17
    1818! Reference:
    1919!   Uses method of Liebe (1985)
    20 !
     20
    2121! Created:
    2222!   12/09/05  John Haynes (haynes@atmos.colostate.edu)
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/icarus.F

    r5095 r5099  
    4141! Corporation.
    4242! All rights reserved.
    43 !
     43
    4444! Redistribution and use in source and binary forms, with or without
    4545! modification, are permitted provided that the
    4646! following conditions are met:
    47 !
     47
    4848!     * Redistributions of source code must retain the above
    4949!       copyright  notice, this list of conditions and the following
     
    5858!       derived from this software without specific prior written
    5959!       permission.
    60 !
     60
    6161! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
    6262! "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
     
    7070! (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
    7171! OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 
    72 !
     72
    7373! *****************************COPYRIGHT*******************************
    7474! *****************************COPYRIGHT*******************************
     
    153153                                 ! with interpolated temperature equal to the radiance
    154154                                 ! determined cloud-top temperature
    155                                  !
     155
    156156                                 ! 1 = find the *lowest* altitude (highest pressure) level
    157157                                 ! with interpolated temperature equal to the radiance
    158158                                 ! determined cloud-top temperature
    159                                  !
     159
    160160                                 ! 2 = find the *highest* altitude (lowest pressure) level
    161161                                 ! with interpolated temperature equal to the radiance
    162162                                 ! determined cloud-top temperature
    163                                  !
     163
    164164                                 ! ONLY APPLICABLE IF top_height EQUALS 1 or 3
    165165                                 !                               !
    166166                                 ! 1 = old setting: matches all versions of
    167167                                 ! ISCCP simulator with versions numbers 3.5.1 and lower
    168                                  !
     168
    169169                                 ! 2 = default setting: for version numbers 4.0 and higher
    170 !
     170
    171171!     The following input variables are used only if top_height = 1 or top_height = 3
    172 !
     172
    173173      REAL skt(npoints)                 !  skin Temperature (K)
    174174      REAL emsfc_lw                     !  10.5 micron emissivity of surface (fraction)                                           
     
    222222     
    223223      REAL boxptop(npoints,ncol)        !  cloud top pressure (mb) in each column
    224                              
    225                                                                                          
    226 !
     224
    227225!     ------
    228226!     Working variables added when program updated to mimic Mark Webb's PV-Wave code
     
    456454      enddo
    457455
    458 !
    459456!     ---------------------------------------------------!
    460457
    461      
    462 !
    463458!     ---------------------------------------------------!
    464459!     COMPUTE CLOUD OPTICAL DEPTH FOR EACH COLUMN and
     
    502497              enddo
    503498          endif
    504 !
     499
    505500!     ---------------------------------------------------!
    506501
    507 
    508 
    509 !     
    510502!     ---------------------------------------------------!
    511503!     COMPUTE INFRARED BRIGHTNESS TEMPERUATRES
    512504!     AND CLOUD TOP TEMPERATURE SATELLITE SHOULD SEE
    513 !
     505
    514506!     again this is only done if top_height = 1 or 3
    515 !
     507
    516508!     fluxtop is the 10.5 micron radiance at the top of the
    517509!              atmosphere
     
    525517
    526518        !----------------------------------------------------------------------
    527         !   
     519
    528520        !             DO CLEAR SKY RADIANCE CALCULATION FIRST
    529         !
     521
    530522        !compute water vapor continuum emissivity
    531523        !this treatment follows Schwarkzopf and Ramasamy
     
    640632        enddo
    641633      endif
    642    
    643 
    644         !
     634
    645635        !           END OF CLEAR SKY CALCULATION
    646         !
     636
    647637        !----------------------------------------------------------------
    648638
     
    795785        !clouds which have partial emissivity and need the
    796786        !adjustment performed in this section
    797         !
     787
    798788      !If it turns out that the cloud brightness temperature
    799789      !is greater than 260K, then the liquid cloud conversion
    800790        !factor of 2.56 is used.
    801       !
     791
    802792        !Note that this is discussed on pages 85-87 of
    803793        !the ISCCP D level documentation (Rossow et al. 1996)
     
    915905!     ---------------------------------------------------!
    916906
    917 !     
    918907!     ---------------------------------------------------!
    919908!     DETERMINE CLOUD TOP PRESSURE
    920 !
     909
    921910!     again the 2 methods differ according to whether
    922911!     or not you use the physical cloud top pressure (top_height = 2)
    923912!     or the radiatively determined cloud top pressure (top_height = 1 or 3)
    924 !
    925913
    926914      !compute cloud top pressure
     
    1004992
    100599330    continue
    1006              
    1007 !
    1008 !
     994
     995
    1009996!     ---------------------------------------------------!
    1010997
    1011 
    1012 !     
    1013998!     ---------------------------------------------------!
    1014999!     DETERMINE ISCCP CLOUD TYPE FREQUENCIES
    1015 !
     1000
    10161001!     Now that ptop and tau have been determined,
    10171002!     determine amount of each of the 49 ISCCP cloud
    10181003!     types
    1019 !
     1004
    10201005!     Also compute grid box mean cloud top pressure and
    10211006!     optical thickness.  The mean cloud top pressure and
     
    10271012!     then converting the average albedo back to a mean
    10281013!     optical thickness. 
    1029 !
    10301014
    10311015      !compute isccp frequencies
     
    11871171        end if
    11881172      enddo ! j
    1189 !
     1173
    11901174!     ---------------------------------------------------!
    11911175
    11921176!     ---------------------------------------------------!
    11931177!     OPTIONAL PRINTOUT OF DATA TO CHECK PROGRAM
    1194 !
     1178
    11951179      if (debugcol.ne.0) then
    1196 !     
     1180
    11971181         do j=1,npoints,debugcol
    11981182
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/isccp_cloud_types.F

    r2432 r5099  
    4141! (c) British Crown Copyright 2009, the Met Office.
    4242! All rights reserved.
    43 !
     43
    4444! Redistribution and use in source and binary forms, with or without
    4545! modification, are permitted provided that the
    4646! following conditions are met:
    47 !
     47
    4848!     * Redistributions of source code must retain the above
    4949!       copyright  notice, this list of conditions and the following
     
    5757!       derived from this software without specific prior written
    5858!       permission.
    59 !
     59
    6060! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
    6161! "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
     
    6969! (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
    7070! OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 
    71 !
     71
    7272! *****************************COPYRIGHT*******************************
    7373! *****************************COPYRIGHT*******************************
     
    161161                                 ! with interpolated temperature equal to the radiance
    162162                                 ! determined cloud-top temperature
    163                                  !
     163
    164164                                 ! 1 = find the *lowest* altitude (highest pressure) level
    165165                                 ! with interpolated temperature equal to the radiance
    166166                                 ! determined cloud-top temperature
    167                                  !
     167
    168168                                 ! 2 = find the *highest* altitude (lowest pressure) level
    169169                                 ! with interpolated temperature equal to the radiance
    170170                                 ! determined cloud-top temperature
    171                                  !
     171
    172172                                 ! ONLY APPLICABLE IF top_height EQUALS 1 or 3
    173                                  !
     173
    174174                                 ! 1 = old setting: matches all versions of
    175175                                 ! ISCCP simulator with versions numbers 3.5.1 and lower
    176                                  !
     176
    177177                                 ! 2 = default setting: for version numbers 4.0 and higher 
    178 !
     178
    179179!     The following input variables are used only if top_height = 1 or top_height = 3
    180 !
     180
    181181      REAL skt(npoints)                 !  skin Temperature (K)
    182182      REAL emsfc_lw                     !  10.5 micron emissivity of surface (fraction)                                           
     
    230230     
    231231      REAL boxptop(npoints,ncol)        !  cloud top pressure (mb) in each column
    232                              
    233                                                                                          
    234 !
     232
    235233!     ------
    236234!     Working variables added when program updated to mimic Mark Webb's PV-Wave code
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/lidar_simulator.F90

    r5095 r5099  
    33! $Revision: 88 $, $Date: 2013-11-13 15:08:38 +0100 (mer. 13 nov. 2013) $
    44! $URL: http://cfmip-obs-sim.googlecode.com/svn/stable/v1.4.0/actsim/lidar_simulator.F90 $
    5 !
     5
    66! Redistribution and use in source and binary forms, with or without modification, are permitted
    77! provided that the following conditions are met:
    8 !
     8
    99!     * Redistributions of source code must retain the above copyright notice, this list
    1010!       of conditions and the following disclaimer.
     
    1515!       contributors may be used to endorse or promote products derived from this software without
    1616!       specific prior written permission.
    17 !
     17
    1818! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR
    1919! IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
     
    3131                , ls_radliq, ls_radice, cv_radliq, cv_radice &
    3232                , ice_type, pmol, pnorm, pnorm_perp_tot,tautot, refl )
    33 !
     33
    3434!---------------------------------------------------------------------------------
    3535! Purpose: To compute lidar signal from model-simulated profiles of cloud water
    3636!          and cloud fraction in each sub-column of each model gridbox.
    37 !
     37
    3838! References:
    3939! Chepfer H., S. Bony, D. Winker, M. Chiriaco, J.-L. Dufresne, G. Seze (2008),
    4040! Use of CALIPSO lidar observations to evaluate the cloudiness simulated by a
    4141! climate model, Geophys. Res. Lett., 35, L15704, doi:10.1029/2008GL034207.
    42 !
     42
    4343! Previous references:
    4444! Chiriaco et al, MWR, 2006; Chepfer et al., MWR, 2007
    45 !
     45
    4646! Contacts: Helene Chepfer (chepfer@lmd.polytechnique.fr), Sandrine Bony (bony@lmd.jussieu.fr)
    47 !
     47
    4848! May 2007: ActSim code of M. Chiriaco and H. Chepfer rewritten by S. Bony
    49 !
     49
    5050! May 2008, H. Chepfer:
    5151! - Units of pressure inputs: Pa
    5252! - Non Spherical particles : LS Ice NS coefficients, CONV Ice NS coefficients
    5353! - New input: ice_type (0=ice-spheres ; 1=ice-non-spherical)
    54 !
     54
    5555! June 2008, A. Bodas-Salcedo:
    5656! - Ported to Fortran 90 and optimisation changes
    57 !
     57
    5858! August 2008, J-L Dufresne:
    5959! - Optimisation changes (sum instructions suppressed)
    60 !
     60
    6161! October 2008, S. Bony,  H. Chepfer and J-L. Dufresne : 
    6262! - Interface with COSP v2.0:
     
    6565!      depolarisation diagnostic removed
    6666!      parasol (polder) reflectances (for 5 different solar zenith angles) added
    67 !
     67
    6868! December 2008, S. Bony,  H. Chepfer and J-L. Dufresne :
    6969! - Modification of the integration of the lidar equation.
    7070! - change the cloud detection threshold
    71 !
     71
    7272! April 2008, A. Bodas-Salcedo:
    7373! - Bug fix in computation of pmol and pnorm of upper layer
    74 !
     74
    7575! April 2008, J-L. Dufresne
    7676! - Bug fix in computation of pmol and pnorm, thanks to Masaki Satoh: a factor 2
    7777! was missing. This affects the ATB values but not the cloud fraction.
    78 !
     78
    7979! January 2013, G. Cesana and H. Chepfer:
    8080! - Add the perpendicular component of the backscattered signal (pnorm_perp_tot) in the arguments
     
    8383! Reference: Cesana G. and H. Chepfer (2013): Evaluation of the cloud water phase
    8484! in a climate model using CALIPSO-GOCCP, J. Geophys. Res., doi: 10.1002/jgrd.50376
    85 !
     85
    8686!---------------------------------------------------------------------------------
    87 !
     87
    8888! Inputs:
    8989!  npoints  : number of horizontal points
     
    106106!  ice_type : ice particle shape hypothesis (ice_type=0 for spheres, ice_type=1
    107107!             for non spherical particles)
    108 !
     108
    109109! Outputs:
    110110!  pmol : molecular attenuated backscatter lidar signal power (m^-1.sr^-1)
     
    113113!  tautot: optical thickess integrated from top to level z
    114114!  refl : parasol(polder) reflectance
    115 !
     115
    116116! Version 1.0 (June 2007)
    117117! Version 1.1 (May 2008)
     
    228228! as a function of the ATB for ice or liquid cloud particles derived from CALIPSO-GOCCP
    229229! observations at 120m vertical grid (Cesana and Chepfer, JGR, 2013).
    230 !
     230
    231231! Relationship between ATBice and ATBperp,ice for ice particles
    232232!  ATBperp,ice = Alpha*ATBice
     
    414414
    415415! Total signal (molecular + particules):
    416 !
    417 !
     416
     417
    418418! For performance reason on vector computers, the 2 following lines should not be used
    419419! and should be replace by the later one.
     
    426426           tautot(:,:) = tautot(:,:)  + tau_part(:,:,i)
    427427      ENDDO ! i
    428 !
     428
    429429!     Upper layer
    430430      pnorm(:,nlev) = betatot(:,nlev) / (2.*tautot(:,nlev)) &
     
    585585!---------------------------------------------------------
    586586!  Parasol/Polder module
    587 !
     587
    588588!  Purpose : Compute reflectance for one particular viewing direction
    589589!  and 5 solar zenith angles (calculation valid only over ocean)
     
    611611
    612612  END SUBROUTINE lidar_simulator
    613 !
     613
    614614!---------------------------------------------------------------------------------
    615 !
     615
    616616  SUBROUTINE parasol(npoints,nrefl,undef  &
    617617                       ,tautot_S_liq,tautot_S_ice  &
     
    621621!          of cloud water and cloud fraction in each sub-column of each model
    622622!          gridbox.
    623 !
    624 !
     623
     624
    625625! December 2008, S. Bony,  H. Chepfer and J-L. Dufresne :
    626626! - optimization for vectorization
    627 !
     627
    628628! Version 2.0 (October 2008)
    629629! Version 2.1 (December 2008)
     
    642642! outputs
    643643    REAL refl(npoints,nrefl)     ! Parasol reflectances
    644 !
     644
    645645! Local variables
    646646    REAL tautot_S(npoints)       ! cloud optical thickness, from TOA to surface
     
    702702    rlumA_mod=0
    703703    rlumB_mod=0
    704 !
     704
    705705    pi = ACOS(-1.0)
    706706    r_norm(:)=1./ COS(pi/180.*tetas(:))
    707 !
     707
    708708    tautot_S_liq(:)=MAX(tautot_S_liq(:),tau(1))
    709709    tautot_S_ice(:)=MAX(tautot_S_ice(:),tau(1))
    710710    tautot_S(:) = tautot_S_ice(:) + tautot_S_liq(:)
    711 !
     711
    712712! relative fraction of the opt. thick due to liquid or ice clouds
    713713    WHERE (tautot_S(:) .GT. 0.)
     
    719719    END WHERE
    720720    tautot_S(:)=MIN(tautot_S(:),tau(nbtau))
    721 !
     721
    722722! Linear interpolation :
    723723
     
    730730      bB(:,ny) = rlumB(:,ny) - aB(:,ny)*tau(ny)
    731731    ENDDO
    732 !
     732
    733733    DO it=1,ntetas
    734734      DO ny=1,nbtau-1
     
    739739      END DO
    740740    END DO
    741 !
     741
    742742    DO it=1,ntetas
    743743      refl(:,it) = frac_taucol_liq(:) * rlumA_mod(:,it) &
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/m_mrgrnk.F90

    r5095 r5099  
    2121! __________________________________________________________
    2222      Real (kind=kdp) :: XVALA, XVALB
    23 !
     23
    2424      Integer, Dimension (SIZE(IRNGT)) :: JWRKT
    2525      Integer :: LMTNA, LMTNC, IRNG1, IRNG2
    2626      Integer :: NVAL, IIND, IWRKD, IWRK, IWRKF, JINDA, IINDA, IINDB
    27 !
     27
    2828      NVAL = Min (SIZE(XDONT), SIZE(IRNGT))
    2929      Select Case (NVAL)
     
    3636         Continue
    3737      End Select
    38 !
     38
    3939!  Fill-in the index array, creating ordered couples
    40 !
     40
    4141      Do IIND = 2, NVAL, 2
    4242         If (XDONT(IIND-1) <= XDONT(IIND)) Then
     
    5151         IRNGT (NVAL) = NVAL
    5252      End If
    53 !
     53
    5454!  We will now have ordered subsets A - B - A - B - ...
    5555!  and merge A and B couples into     C   -   C   - ...
    56 !
     56
    5757      LMTNA = 2
    5858      LMTNC = 4
    59 !
     59
    6060!  First iteration. The length of the ordered subsets goes from 2 to 4
    61 !
     61
    6262      Do
    6363         If (NVAL <= 2) Exit
    64 !
     64
    6565!   Loop on merges of A and B into C
    66 !
     66
    6767         Do IWRKD = 0, NVAL - 1, 4
    6868            If ((IWRKD+4) > NVAL) Then
    6969               If ((IWRKD+2) >= NVAL) Exit
    70 !
     70
    7171!   1 2 3
    72 !
     72
    7373               If (XDONT(IRNGT(IWRKD+2)) <= XDONT(IRNGT(IWRKD+3))) Exit
    74 !
     74
    7575!   1 3 2
    76 !
     76
    7777               If (XDONT(IRNGT(IWRKD+1)) <= XDONT(IRNGT(IWRKD+3))) Then
    7878                  IRNG2 = IRNGT (IWRKD+2)
    7979                  IRNGT (IWRKD+2) = IRNGT (IWRKD+3)
    8080                  IRNGT (IWRKD+3) = IRNG2
    81 !
     81
    8282!   3 1 2
    83 !
     83
    8484               Else
    8585                  IRNG1 = IRNGT (IWRKD+1)
     
    9090               Exit
    9191            End If
    92 !
     92
    9393!   1 2 3 4
    94 !
     94
    9595            If (XDONT(IRNGT(IWRKD+2)) <= XDONT(IRNGT(IWRKD+3))) Cycle
    96 !
     96
    9797!   1 3 x x
    98 !
     98
    9999            If (XDONT(IRNGT(IWRKD+1)) <= XDONT(IRNGT(IWRKD+3))) Then
    100100               IRNG2 = IRNGT (IWRKD+2)
     
    108108                  IRNGT (IWRKD+4) = IRNG2
    109109               End If
    110 !
     110
    111111!   3 x x x
    112 !
     112
    113113            Else
    114114               IRNG1 = IRNGT (IWRKD+1)
     
    133133            End If
    134134         End Do
    135 !
     135
    136136!  The Cs become As and Bs
    137 !
     137
    138138         LMTNA = 4
    139139         Exit
    140140      End Do
    141 !
     141
    142142!  Iteration loop. Each time, the length of the ordered subsets
    143143!  is doubled.
    144 !
     144
    145145      Do
    146146         If (LMTNA >= NVAL) Exit
    147147         IWRKF = 0
    148148         LMTNC = 2 * LMTNC
    149 !
     149
    150150!   Loop on merges of A and B into C
    151 !
     151
    152152         Do
    153153            IWRK = IWRKF
     
    161161            IINDA = 1
    162162            IINDB = JINDA + 1
    163 !
     163
    164164!   Shortcut for the case when the max of A is smaller
    165165!   than the min of B. This line may be activated when the
    166166!   initial set is already close to sorted.
    167 !
     167
    168168!          IF (XDONT(IRNGT(JINDA)) <= XDONT(IRNGT(IINDB))) CYCLE
    169 !
     169
    170170!  One steps in the C subset, that we build in the final rank array
    171 !
     171
    172172!  Make a copy of the rank array for the merge iteration
    173 !
     173
    174174            JWRKT (1:LMTNA) = IRNGT (IWRKD:JINDA)
    175 !
     175
    176176            XVALA = XDONT (JWRKT(IINDA))
    177177            XVALB = XDONT (IRNGT(IINDB))
    178 !
     178
    179179            Do
    180180               IWRK = IWRK + 1
    181 !
     181
    182182!  We still have unprocessed values in both A and B
    183 !
     183
    184184               If (XVALA > XVALB) Then
    185185                  IRNGT (IWRK) = IRNGT (IINDB)
     
    197197                  XVALA = XDONT (JWRKT(IINDA))
    198198               End If
    199 !
     199
    200200            End Do
    201201         End Do
    202 !
     202
    203203!  The Cs become As and Bs
    204 !
     204
    205205         LMTNA = 2 * LMTNA
    206206      End Do
    207 !
     207
    208208      Return
    209 !
     209
    210210End Subroutine D_mrgrnk
    211211
     
    221221! __________________________________________________________
    222222      Integer :: XVALA, XVALB
    223 !
     223
    224224      Integer, Dimension (SIZE(IRNGT)) :: JWRKT
    225225      Integer :: LMTNA, LMTNC, IRNG1, IRNG2
    226226      Integer :: NVAL, IIND, IWRKD, IWRK, IWRKF, JINDA, IINDA, IINDB
    227 !
     227
    228228      NVAL = Min (SIZE(XDONT), SIZE(IRNGT))
    229229      Select Case (NVAL)
     
    236236         Continue
    237237      End Select
    238 !
     238
    239239!  Fill-in the index array, creating ordered couples
    240 !
     240
    241241      Do IIND = 2, NVAL, 2
    242242         If (XDONT(IIND-1) <= XDONT(IIND)) Then
     
    251251         IRNGT (NVAL) = NVAL
    252252      End If
    253 !
     253
    254254!  We will now have ordered subsets A - B - A - B - ...
    255255!  and merge A and B couples into     C   -   C   - ...
    256 !
     256
    257257      LMTNA = 2
    258258      LMTNC = 4
    259 !
     259
    260260!  First iteration. The length of the ordered subsets goes from 2 to 4
    261 !
     261
    262262      Do
    263263         If (NVAL <= 2) Exit
    264 !
     264
    265265!   Loop on merges of A and B into C
    266 !
     266
    267267         Do IWRKD = 0, NVAL - 1, 4
    268268            If ((IWRKD+4) > NVAL) Then
    269269               If ((IWRKD+2) >= NVAL) Exit
    270 !
     270
    271271!   1 2 3
    272 !
     272
    273273               If (XDONT(IRNGT(IWRKD+2)) <= XDONT(IRNGT(IWRKD+3))) Exit
    274 !
     274
    275275!   1 3 2
    276 !
     276
    277277               If (XDONT(IRNGT(IWRKD+1)) <= XDONT(IRNGT(IWRKD+3))) Then
    278278                  IRNG2 = IRNGT (IWRKD+2)
    279279                  IRNGT (IWRKD+2) = IRNGT (IWRKD+3)
    280280                  IRNGT (IWRKD+3) = IRNG2
    281 !
     281
    282282!   3 1 2
    283 !
     283
    284284               Else
    285285                  IRNG1 = IRNGT (IWRKD+1)
     
    290290               Exit
    291291            End If
    292 !
     292
    293293!   1 2 3 4
    294 !
     294
    295295            If (XDONT(IRNGT(IWRKD+2)) <= XDONT(IRNGT(IWRKD+3))) Cycle
    296 !
     296
    297297!   1 3 x x
    298 !
     298
    299299            If (XDONT(IRNGT(IWRKD+1)) <= XDONT(IRNGT(IWRKD+3))) Then
    300300               IRNG2 = IRNGT (IWRKD+2)
     
    308308                  IRNGT (IWRKD+4) = IRNG2
    309309               End If
    310 !
     310
    311311!   3 x x x
    312 !
     312
    313313            Else
    314314               IRNG1 = IRNGT (IWRKD+1)
     
    333333            End If
    334334         End Do
    335 !
     335
    336336!  The Cs become As and Bs
    337 !
     337
    338338         LMTNA = 4
    339339         Exit
    340340      End Do
    341 !
     341
    342342!  Iteration loop. Each time, the length of the ordered subsets
    343343!  is doubled.
    344 !
     344
    345345      Do
    346346         If (LMTNA >= NVAL) Exit
    347347         IWRKF = 0
    348348         LMTNC = 2 * LMTNC
    349 !
     349
    350350!   Loop on merges of A and B into C
    351 !
     351
    352352         Do
    353353            IWRK = IWRKF
     
    361361            IINDA = 1
    362362            IINDB = JINDA + 1
    363 !
     363
    364364!   Shortcut for the case when the max of A is smaller
    365365!   than the min of B. This line may be activated when the
    366366!   initial set is already close to sorted.
    367 !
     367
    368368!          IF (XDONT(IRNGT(JINDA)) <= XDONT(IRNGT(IINDB))) CYCLE
    369 !
     369
    370370!  One steps in the C subset, that we build in the final rank array
    371 !
     371
    372372!  Make a copy of the rank array for the merge iteration
    373 !
     373
    374374            JWRKT (1:LMTNA) = IRNGT (IWRKD:JINDA)
    375 !
     375
    376376            XVALA = XDONT (JWRKT(IINDA))
    377377            XVALB = XDONT (IRNGT(IINDB))
    378 !
     378
    379379            Do
    380380               IWRK = IWRK + 1
    381 !
     381
    382382!  We still have unprocessed values in both A and B
    383 !
     383
    384384               If (XVALA > XVALB) Then
    385385                  IRNGT (IWRK) = IRNGT (IINDB)
     
    397397                  XVALA = XDONT (JWRKT(IINDA))
    398398               End If
    399 !
     399
    400400            End Do
    401401         End Do
    402 !
     402
    403403!  The Cs become As and Bs
    404 !
     404
    405405         LMTNA = 2 * LMTNA
    406406      End Do
    407 !
     407
    408408      Return
    409 !
     409
    410410End Subroutine I_mrgrnk
    411411end module m_mrgrnk
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/math_lib.F90

    r5095 r5099  
    44! Compiled/Modified:
    55!   07/01/06  John Haynes (haynes@atmos.colostate.edu)
    6 !
     6
    77! gamma (function)
    88! path_integral (function)
     
    1919  function gamma(x)
    2020  implicit none
    21 !
     21
    2222! Purpose:
    2323!   Returns the gamma function
    24 !
     24
    2525! Input:
    2626!   [x]   value to compute gamma function of
    27 !
     27
    2828! Returns:
    2929!   gamma(x)
    30 !
     30
    3131! Coded:
    3232!   02/02/06  John Haynes (haynes@atmos.colostate.edu)
     
    100100  use array_lib
    101101  implicit none
    102 !
     102
    103103! Purpose:
    104104!   evalues the integral (f ds) between f(index=i1) and f(index=i2)
    105105!   using the AVINT procedure
    106 !
     106
    107107! Inputs:
    108108!   [f]    functional values
     
    110110!   [i1]   index of lower limit
    111111!   [i2]   index of upper limit
    112 !
     112
    113113! Returns:
    114114!   result of path integral
    115 !
     115
    116116! Notes:
    117117!   [s] may be in forward or reverse numerical order
    118 !
     118
    119119! Requires:
    120120!   mrgrnk package
    121 !
     121
    122122! Created:
    123123!   02/02/06  John Haynes (haynes@atmos.colostate.edu)
     
    162162  subroutine avint ( ftab, xtab, ntab, a_in, b_in, result )
    163163  implicit none
    164 !
     164
    165165! Purpose:
    166166!   estimate the integral of unevenly spaced data
    167 !
     167
    168168! Inputs:
    169169!   [ftab]     functional values
     
    172172!   [a]        lower limit of integration
    173173!   [b]        upper limit of integration
    174 !
     174
    175175! Outputs:
    176176!   [result]   approximate value of integral
    177 !
     177
    178178! Reference:
    179179!   From SLATEC libraries, in public domain
    180 !
     180
    181181!***********************************************************************
    182 !
     182
    183183!  AVINT estimates the integral of unevenly spaced data.
    184 !
     184
    185185!  Discussion:
    186 !
     186
    187187!    The method uses overlapping parabolas and smoothing.
    188 !
     188
    189189!  Modified:
    190 !
     190
    191191!    30 October 2000
    192192!    4 January 2008, A. Bodas-Salcedo. Error control for XTAB taken out of
    193193!                    loop to allow vectorization.
    194 !
     194
    195195!  Reference:
    196 !
     196
    197197!    Philip Davis and Philip Rabinowitz,
    198198!    Methods of Numerical Integration,
    199199!    Blaisdell Publishing, 1967.
    200 !
     200
    201201!    P E Hennion,
    202202!    Algorithm 77,
     
    204204!    Communications of the Association for Computing Machinery,
    205205!    Volume 5, page 96, 1962.
    206 !
     206
    207207!  Parameters:
    208 !
     208
    209209!    Input, real ( kind = 8 ) FTAB(NTAB), the function values,
    210210!    FTAB(I) = F(XTAB(I)).
    211 !
     211
    212212!    Input, real ( kind = 8 ) XTAB(NTAB), the abscissas at which the
    213213!    function values are given.  The XTAB's must be distinct
    214214!    and in ascending order.
    215 !
     215
    216216!    Input, integer NTAB, the number of entries in FTAB and
    217217!    XTAB.  NTAB must be at least 3.
    218 !
     218
    219219!    Input, real ( kind = 8 ) A, the lower limit of integration.  A should
    220220!    be, but need not be, near one endpoint of the interval
    221221!    (X(1), X(NTAB)).
    222 !
     222
    223223!    Input, real ( kind = 8 ) B, the upper limit of integration.  B should
    224224!    be, but need not be, near one endpoint of the interval
    225225!    (X(1), X(NTAB)).
    226 !
     226
    227227!    Output, real ( kind = 8 ) RESULT, the approximate value of the integral.
    228228
     
    293293    return
    294294  end if
    295 !
     295
    296296!  If B < A, temporarily switch A and B, and store sign.
    297 !
     297
    298298  if ( b < a ) then
    299299    syl = b
     
    305305    ind = 1
    306306  end if
    307 !
     307
    308308!  Bracket A and B between XTAB(ILO) and XTAB(IHI).
    309 !
     309
    310310  ilo = 1
    311311  ihi = ntab
     
    330330  ihi = min ( ihi, ntab - 1 )
    331331  ihi = max ( ilo, ihi - 1 )
    332 !
     332
    333333!  Carry out approximate integration from XTAB(ILO) to XTAB(IHI).
    334 !
     334
    335335  sum1 = 0.0D+00
    336336 
     
    380380        + cb * 0.5D+00 * ( b**2 - syl**2 ) &
    381381        + cc * ( b - syl )
    382 !
     382
    383383!  Restore original values of A and B, reverse sign of integral
    384384!  because of earlier switch.
    385 !
     385
    386386  if ( ind /= 1 ) then
    387387    ind = 1
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/mod_cosp.F90

    r5095 r5099  
    11! (c) British Crown Copyright 2008, the Met Office.
    22! All rights reserved.
    3 !
     3
    44! Redistribution and use in source and binary forms, with or without modification, are permitted
    55! provided that the following conditions are met:
    6 !
     6
    77!     * Redistributions of source code must retain the above copyright notice, this list
    88!       of conditions and the following disclaimer.
     
    1313!       to endorse or promote products derived from this software without specific prior written
    1414!       permission.
    15 !
     15
    1616! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR
    1717! IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/mod_cosp_constants.F90

    r4785 r5099  
    11! (c) British Crown Copyright 2008, the Met Office.
    22! All rights reserved.
    3 !
     3
    44! Redistribution and use in source and binary forms, with or without modification, are permitted
    55! provided that the following conditions are met:
    6 !
     6
    77!     * Redistributions of source code must retain the above copyright notice, this list
    88!       of conditions and the following disclaimer.
     
    1313!       to endorse or promote products derived from this software without specific prior written
    1414!       permission.
    15 !
     15
    1616! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR
    1717! IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
     
    2323! OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
    2424
    25 !
    2625! History:
    2726! Jul 2007 - A. Bodas-Salcedo - Initial version
     
    2928! Oct 2008 - H. Chepfer       - Added PARASOL_NREFL
    3029! Jun 2010 - R. Marchand      - Modified to support quickbeam V3, added ifdef for hydrometeor definitions
    31 !
    32 !
    33 !
     30
     31
    3432
    3533!!#INCLUDE "cosp_defs.h"
     
    143141                                             shape=(/2,MISR_N_CTH/))
    144142
    145 
    146     !
    147143    ! The following code was modifed by Roj with implementation of quickbeam V3
    148144    !   (1) use ifdef to support more than one microphyscis scheme
    149145    !   (2) added constants  microphysic_scheme_name, LOAD_scale_LUTs, and SAVE_scale_LUTs
    150     !
    151146
    152147    ! directory where LUTs will be stored
     
    155150#ifdef MMF_V3_SINGLE_MOMENT
    156151
    157     !       
    158152    !  Table hclass for quickbeam to support one-moment (bulk) microphysics scheme used by MMF V3.0 & V3.5
    159     !
    160 
    161     !
     153
    162154    ! NOTE:  if ANY value in this section of code is changed, the existing LUT
    163155    !        (i.e., the associated *.dat file) MUST be deleted so that a NEW
    164156    !        LUT will be created !!!
    165     !
     157
    166158    character*120 :: RADAR_SIM_MICROPHYSICS_SCHEME_NAME = 'MMF_v3_single_moment'
    167159
     
    190182
    191183    ! NOTES on HCLASS variables
    192     !
     184
    193185    ! TYPE - Set to
    194186    ! 1 for modified gamma distribution,
     
    217209    ! P1, P2, P3 - are default distribution parameters that depend on the type
    218210    ! of distribution (see quickmbeam documentation for more information)
    219     !
     211
    220212    ! Modified Gamma (must set P3 and one of P1 or P2)
    221213    ! P1 - Set to the total particle number concentration Nt /rho_a (kg-1 ), where
     
    223215    ! P2 - Set to the particle mean diameter D (micron).
    224216    ! P3 - Set to the distribution width nu.
    225     !
     217
    226218    ! Exponetial (set one of)
    227219    ! P1 - Set to a constant intercept parameter N0 (m-4).
    228220    ! P2 - Set to a constant lambda (micron-1).
    229     !
     221
    230222    ! Power Law
    231223    ! P1 - Set this to the value of a constant power law parameter br
    232     !
     224
    233225    ! Monodisperse
    234226    ! P1 - Set to a constant diameter D0 (micron) = Re.
    235     !
     227
    236228    ! Log-normal (must set P3 and one of P1 or P2)
    237229    ! P1 - Set to the total particle number concentration Nt /rho_a (kg-1 )
    238230    ! P2 - Set to the geometric mean particle radius rg (micron).
    239231    ! P3 - Set to the natural logarithm of the geometric standard deviation.
    240     !
    241 
    242232
    243233    real,dimension(N_HYDRO) :: N_ax,N_bx,alpha_x,c_x,d_x,g_x,a_x,b_x,gamma_1,gamma_2,gamma_3,gamma_4
     
    265255#ifdef MMF_V3p5_TWO_MOMENT
    266256
    267     !
    268257    !  Table hclass for quickbeam to support two-moment "morrison" microphysics scheme used by V3.5 (SAM 6.8)
    269     !
     258
    270259    !  This Number concentriation Np in [1/kg] MUST be input to COSP/radar simulator
    271     !
     260
    272261    !  NOTE:  Be sure to check that the ice-density (rho) set it this tables matches what you used
    273     !
    274 
    275     !
     262
    276263    ! NOTE:  if ANY value in this section of code is changed, the existing LUT
    277264    !        (i.e., the associated *.dat file) MUST be deleted so that a NEW
    278265    !        LUT will be created !!!
    279     !
     266
    280267    character*120 :: RADAR_SIM_MICROPHYSICS_SCHEME_NAME = 'MMF_v3.5_two_moment'
    281268
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/mod_cosp_isccp_simulator.F90

    r3233 r5099  
    33! $Revision: 23 $, $Date: 2011-03-31 15:41:37 +0200 (jeu. 31 mars 2011) $
    44! $URL: http://cfmip-obs-sim.googlecode.com/svn/stable/v1.4.0/cosp_isccp_simulator.F90 $
    5 !
     5
    66! Redistribution and use in source and binary forms, with or without modification, are permitted
    77! provided that the following conditions are met:
    8 !
     8
    99!     * Redistributions of source code must retain the above copyright notice, this list
    1010!       of conditions and the following disclaimer.
     
    1515!       to endorse or promote products derived from this software without specific prior written
    1616!       permission.
    17 !
     17
    1818! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR
    1919! IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/mod_cosp_lidar.F90

    r3233 r5099  
    33! $Revision: 88 $, $Date: 2013-11-13 15:08:38 +0100 (mer. 13 nov. 2013) $
    44! $URL: http://cfmip-obs-sim.googlecode.com/svn/stable/v1.4.0/cosp_lidar.F90 $
    5 !
     5
    66! Redistribution and use in source and binary forms, with or without modification, are permitted
    77! provided that the following conditions are met:
    8 !
     8
    99!     * Redistributions of source code must retain the above copyright notice, this list
    1010!       of conditions and the following disclaimer.
     
    1515!       to endorse or promote products derived from this software without specific prior written
    1616!       permission.
    17 !
     17
    1818! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR
    1919! IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
     
    2525! OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
    2626
    27 !
    2827! History:
    2928! Jul 2007 - A. Bodas-Salcedo - Initial version
     
    3231!                               frac_out changed in sgx%frac_out)
    3332! Jun 2011 - G. Cesana        - Added betaperp_tot argument
    34 !
    35 !
     33
     34
    3635MODULE MOD_COSP_LIDAR
    3736  USE MOD_COSP_CONSTANTS
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/mod_cosp_misr_simulator.F90

    r3233 r5099  
    33! $Revision: 23 $, $Date: 2011-03-31 15:41:37 +0200 (jeu. 31 mars 2011) $
    44! $URL: http://cfmip-obs-sim.googlecode.com/svn/stable/v1.4.0/cosp_misr_simulator.F90 $
    5 !
     5
    66! Redistribution and use in source and binary forms, with or without modification, are permitted
    77! provided that the following conditions are met:
    8 !
     8
    99!     * Redistributions of source code must retain the above copyright notice, this list
    1010!       of conditions and the following disclaimer.
     
    1515!       to endorse or promote products derived from this software without specific prior written
    1616!       permission.
    17 !
     17
    1818! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR
    1919! IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
     
    2525! OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
    2626
    27 !
    2827! History:
    2928! Nov 2008 - A. Bodas-Salcedo - Initial version
    30 !
    31 !
     29
    3230
    3331MODULE MOD_COSP_MISR_SIMULATOR
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/mod_cosp_modis_simulator.F90

    r5095 r5099  
    44! $Revision: 88 $, $Date: 2013-11-13 07:08:38 -0700 (Wed, 13 Nov 2013) $
    55! $URL: http://cfmip-obs-sim.googlecode.com/svn/stable/v1.4.0/cosp_modis_simulator.F90 $
    6 !
     6
    77! Redistribution and use in source and binary forms, with or without modification, are permitted
    88! provided that the following conditions are met:
    9 !
     9
    1010!     * Redistributions of source code must retain the above copyright notice, this list
    1111!       of conditions and the following disclaimer.
     
    1616!       to endorse or promote products derived from this software without specific prior written
    1717!       permission.
    18 !
     18
    1919! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR
    2020! IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
     
    2525! IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT
    2626! OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
    27 !
    28 
    29 !
     27
    3028! History:
    3129!   May 2009 - Robert Pincus - Initial version
    3230!   Dec 2009 - Robert Pincus - Tiny revisions
    33 !
     31
    3432MODULE MOD_COSP_Modis_Simulator
    3533  USE MOD_COSP_CONSTANTS
     
    4442  !------------------------------------------------------------------------------------------------
    4543  ! Public type
    46   !
     44
    4745  ! Summary statistics from MODIS retrievals
    4846  type COSP_MODIS
    4947     ! Dimensions
    5048     integer :: Npoints   ! Number of gridpoints
    51      
    52      !
     49
    5350     ! Grid means; dimension nPoints
    54      !
     51
    5552     real, dimension(:),       pointer :: &
    5653       Cloud_Fraction_Total_Mean,       Cloud_Fraction_Water_Mean,       Cloud_Fraction_Ice_Mean,       &
     
    6158       Cloud_Top_Pressure_Total_Mean,                                                                   &
    6259                                        Liquid_Water_Path_Mean,          Ice_Water_Path_Mean
    63      !
     60
    6461     ! Also need the ISCCP-type optical thickness/cloud top pressure histogram
    65      !
     62
    6663     real, dimension(:, :, :), pointer :: Optical_Thickness_vs_Cloud_Top_Pressure
    6764     real, dimension(:, :, :), pointer :: Optical_Thickness_vs_ReffICE
     
    126123    integer, dimension(count(gridBox%sunlit(:) <= 0)) :: notSunlit
    127124    ! ------------------------------------------------------------
    128    
    129     !
     125
    130126    ! Are there any sunlit points?
    131     !
     127
    132128    nSunlit = count(gridBox%sunlit(:) > 0)
    133129    if(nSunlit > 0) then
     
    135131      nPoints  = gridBox%Npoints
    136132      nSubCols = subCols%Ncolumns
    137       !
     133
    138134      ! This is a vector index indicating which points are sunlit
    139       !
     135
    140136      sunlit(:)    = pack((/ (i, i = 1, nPoints ) /), mask =       gridBox%sunlit(:) > 0)
    141137      notSunlit(:) = pack((/ (i, i = 1, nPoints ) /), mask = .not. gridBox%sunlit(:) > 0)
    142                
    143       !
     138
    144139      ! Copy needed quantities, reversing vertical order and removing points with no sunlight
    145       !
     140
    146141      pressureLevels(:, 1) = 0.0 ! Top of model, following ISCCP sim
    147142      temperature(:, :)     = gridBox%T (sunlit(:), nLevels:1:-1)
    148143      pressureLayers(:, :)  = gridBox%p (sunlit(:), nLevels:1:-1)
    149144      pressureLevels(:, 2:) = gridBox%ph(sunlit(:), nLevels:1:-1)
    150      
    151       !
     145
    152146      ! Subcolumn properties - first stratiform cloud...
    153       !
     147
    154148      where(subCols%frac_out(sunlit(:), :, :) == I_LSC)
    155149        !opticalThickness(:, :, :) = &
     
    180174      end do
    181175
    182       !
    183176      ! .. then add convective cloud...
    184       !
     177
    185178      where(subCols%frac_out(sunlit(:), :, :) == I_CVC)
    186179        !opticalThickness(:, :, :) = &
     
    201194      end do
    202195
    203       !
    204196      ! Reverse vertical order
    205       !
     197
    206198      opticalThickness(:, :, :)  = opticalThickness(:, :, nLevels:1:-1)
    207199      cloudWater      (:, :, :)  = cloudWater      (:, :, nLevels:1:-1)
     
    242234                        jointHistogram,jointHistogram2,jointHistogram3)
    243235      ! DJS2015: END
    244      
    245       !
     236
    246237      ! Copy results into COSP structure
    247       !
     238
    248239      modisSim%Cloud_Fraction_Total_Mean(sunlit(:)) = cfTotal(:)
    249240      modisSim%Cloud_Fraction_Water_Mean(sunlit(:)) = cfLiquid
     
    273264      modisSim%Optical_Thickness_vs_ReffICE(sunlit(:),2:numModisTauBins+1,:)              = jointHistogram2(:, :, :)
    274265      modisSim%Optical_Thickness_vs_ReffLIQ(sunlit(:),2:numModisTauBins+1,:)              = jointHistogram3(:, :, :)
    275       !
     266
    276267      ! Reorder pressure bins in joint histogram to go from surface to TOA
    277       !
     268
    278269      modisSim%Optical_Thickness_vs_Cloud_Top_Pressure(:,:,:) = modisSim%Optical_Thickness_vs_Cloud_Top_Pressure(:, :, numModisPressureBins:1:-1)
    279270      if(nSunlit < nPoints) then
    280         !
     271
    281272        ! Where it's night and we haven't done the retrievals the values are undefined
    282         !
     273
    283274        modisSim%Cloud_Fraction_Total_Mean(notSunlit(:)) = R_UNDEF
    284275        modisSim%Cloud_Fraction_Water_Mean(notSunlit(:)) = R_UNDEF
     
    310301      end if
    311302    else
    312       !
     303
    313304      ! It's nightime everywhere - everything is undefined
    314       !
     305
    315306      modisSim%Cloud_Fraction_Total_Mean(:) = R_UNDEF
    316307      modisSim%Cloud_Fraction_Water_Mean(:) = R_UNDEF
     
    350341    integer,           intent(in)  :: Npoints  ! Number of sampled points
    351342    type(cosp_MODIS),  intent(out) :: x
    352     !
     343
    353344    ! Allocate minumum storage if simulator not used
    354     !
     345
    355346    if (cfg%LMODIS_sim) then
    356347      x%nPoints  = nPoints
     
    398389  SUBROUTINE FREE_COSP_MODIS(x)
    399390    type(cosp_MODIS),intent(inout) :: x
    400     !
     391
    401392    ! Free space used by cosp_modis variable.
    402     !
    403    
     393
    404394    if(associated(x%Cloud_Fraction_Total_Mean)) deallocate(x%Cloud_Fraction_Total_Mean)
    405395    if(associated(x%Cloud_Fraction_Water_Mean)) deallocate(x%Cloud_Fraction_Water_Mean)
     
    439429    type(cosp_modis),      intent(in   ) :: orig
    440430    type(cosp_modis),      intent(  out) :: copy
    441     !
     431
    442432    ! Copy a set of grid points from one cosp_modis variable to another.
    443433    !   Should test to be sure ix and iy refer to the same number of grid points
    444     !
     434
    445435    integer :: orig_start, orig_end, copy_start, copy_end
    446436   
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/mod_cosp_radar.F90

    r5095 r5099  
    11! (c) British Crown Copyright 2008, the Met Office.
    22! All rights reserved.
    3 !
     3
    44! Redistribution and use in source and binary forms, with or without modification, are permitted
    55! provided that the following conditions are met:
    6 !
     6
    77!     * Redistributions of source code must retain the above copyright notice, this list
    88!       of conditions and the following disclaimer.
     
    1313!       to endorse or promote products derived from this software without specific prior written
    1414!       permission.
    15 !
     15
    1616! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR
    1717! IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/mod_cosp_simulator.F90

    r4619 r5099  
    33! $Revision: 88 $, $Date: 2013-11-13 15:08:38 +0100 (mer. 13 nov. 2013) $
    44! $URL: http://cfmip-obs-sim.googlecode.com/svn/stable/v1.4.0/cosp_simulator.F90 $
    5 !
     5
    66! Redistribution and use in source and binary forms, with or without modification, are permitted
    77! provided that the following conditions are met:
    8 !
     8
    99!     * Redistributions of source code must retain the above copyright notice, this list
    1010!       of conditions and the following disclaimer.
     
    1515!       to endorse or promote products derived from this software without specific prior written
    1616!       permission.
    17 !
     17
    1818! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR
    1919! IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
     
    2525! OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
    2626
    27 !
    2827! History:
    2928! Jul 2007 - A. Bodas-Salcedo - Initial version
    3029! Jan 2013 - G. Cesana - Add new variables linked to the lidar cloud phase
    31 !
    3230
    3331#include "cosp_defs.h"
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/mod_cosp_stats.F90

    r4619 r5099  
    33! $Revision: 88 $, $Date: 2013-11-13 15:08:38 +0100 (mer. 13 nov. 2013) $
    44! $URL: http://cfmip-obs-sim.googlecode.com/svn/stable/v1.4.0/cosp_stats.F90 $
    5 !
     5
    66! Redistribution and use in source and binary forms, with or without modification, are permitted
    77! provided that the following conditions are met:
    8 !
     8
    99!     * Redistributions of source code must retain the above copyright notice, this list
    1010!       of conditions and the following disclaimer.
     
    1515!       to endorse or promote products derived from this software without specific prior written
    1616!       permission.
    17 !
     17
    1818! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR
    1919! IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
     
    2525! OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
    2626
    27 !
    2827! History:
    2928! Jul 2007 - A. Bodas-Salcedo - Initial version
     
    3433! Jan 2013 - G. Cesana        - Added betaperp and temperature arguments
    3534!                             - Added phase 3D/3Dtemperature/Map output variables in diag_lidar
    36 !
    37 !
     35
     36
    3837#include "cosp_defs.h"
    3938MODULE MOD_COSP_STATS
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/mod_cosp_types.F90

    r5095 r5099  
    11! (c) British Crown Copyright 2008, the Met Office.
    22! All rights reserved.
    3 !
     3
    44! Redistribution and use in source and binary forms, with or without modification, are permitted
    55! provided that the following conditions are met:
    6 !
     6
    77!     * Redistributions of source code must retain the above copyright notice, this list
    88!       of conditions and the following disclaimer.
     
    1313!       to endorse or promote products derived from this software without specific prior written
    1414!       permission.
    15 !
     15
    1616! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR
    1717! IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/mod_cosp_utils.F90

    r3241 r5099  
    33! $Revision: 23 $, $Date: 2011-03-31 15:41:37 +0200 (jeu. 31 mars 2011) $
    44! $URL: http://cfmip-obs-sim.googlecode.com/svn/stable/v1.4.0/cosp_utils.F90 $
    5 !
     5
    66! Redistribution and use in source and binary forms, with or without modification, are permitted
    77! provided that the following conditions are met:
    8 !
     8
    99!     * Redistributions of source code must retain the above copyright notice, this list
    1010!       of conditions and the following disclaimer.
     
    1515!       to endorse or promote products derived from this software without specific prior written
    1616!       permission.
    17 !
     17
    1818! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR
    1919! IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
     
    2525! OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
    2626
    27 !
    2827! History:
    2928! Jul 2007 - A. Bodas-Salcedo - Initial version
    30 !
    3129
    3230MODULE MOD_COSP_UTILS
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/mod_llnl_stats.F90

    r5095 r5099  
    33! $Revision: 88 $, $Date: 2013-11-13 15:08:38 +0100 (mer. 13 nov. 2013) $
    44! $URL: http://cfmip-obs-sim.googlecode.com/svn/stable/v1.4.0/llnl/llnl_stats.F90 $
    5 !
     5
    66! Redistribution and use in source and binary forms, with or without modification, are permitted
    77! provided that the following conditions are met:
    8 !
     8
    99!     * Redistributions of source code must retain the above copyright notice, this list
    1010!       of conditions and the following disclaimer.
     
    1515!       nor the names of its contributors may be used to endorse or promote products derived from
    1616!       this software without specific prior written permission.
    17 !
     17
    1818! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR
    1919! IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
     
    2424! IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT
    2525! OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
    26 !
     26
    2727! History
    28 !
     28
    2929! Jan 2013 - G. Cesana        - Added betaperp_tot and temp_tot arguments
    30 !
    31 
    3230
    3331MODULE MOD_LLNL_STATS
     
    6260   ! bmin: mimumum value of first bin
    6361   ! bwidth: bin width
    64    !
     62
    6563   ! Output: 2D histogram on each horizontal point (Npoints,Nbins,Nlevels)
    6664   
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/mod_lmd_ipsl_stats.F90

    r5095 r5099  
    33! $Revision: 88 $, $Date: 2013-11-13 15:08:38 +0100 (mer. 13 nov. 2013) $
    44! $URL: http://cfmip-obs-sim.googlecode.com/svn/stable/v1.4.0/actsim/lmd_ipsl_stats.F90 $
    5 !
     5
    66! Redistribution and use in source and binary forms, with or without modification, are permitted
    77! provided that the following conditions are met:
    8 !
     8
    99!     * Redistributions of source code must retain the above copyright notice, this list
    1010!       of conditions and the following disclaimer.
     
    1515!       contributors may be used to endorse or promote products derived from this software without
    1616!       specific prior written permission.
    17 !
     17
    1818! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR
    1919! IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
     
    3838                  ,cfad2,srbval,ncat,ntype,lidarcld,lidarcldtype,lidarcldphase,cldlayer & !OPAQ
    3939                  ,cldtype,cldlayerphase,lidarcldtmp,parasolrefl,vgrid_z,profSR)          !OPAQ !TIBO
    40 !
     40
    4141! -----------------------------------------------------------------------------------
    4242! Lidar outputs :
    43 !
     43
    4444! Diagnose cloud fraction (3D cloud fraction + low/middle/high/total cloud fraction)
    4545! and phase cloud fraction (3D, low/mid/high/total and 3D temperature)
     
    4747!      +
    4848! Compute CFADs of lidar scattering ratio SR and of depolarization index
    49 !
     49
    5050! Authors: Sandrine Bony and Helene Chepfer (LMD/IPSL, CNRS, UPMC, France).
    51 !
     51
    5252! December 2008, S. Bony,  H. Chepfer and J-L. Dufresne :
    5353! - change of the cloud detection threshold S_cld from 3 to 5, for better
     
    6262! June 2010, T. Yokohata, T. Nishimura and K. Ogochi
    6363! Optimisation of COSP_CFAD_SR
    64 !
     64
    6565! January 2013, G. Cesana, H. Chepfer:
    6666! - Add the perpendicular component of the backscattered signal (pnorm_perp) in the arguments
     
    7373! Reference: Cesana G. and H. Chepfer (2013): Evaluation of the cloud water phase
    7474! in a climate model using CALIPSO-GOCCP, J. Geophys. Res., doi: 10.1002/jgrd.50376
    75 !
     75
    7676! ------------------------------------------------------------------------------------
    7777
     
    133133      real, parameter  ::  S_cld_att = 30. ! New threshold for undefine cloud phase detection   
    134134
    135 
    136 !
    137135! c -------------------------------------------------------
    138136! c 0- Initializations
    139137! c -------------------------------------------------------
    140 !
     138
    141139!  Should be modified in future version
    142140      xmax=undef-1.0
     
    174172! c -------------------------------------------------------
    175173      if (ok_lidar_cfad) then
    176 !
     174
    177175! c CFADs of subgrid-scale lidar scattering ratios :
    178176! c -------------------------------------------------------
     
    224222! S_att: Threshold for full attenuation
    225223! S_clr: Threshold for clear-sky layer
    226 !
     224
    227225!--- Input-Outout arguments
    228226! x: variable to process (Npoints,Ncolumns,Nlevels), mofified where saturation occurs
    229 !
     227
    230228! -- Output arguments
    231229! srbval : values of the histogram bins
     
    528526          enddo
    529527          endif
    530 !
    531528
    532529          iz=1
     
    600597
    601598!____________________________________________________________________________________________________
    602 !
     599
    603600!4.1.a Ice: ATBperp above the phase discrimination line
    604601!____________________________________________________________________________________________________
    605 !
     602
    606603           if( (ATBperp(i,ncol,nlev)-ATBperp_tmp).ge.0. )then   ! Ice clouds
    607604             ! ICE with temperature above 273,15°K = Liquid (false ice)
     
    645642
    646643!____________________________________________________________________________________________________
    647 !
     644
    648645! 4.1.b Liquid: ATBperp below the phase discrimination line
    649646!____________________________________________________________________________________________________
    650 !
     647
    651648             else                                        ! Liquid clouds
    652649              ! Liquid with temperature above 231,15°K
     
    710707                          ATB(i,ncol,nlev)*epsilon50 + zeta50
    711708!____________________________________________________________________________________________________
    712 !
     709
    713710! 4.2.a Ice: ATBperp above the phase discrimination line
    714711!____________________________________________________________________________________________________
    715 !
     712
    716713            ! ICE with temperature above 273,15°K = Liquid (false ice)
    717714          if( (ATBperp(i,ncol,nlev)-ATBperp_tmp).ge.0. )then   ! Ice clouds
     
    756753
    757754!____________________________________________________________________________________________________
    758 !
     755
    759756! 4.2.b Liquid: ATBperp below the phase discrimination line
    760757!____________________________________________________________________________________________________
    761 !
     758
    762759          else 
    763760             ! Liquid with temperature above 231,15°K
     
    816813
    817814!____________________________________________________________________________________________________
    818 !
     815
    819816! Undefined phase: For a cloud located below another cloud with SR>30
    820817! see Cesana and Chepfer 2013 Sect.III.2
    821818!____________________________________________________________________________________________________
    822 !
     819
    823820if(toplvlsat.ne.0)then         
    824821      do nlev=toplvlsat,1,-1
     
    849846
    850847!____________________________________________________________________________________________________
    851 !
     848
    852849! Computation of final cloud phase diagnosis
    853850!____________________________________________________________________________________________________
    854 !
    855851
    856852! Compute the Ice percentage in cloud = ice/(ice+liq) as a function
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/mod_modis_sim.F90

    r5095 r5099  
    44! $Revision: 88 $, $Date: 2013-11-13 07:08:38 -0700 (Wed, 13 Nov 2013) $
    55! $URL: http://cfmip-obs-sim.googlecode.com/svn/stable/v1.4.0/MODIS_simulator/modis_simulator.F90 $
    6 !
     6
    77! Redistribution and use in source and binary forms, with or without modification, are permitted
    88! provided that the following conditions are met:
    9 !
     9
    1010!     * Redistributions of source code must retain the above copyright notice, this list
    1111!       of conditions and the following disclaimer.
     
    1616!       to endorse or promote products derived from this software without specific prior written
    1717!       permission.
    18 !
     18
    1919! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR
    2020! IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
     
    2525! IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT
    2626! OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
    27 !
    28 
    29 !
     27
    3028! History:
    3129!   May 2009 - Robert Pincus - Initial version
     
    3432!   November 2009 - Robert Pincus - Bux fixes and speed-ups after experience with Rick Hemler using AM2 (GFDL)
    3533!   January 2010 - Robert Pincus - Added high, middle, low cloud fractions
    36 !
    37 
    38 !
     34
    3935! Notes on using the MODIS simulator:
    4036!  *) You may provide either layer-by-layer values of optical thickness at 0.67 and 2.1 microns, or
     
    4743!     libraries (parameters re_water_min, re_ice_min, etc.) the simulator returns missing values (re_fill)
    4844
    49 !
    5045! When error conditions are encountered this code calls the function complain_and_die, supplied at the
    5146!   bottom of this module. Users probably want to replace this with something more graceful.
    52 !
     47
    5348module mod_modis_sim
    5449  USE MOD_COSP_TYPES, only: R_UNDEF, ok_debug_cosp
     
    5651  ! ------------------------------
    5752  ! Algorithmic parameters
    58   !
    59  
     53
    6054  real, parameter :: ice_density          = 0.93               ! liquid density is 1. 
    61   !
     55
    6256  ! Retrieval parameters
    63   !
     57
    6458  real, parameter :: min_OpticalThickness = 0.3,             & ! Minimum detectable optical thickness
    6559                     CO2Slicing_PressureLimit = 700. * 100., & ! Cloud with higher pressures use thermal methods, units Pa
     
    7367 
    7468  logical, parameter :: useSimpleReScheme = .false.
    75   !
     69
    7670  ! These are the limits of the libraries for the MODIS collection 5 algorithms
    7771  !   They are also the limits used in the fits for g and w0
    78   !
     72
    7973  real,    parameter :: re_water_min= 4., re_water_max= 30., re_ice_min= 5., re_ice_max= 90.
    8074  integer, parameter :: num_trial_res = 15             ! increase to make the linear pseudo-retrieval of size more accurate
     
    8276!  logical, parameter :: use_two_re_iterations = .false. ! do two retrieval iterations?
    8377! DJS2015 END 
    84   !
     78
    8579  ! Precompute near-IR optical params vs size for retrieval scheme
    86   !
     80
    8781  integer, private :: i
    8882  real, dimension(num_trial_res), parameter :: &
     
    9488  ! ------------------------------
    9589  ! Bin boundaries for the joint optical thickness/cloud top pressure histogram
    96   !
     90
    9791  integer, parameter :: numTauHistogramBins = 6, numPressureHistogramBins = 7
    9892
     
    10397    pressureHistogramBoundaries = (/ 0., 18000., 31000., 44000., 56000., 68000., 80000., 1000000. /)
    10498  real, parameter :: highCloudPressureLimit = 440. * 100., lowCloudPressureLimit = 680. * 100.
    105   !
     99
    106100  ! For output - nominal bin centers and  bin boundaries. On output pressure bins are highest to lowest.
    107   !
     101
    108102  integer, private :: k, l
    109103  real, parameter, dimension(2, numTauHistogramBins) ::   &
     
    155149  !     subroutine modis_L2_simulator_oneTau, or
    156150  !  2) Provide ice and liquid optical depths in each layer
    157   !
     151
    158152  interface modis_L2_simulator
    159153    module procedure modis_L2_simulator_oneTau, modis_L2_simulator_twoTaus
     
    162156  !------------------------------------------------------------------------------------------------
    163157  ! MODIS simulator using specified liquid and ice optical thickness in each layer
    164   !
     158
    165159  !   Note: this simulator operates on all points; to match MODIS itself night-time
    166160  !     points should be excluded
    167   !
     161
    168162  !   Note: the simulator requires as input the optical thickness and cloud top pressure
    169163  !     derived from the ISCCP simulator run with parameter top_height = 1.
     
    171165  !     and reverts to a thermal algorithm much like ISCCP's. Rather than replicate that
    172166  !     alogrithm in this simulator we simply report the values from the ISCCP simulator.
    173   !
     167
    174168  subroutine modis_L2_simulator_twoTaus(                                       &
    175169                                temp, pressureLayers, pressureLevels,          &
     
    212206    nSubcols = size(liquid_opticalThickness, 1)
    213207    nLevels  = size(liquid_opticalThickness, 2)
    214  
    215     !
     208
    216209    ! Initial error checks
    217     !   
     210
    218211    if(any((/ size(ice_opticalThickness, 1), size(waterSize, 1), size(iceSize, 1), &
    219212              size(isccpTau), size(isccpCloudTopPressure),              &
     
    234227             
    235228    ! ---------------------------------------------------
    236     !
     229
    237230    ! Compute the total optical thickness and the proportion due to liquid in each cell
    238     !
     231
    239232    where(liquid_opticalThickness(:, :) + ice_opticalThickness(:, :) > 0.)
    240233      tauLiquidFraction(:, :) = liquid_opticalThickness(:, :)/(liquid_opticalThickness(:, :) + ice_opticalThickness(:, :))
     
    243236    end  where
    244237    tauTotal(:, :) = liquid_opticalThickness(:, :) + ice_opticalThickness(:, :)
    245    
    246     !
     238
    247239    ! Optical depth retrieval
    248240    !   This is simply a sum over the optical thickness in each layer
    249241    !   It should agree with the ISCCP values after min values have been excluded
    250     !
     242
    251243    retrievedTau(:) = sum(tauTotal(:, :), dim = 2)
    252244
    253     !
    254245    ! Cloud detection - does optical thickness exceed detection threshold?
    255     !
     246
    256247    cloudMask = retrievedTau(:) >= min_OpticalThickness
    257    
    258     !
     248
    259249    ! Initialize initial estimates for size retrievals
    260     !
     250
    261251    if(any(cloudMask) .and. .not. useSimpleReScheme) then
    262252      g_w(:)  = get_g_nir(  phaseIsLiquid, trial_re_w(:))
     
    268258    do i = 1, nSubCols
    269259      if(cloudMask(i)) then
    270         !
     260
    271261        ! Cloud top pressure determination
    272262        !   MODIS uses CO2 slicing for clouds with tops above about 700 mb and thermal methods for clouds
     
    276266        ! This assumes linear variation in p between levels. Linear in ln(p) is probably better,
    277267        !   though we'd need to deal with the lowest pressure gracefully.
    278         !
     268
    279269        retrievedCloudTopPressure(i) = cloud_top_pressure((/ 0., tauTotal(i, :) /), &
    280270                                                          pressureLevels,           &
    281271                                                          CO2Slicing_TauLimit) 
    282        
    283        
    284         !
     272
    285273        ! Phase determination - determine fraction of total tau that's liquid
    286274        ! When ice and water contribute about equally to the extinction we can't tell
    287275        !   what the phase is
    288         !
     276
    289277        integratedLiquidFraction = weight_by_extinction(tauTotal(i, :),          &
    290278                                                        tauLiquidFraction(i, :), &
     
    297285          retrievedPhase(i) = phaseIsUndetermined
    298286        end if
    299        
    300         !
     287
    301288        ! Size determination
    302         !
     289
    303290        if(useSimpleReScheme) then
    304291          !   This is the extinction-weighted size considering only the phase we've chosen
    305           !
     292
    306293          if(retrievedPhase(i) == phaseIsIce) then
    307294            retrievedSize(i) = weight_by_extinction(ice_opticalThickness(i, :),  &
     
    323310        end if
    324311      else
    325         !
     312
    326313        ! Values when we don't think there's a cloud.
    327         !
     314
    328315        retrievedCloudTopPressure(i) = R_UNDEF
    329316        retrievedPhase(i) = phaseIsNone
     
    337324    !   mimics what MODIS does to first order.
    338325    !   Of course, ISCCP cloud top pressures are in mb.
    339     !   
     326
    340327    where(cloudMask(:) .and. retrievedCloudTopPressure(:) > CO2Slicing_PressureLimit) &
    341328      retrievedCloudTopPressure(:) = isccpCloudTopPressure * 100.
     
    343330  end subroutine modis_L2_simulator_twoTaus
    344331  !------------------------------------------------------------------------------------------------
    345   !
     332
    346333  ! MODIS simulator: provide a single optical thickness and the cloud ice and liquid contents;
    347334  !   we'll partition this into ice and liquid optical thickness and call the full MODIS simulator
    348   !
     335
    349336  subroutine modis_L2_simulator_oneTau(                                         &
    350337                                temp, pressureLayers, pressureLevels,           &
     
    399386        tauLiquidFraction(:, :) = 0.
    400387      elsewhere
    401         !
     388
    402389        ! Geometic optics limit - tau as LWP/re  (proportional to LWC/re)
    403         !
     390
    404391! Modif AI 02 2018
    405392        tauLiquidFraction(:, :) = (cloudWater(:, :)/max(waterSize(:, :), seuil) ) / &
     
    682669                                        Liquid_Water_Path_Mean,          Ice_Water_Path_Mean,           &   
    683670       Optical_Thickness_vs_Cloud_Top_Pressure)
    684     !
     671
    685672    ! Inputs; dimension nPoints, nSubcols
    686     !
     673
    687674    integer, dimension(:, :),   intent(in)  :: phase
    688675    real,    dimension(:, :),   intent(in)  :: cloud_top_pressure, optical_thickness, particle_size
    689     !
     676
    690677    ! Outputs; dimension nPoints
    691     !
     678
    692679    real,    dimension(:),      intent(out) :: &
    693680       Cloud_Fraction_Total_Mean,       Cloud_Fraction_Water_Mean,       Cloud_Fraction_Ice_Mean,       &
     
    702689    ! ---------------------------
    703690    ! Local variables
    704     !
     691
    705692    real, parameter :: LWP_conversion = 2./3. * 1000. ! MKS units 
    706693    integer :: i, j
     
    714701    nPoints  = size(phase, 1)
    715702    nSubcols = size(phase, 2)
    716     !
     703
    717704    ! Array conformance checks
    718     !
     705
    719706    if(any( (/ size(cloud_top_pressure, 1), size(optical_thickness, 1), size(particle_size, 1),                                &
    720707               size(Cloud_Fraction_Total_Mean),       size(Cloud_Fraction_Water_Mean),       size(Cloud_Fraction_Ice_Mean),    &
     
    728715    if(any( (/ size(cloud_top_pressure, 2), size(optical_thickness, 2), size(particle_size, 2) /)  /= nSubcols)) &
    729716      call complain_and_die("Some L3 arrays have wrong number of subcolumns")
    730    
    731    
    732     !
     717
    733718    ! Include only those pixels with successful retrievals in the statistics
    734     !
     719
    735720    validRetrievalMask(:, :) = particle_size(:, :) > 0.
    736721    cloudMask      = phase(:, :) /= phaseIsNone   .and. validRetrievalMask(:, :)
    737722    waterCloudMask = phase(:, :) == phaseIsLiquid .and. validRetrievalMask(:, :)
    738723    iceCloudMask   = phase(:, :) == phaseIsIce    .and. validRetrievalMask(:, :)
    739     !
     724
    740725    ! Use these as pixel counts at first
    741     !
     726
    742727    Cloud_Fraction_Total_Mean(:) = real(count(cloudMask,      dim = 2))
    743728    Cloud_Fraction_Water_Mean(:) = real(count(waterCloudMask, dim = 2))
     
    747732    Cloud_Fraction_Low_Mean(:)  = real(count(cloudMask .and. cloud_top_pressure >  lowCloudPressureLimit,  dim = 2))
    748733    Cloud_Fraction_Mid_Mean(:)  = Cloud_Fraction_Total_Mean(:) - Cloud_Fraction_High_Mean(:) - Cloud_Fraction_Low_Mean(:)
    749    
    750     !
     734
    751735    ! Don't want to divide by 0, even though the sums will be 0 where the pixel counts are 0.
    752     !
     736
    753737    where (Cloud_Fraction_Total_Mean == 0) Cloud_Fraction_Total_Mean = -1.
    754738    where (Cloud_Fraction_Water_Mean == 0) Cloud_Fraction_Water_Mean = -1.
     
    780764                             / Cloud_Fraction_Ice_Mean(:)
    781765
    782     !
    783766    ! Normalize pixel counts to fraction
    784767    !   The first three cloud fractions have been set to -1 in cloud-free areas, so set those places to 0.
    785     !
     768
    786769    Cloud_Fraction_Total_Mean(:) = max(0., Cloud_Fraction_Total_Mean(:)/nSubcols)
    787770    Cloud_Fraction_Water_Mean(:) = max(0., Cloud_Fraction_Water_Mean(:)/nSubcols)
     
    794777    ! ----
    795778    ! Joint histogram
    796     !
     779
    797780    do i = 1, numTauHistogramBins
    798781      where(cloudMask(:, :))
     
    826809    real,               intent(in) :: tauLimit
    827810    real                           :: cloud_top_pressure
    828     !
     811
    829812    ! Find the extinction-weighted pressure. Assume that pressure varies linearly between
    830813    !   layers and use the trapezoidal rule.
    831     !
    832    
     814
    833815    real :: deltaX, totalTau, totalProduct
    834816    integer :: i
     
    839821        deltaX = tauLimit - totalTau
    840822        totalTau = totalTau + deltaX
    841         !
     823
    842824        ! Result for trapezoidal rule when you take less than a full step
    843825        !   tauIncrement is a layer-integrated value
    844         !
     826
    845827        totalProduct = totalProduct           &
    846828                     + pressure(i-1) * deltaX &
     
    859841    real,               intent(in) :: tauLimit
    860842    real                           :: weight_by_extinction
    861     !
     843
    862844    ! Find the extinction-weighted value of f(tau), assuming constant f within each layer
    863     !
    864    
     845
    865846    real    :: deltaX, totalTau, totalProduct
    866847    integer :: i
     
    892873    ice_g(:)    = get_g_nir(  phaseIsIce,    ice_size)
    893874    ice_w0(:)   = get_ssa_nir(phaseIsIce,    ice_size)
    894     !
     875
    895876    ! Combine ice and water optical properties
    896     !
     877
    897878    g(:) = 0; w0(:) = 0.
    898879    tau(:) = ice_tau(:) + water_tau(:)
     
    913894      real,    intent(in) :: tau, obs_Refl_nir
    914895      real                :: retrieve_re
    915       !
     896
    916897      ! Finds the re that produces the minimum mis-match between predicted and observed reflectance in
    917898      !   MODIS band 7 (near IR)
     
    920901      !  two-stream for layer reflectance and transmittance as a function of optical thickness tau, g, and w0
    921902      !  adding-doubling for total reflectance
    922       ! 
    923       !
    924       !
     903
     904
     905
    925906      ! Local variables
    926       !
     907
    927908      real, parameter :: min_distance_to_boundary = 0.01
    928909      real    :: re_min, re_max, delta_re
     
    946927        w0(:)  = w0_i(:)
    947928      end if
    948       !
     929
    949930      ! 1st attempt at index: w/coarse re resolution
    950       !
     931
    951932      predicted_Refl_nir(:) = two_stream_reflectance(tau, g(:), w0(:))
    952933      retrieve_re = interpolate_to_min(trial_re(:), predicted_Refl_nir(:), obs_Refl_nir)
    953       !
     934
    954935      ! If first retrieval works, can try 2nd iteration using greater re resolution
    955       !
     936
    956937! DJS2015: Remove unused piece of code     
    957938!      if(use_two_re_iterations .and. retrieve_re > 0.) then
     
    959940!        re_max = retrieve_re + delta_re
    960941!        delta_re = (re_max - re_min)/real(num_trial_res-1)
    961 
     942
    962943!        trial_re(:) = re_min + delta_re * (/ (i - 1, i = 1, num_trial_res) /)
    963944!        g(:)  = get_g_nir(  phase, trial_re(:))
     
    977958    real,               intent(in) :: yobs
    978959    real                           :: interpolate_to_min
    979     !
     960
    980961    ! Given a set of values of y as y(x), find the value of x that minimizes abs(y - yobs)
    981962    !   y must be monotonic in x
    982     !
     963
    983964    real, dimension(size(x)) :: diff
    984965    integer                  :: nPoints, minDiffLoc, lowerBound, upperBound
     
    1005986   
    1006987    if(diff(lowerBound) * diff(upperBound) < 0) then     
    1007       !
     988
    1008989      ! Interpolate the root position linearly if we bracket the root
    1009       !
     990
    1010991      interpolate_to_min = x(upperBound) - &
    1011992                           diff(upperBound) * (x(upperBound) - x(lowerBound)) / (diff(upperBound) - diff(lowerBound))
     
    10201001  ! --------------------------------------------
    10211002  elemental function get_g_nir (phase, re)
    1022     !
     1003
    10231004    ! Polynomial fit for asummetry parameter g in MODIS band 7 (near IR) as a function
    10241005    !   of size for ice and water
    10251006    ! Fits from Steve Platnick
    1026     !
    10271007
    10281008    integer, intent(in) :: phase
     
    10561036        real,    intent(in) :: re
    10571037        real                :: get_ssa_nir
    1058         !
     1038
    10591039        ! Polynomial fit for single scattering albedo in MODIS band 7 (near IR) as a function
    10601040        !   of size for ice and water
    10611041        ! Fits from Steve Platnick
    1062         !
     1042
    10631043        real, dimension(4), parameter :: ice_coefficients   = (/ 0.9625, -1.8069e-2, 3.3281e-4,-2.2865e-6/)
    10641044        real, dimension(3), parameter :: water_coefficients = (/ 1.0044, -1.1397e-2, 1.3300e-4 /)
     
    11071087    integer                               :: i
    11081088    ! ---------------------------------------
    1109     !
     1089
    11101090    ! This wrapper reports reflectance only and strips out non-cloudy elements from the calculation
    1111     !
     1091
    11121092    cloudMask = tau(:) > 0.
    11131093    cloudIndicies = pack((/ (i, i = 1, size(tau)) /), mask = cloudMask)
     
    11251105    real, intent(in)  :: tauint, gint, w0int
    11261106    real, intent(out) :: ref, tra
    1127     !
     1107
    11281108    ! Compute reflectance in a single layer using the two stream approximation
    11291109    !   The code itself is from Lazaros Oreopoulos via Steve Platnick
    1130     !
     1110
    11311111    ! ------------------------
    11321112    ! Local variables
     
    11371117    real :: tau, w0, g, f, gamma1, gamma2, gamma3, gamma4, &
    11381118            rh, a1, a2, rk, r1, r2, r3, r4, r5, t1, t2, t3, t4, t5, beta, e1, e2, ef1, ef2, den, th
    1139     !
     1119
    11401120    ! Compute reflectance and transmittance in a single layer using the two stream approximation
    11411121    !   The code itself is from Lazaros Oreopoulos via Steve Platnick
    1142     !
     1122
    11431123    f   = gint**2
    11441124    tau = (1 - w0int * f) * tauint
     
    12041184    real, intent(in) :: tauint, gint, w0int
    12051185    real             :: two_stream_reflectance
    1206     !
     1186
    12071187    ! Compute reflectance in a single layer using the two stream approximation
    12081188    !   The code itself is from Lazaros Oreopoulos via Steve Platnick
    1209     !
     1189
    12101190    ! ------------------------
    12111191    ! Local variables
     
    12791259      real,    dimension(:), intent(in)  :: Refl,     Tran
    12801260      real,                  intent(out) :: Refl_tot, Tran_tot
    1281       !
     1261
    12821262      ! Use adding/doubling formulas to compute total reflectance and transmittance from layer values
    1283       !
    1284      
     1263
    12851264      integer :: i
    12861265      real, dimension(size(Refl)) :: Refl_cumulative, Tran_cumulative
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/optics_lib.F90

    r5095 r5099  
    44! Compiled/Modified:
    55!   07/01/06  John Haynes (haynes@atmos.colostate.edu)
    6 !
     6
    77! m_wat (subroutine)
    88! m_ice (subroutine)
     
    1919  subroutine m_wat(freq, tk, n_r, n_i)
    2020  implicit none
    21 
     21
    2222! Purpose:
    2323!   compute complex index of refraction of liquid water
    24 !
     24
    2525! Inputs:
    2626!   [freq]    frequency (GHz)
    2727!   [tk]       temperature (K)
    28 !
     28
    2929! Outputs:
    3030!   [n_r]     real part index of refraction
    3131!   [n_i]     imaginary part index of refraction
    32 !
     32
    3333! Reference:
    3434!   Based on the work of Ray (1972)
    35 !
     35
    3636! Coded:
    3737!   03/22/05  John Haynes (haynes@atmos.colostate.edu)
     
    8383  subroutine m_ice(freq,t,n_r,n_i)
    8484  implicit none
    85 !
     85
    8686! Purpose:
    8787!   compute complex index of refraction of ice
    88 !
     88
    8989! Inputs:
    9090!   [freq]    frequency (GHz)
    9191!   [t]       temperature (K)
    92 !
     92
    9393! Outputs:
    9494!   [n_r]     real part index of refraction
    9595!   [n_i]     imaginary part index of refraction
    96 !
     96
    9797! Reference:
    9898!    Fortran 90 port from IDL of REFICE by Stephen G. Warren
    99 !
     99
    100100! Modified:
    101101!   05/31/05  John Haynes (haynes@atmos.colostate.edu)
     
    121121! Allowable wavelength range extends from 0.045 microns to 8.6 meter
    122122! temperature dependence only considered beyond 167 microns.
    123 !
     123
    124124! interpolation is done     n_r  vs. log(xlam)
    125125!                           n_r  vs.        t
    126126!                       log(n_i) vs. log(xlam)
    127127!                       log(n_i) vs.        t
    128 !
     128
    129129! Stephen G. Warren - 1983
    130130! Dept. of Atmospheric Sciences
    131131! University of Washington
    132132! Seattle, Wa  98195
    133 !
     133
    134134! Based on
    135 !
     135
    136136!    Warren,S.G.,1984.
    137137!    Optical constants of ice from the ultraviolet to the microwave.
    138138!    Applied Optics,23,1206-1225
    139 !
     139
    140140! Reference temperatures are -1.0,-5.0,-20.0, and -60.0 deg C
    141141 
     
    576576! subroutine MIEINT
    577577! ----------------------------------------------------------------------------
    578 !
     578
    579579!     General purpose Mie scattering routine for single particles
    580580!     Author: R Grainger 1990
     
    583583!     code to ensure correct calculation of backscatter coeficient
    584584!     Options/Extend_Source
    585 !
     585
    586586      Subroutine MieInt(Dx, SCm, Inp, Dqv, Dqxt, Dqsc, Dbsc, Dg, Xs1, Xs2, DPh, Error)
    587587
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/pf_to_mr.F

    r5095 r5099  
    33! $Revision: 23 $, $Date: 2011-03-31 15:41:37 +0200 (jeu. 31 mars 2011) $
    44! $URL: http://cfmip-obs-sim.googlecode.com/svn/stable/v1.4.0/llnl/pf_to_mr.f $
    5 !
     5
    66! Redistribution and use in source and binary forms, with or without modification, are permitted
    77! provided that the following conditions are met:
    8 !
     8
    99!     * Redistributions of source code must retain the above copyright notice, this list
    1010!       of conditions and the following disclaimer.
     
    1515!       nor the names of its contributors may be used to endorse or promote products derived from
    1616!       this software without specific prior written permission.
    17 !
     17
    1818! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR
    1919! IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/prec_scops.F

    r5095 r5099  
    33! $Revision: 23 $, $Date: 2011-03-31 15:41:37 +0200 (jeu. 31 mars 2011) $
    44! $URL: http://cfmip-obs-sim.googlecode.com/svn/stable/v1.4.0/llnl/prec_scops.f $
    5 !
     5
    66! Redistribution and use in source and binary forms, with or without modification, are permitted
    77! provided that the following conditions are met:
    8 !
     8
    99!     * Redistributions of source code must retain the above copyright notice, this list
    1010!       of conditions and the following disclaimer.
     
    1515!       nor the names of its contributors may be used to endorse or promote products derived from
    1616!       this software without specific prior written permission.
    17 !
     17
    1818! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR
    1919! IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/radar_simulator.F90

    r5095 r5099  
    1717 
    1818! Purpose:
    19 !
     19
    2020!   Simulates a vertical profile of radar reflectivity
    2121!   Originally Part of QuickBeam v1.04 by John Haynes & Roger Marchand.
    2222!   but has been substantially modified since that time by
    2323!   Laura Fowler and Roger Marchand (see modifications below).
    24 !
     24
    2525! Inputs:
    26 !
     26
    2727!   [hp]              structure that defines hydrometeor types and other radar properties
    28 !
     28
    2929!   [nprof]           number of hydrometeor profiles
    3030!   [ngate]           number of vertical layers
    31 !
     31
    3232!   [undef]           missing data value
    3333!   (The following 5 arrays must be in order from closest to the radar
    3434!    to farthest...)
    35 !
     35
    3636!   [hgt_matrix]      height of hydrometeors (km)
    3737!   [p_matrix]        pressure profile (hPa)
    3838!   [t_matrix]        temperature profile (K)
    3939!   [rh_matrix]       relative humidity profile (%) -- only needed if gaseous aborption calculated.
    40 !
     40
    4141!   [hm_matrix]       table of hydrometeor mixing rations (g/kg)
    4242!   [re_matrix]       table of hydrometeor effective radii.  0 ==> use defaults. (units=microns)   
    4343!   [Np_matrix]       table of hydrometeor number concentration.  0 ==> use defaults. (units = 1/kg)
    44 !
     44
    4545! Outputs:
    46 !
     46
    4747!   [Ze_non]          radar reflectivity without attenuation (dBZ)
    4848!   [Ze_ray]          Rayleigh reflectivity (dBZ)
     
    5050!   [g_atten_to_vol]  gaseous atteunation, radar to vol (dB)
    5151!   [dBZe]            effective radar reflectivity factor (dBZ)
    52 !
     52
    5353! Optional:
    5454!   [g_to_vol_in]     integrated atten due to gases, r>v (dB).
     
    5757!   [g_to_vol_out]    integrated atten due to gases, r>v (dB).
    5858!                     If present then gaseous absorption for each profile is returned here.
    59 !
     59
    6060! Created:
    6161!   11/28/2005  John Haynes (haynes@atmos.colostate.edu)
    62 !
     62
    6363! Modified:
    6464!   09/2006  placed into subroutine form (Roger Marchand,JMH)
     
    6666!   01/2008  'Do while' to determine if hydrometeor(s) present in volume
    6767!             changed for vectorization purposes (A. Bodas-Salcedo)
    68 !
     68
    6969!   07/2010  V3.0 ... Modified to load or save scale factors to disk as a Look-Up Table (LUT)
    7070!  ... All hydrometeor and radar simulator properties now included in hp structure
     
    7373!     Also ... Support of Morrison 2-moment style microphyscis (Np_matrix) added
    7474!  ... Changes implement by Roj Marchand following work by Laura Fowler
    75 !
     75
    7676!   10/2011  Modified ngate loop to go in either direction depending on flag
    7777!     hp%radar_at_layer_one.  This affects the direction in which attenuation is summed.
    78 !
     78
    7979!     Also removed called to AVINT for gas and hydrometeor attenuation and replaced with simple
    8080!     summation. (Roger Marchand)
    81 !
    82 !
     81
     82
    8383! ----- INPUTS ----- 
    8484 
     
    153153  g_to_vol_out_present = present(g_to_vol_out)
    154154
    155   !
    156155  ! load scaling matricies from disk -- but only the first time this subroutine is called
    157   !
     156
    158157  if(hp%load_scale_LUTs) then
    159158    call load_scale_LUTs(hp)
     
    223222          if(Re.gt.0) then
    224223            ! determine index in to scale LUT
    225             !
     224
    226225            ! distance between Re points (defined by "base" and "step") for
    227226            ! each interval of size Re_BIN_LENGTH
     
    278277                ! alternative is to comment out above two lines and use the following block
    279278                ! MG Mie approach - adjust density of sphere with D = D_characteristic to match particle density
    280                 !
     279
    281280                ! hp%rho_eff(tp,1:ns,iRe_type) = (6/pi)*hp%apm(tp)*(Di*1E-6)**(hp%bpm(tp)-3)   !MG Mie approach
    282281
     
    373372
    374373      !     :: attenuation due to hydrometeors between radar and volume
    375       !
     374
    376375      ! NOTE old scheme integrates attenuation only for the layers ABOVE
    377376      ! the current layer ... i.e. 1 to k-1 rather than 1 to k ...
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/radar_simulator_init.F90

    r5095 r5099  
    1212 
    1313! Purpose:
    14 !
     14
    1515!   Initialize variables used by the radar simulator.
    1616!   Part of QuickBeam v3.0 by John Haynes and Roj Marchand
    17 !   
    18 !
     17
     18
    1919! Inputs: 
    2020!   []   from data in hydrometeor class input
    21 !
     21
    2222!   [freq]            radar frequency (GHz)
    23 !
     23
    2424!   [k2]              |K|^2, the dielectric constant, set to -1 to use the
    2525!                     frequency dependent default
    26 !
     26
    2727!   [use_gas_abs]     1=do gaseous abs calcs, 0=no gasesous absorbtion calculated,
    2828!                     2=calculate (and use) absorption for first profile on all profiles
    29 !
     29
    3030!   [undef]           mising data value
    3131!   [nhclass]         number of hydrometeor types
    32 !
     32
    3333!   For each hydrometero type:
    3434!       hclass_type     Type of distribution (see quickbeam documentation)
    3535!       hclass_phase            1==ice, 0=liquid
    36 !
     36
    3737!       hclass_dmin         minimum diameter allowed is drop size distribution N(D<Dmin)=0
    3838!       hclass_dmax         maximum diameter allowed is drop size distribution N(D>Dmax)=0
    39 !
     39
    4040!   hclass_apm,hclass_bpm   Density of partical apm*D^bpm or constant = rho
    4141!   hclass_rho,
    42 !   
     42
    4343!       hclass_p1,hclass_p2,    Default values of DSD parameters (see quickbeam documentation)
    4444!   hclass_p3,
    45 !
     45
    4646!   load_scale_LUTs_flag    Flag, load scale factors Look Up Table from file at start of run
    4747!   update_scale_LUTs_flag  Flag, save new scale factors calculated during this run to LUT
    4848!   LUT_file_name       Name of file containing LUT
    49 !
     49
    5050! Outputs:
    5151!   [hp]            structure that define hydrometeor types
    52 !
     52
    5353! Modified:
    5454!   08/23/2006  placed into subroutine form (Roger Marchand)
     
    7575  integer :: i,j
    7676  real*8  :: delt, deltp
    77        
    78     !
     77
    7978    ! set radar simulation properites
    80     !
     79
    8180    hp%freq=freq
    8281    hp%k2=k2
     
    8887    hp%update_scale_LUTs=update_scale_LUTs_flag
    8988    hp%scale_LUT_file_name=LUT_file_name
    90    
    91     !
     89
    9290        ! Store settings for hydrometeor types in hp (class_parameter) structure.   
    93         !
     91
    9492        do i = 1,nhclass
    9593        hp%dtype(i) = hclass_type(i)
     
    104102        hp%p3(i)    = hclass_p3(i)
    105103        enddo
    106    
    107         !
     104
    108105        ! initialize scaling array
    109         !
     106
    110107        hp%N_scale_flag = .false.
    111108        hp%fc = undef
     
    116113        hp%Zr_scaled = 0.0
    117114        hp%kr_scaled = 0.0
    118  
    119         !
     115
    120116    ! set up Re bin "structure" for z_scaling
    121         !
     117
    122118    hp%base_list(1)=0;
    123119    do j=1,Re_MAX_BIN
     
    131127    enddo
    132128
    133     !
    134129    ! set up Temperature bin structure used for z scaling
    135     !
     130
    136131    do i=1,cnt_ice
    137132        hp%mt_tti(i)=(i-1)*5-90 + 273.15
     
    141136        hp%mt_ttl(i)=(i-1)*5-60 + 273.15
    142137    enddo
    143  
    144     !
     138
    145139        ! define array discrete diameters used in mie calculations
    146         !
     140
    147141        delt = (log(dmax)-log(dmin))/(nd-1)
    148142        deltp = exp(delt)
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/scale_luts_io.F90

    r3233 r5099  
    11  ! scale_LUT_io:  Contains subroutines to load and save scaling Look Up Tables (LUTs) to a file
    2   !
     2
    33  ! June 2010   Written by Roj Marchand
    44 
     
    1616    logical :: LUT_file_exists
    1717    integer :: i,j,k,ind
    18    
    19     !
     18
    2019    ! load scale LUT from file
    21     !
     20
    2221    inquire(file=trim(hp%scale_LUT_file_name) // '_radar_Z_scale_LUT.dat', &
    2322        exist=LUT_file_exists)
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/scops.F

    r5095 r5099  
    88! $Revision: 23 $, $Date: 2011-03-31 15:41:37 +0200 (jeu. 31 mars 2011) $
    99! $URL: http://cfmip-obs-sim.googlecode.com/svn/stable/v1.4.0/icarus-scops-4.1-bsd/scops.f $
    10 !
     10
    1111! Redistribution and use in source and binary forms, with or without
    1212! modification, are permitted provided that the
    1313! following conditions are met:
    14 !
     14
    1515!     * Redistributions of source code must retain the above
    1616!       copyright  notice, this list of conditions and the following
     
    2424!       derived from this software without specific prior written
    2525!       permission.
    26 !
     26
    2727! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
    2828! "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
     
    3636! (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
    3737! OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 
    38 !
     38
    3939! *****************************COPYRIGHT*******************************
    4040! *****************************COPYRIGHT*******************************
  • LMDZ6/branches/Amaury_dev/libf/phylmd/cosp/zeff.F90

    r5095 r5099  
    1010!   Part of QuickBeam v1.03 by John Haynes
    1111!   http://reef.atmos.colostate.edu/haynes/radarsim
    12 !
     12
    1313! Inputs:
    1414!   [freq]      radar frequency (GHz)
     
    2323!   [qs]        ... efficiencies; otherwise set to -1
    2424!   [rho_e]     medium effective density (kg m^-3) (-1 = pure)
    25 !
     25
    2626! Outputs:
    2727!   [z_eff]     unattenuated effective reflectivity factor (mm^6/m^3)
    2828!   [z_ray]     reflectivity factor, Rayleigh only (mm^6/m^3)
    2929!   [kr]        attenuation coefficient (db km^-1)
    30 !
     30
    3131! Created:
    3232!   11/28/05  John Haynes (haynes@atmos.colostate.edu)
Note: See TracChangeset for help on using the changeset viewer.