source: CONFIG/publications/ICOLMDZORINCA_CO2_Transport_GMD_2023/INCA/build/ppsrc/INCA_SRC/diurnal_geom.f90 @ 6610

Last change on this file since 6610 was 6610, checked in by acosce, 10 months ago

INCA used for ICOLMDZORINCA_CO2_Transport_GMD_2023

File size: 6.9 KB
Line 
1
2
3
4
5
6
7
8
9
10
11
12!$Id: diurnal_geom.F90 123 2009-03-27 10:38:52Z acosce $
13!! =========================================================================
14!! INCA - INteraction with Chemistry and Aerosols
15!!
16!! Copyright Laboratoire des Sciences du Climat et de l'Environnement (LSCE)
17!!           Unite mixte CEA-CNRS-UVSQ
18!!
19!! Contributors to this INCA subroutine:
20!!
21!! Stacy Walters, NCAR, stacy@ucar.edu
22!!
23!! Anne Cozic, LSCE, anne.cozic@cea.fr
24!! Yann Meurdesoif, LSCE, yann.meurdesoif@cea.fr
25!!
26!! This software is a computer program whose purpose is to simulate the
27!! atmospheric gas phase and aerosol composition. The model is designed to be
28!! used within a transport model or a general circulation model. This version
29!! of INCA was designed to be coupled to the LMDz GCM. LMDz-INCA accounts
30!! for emissions, transport (resolved and sub-grid scale), photochemical
31!! transformations, and scavenging (dry deposition and washout) of chemical
32!! species and aerosols interactively in the GCM. Several versions of the INCA
33!! model are currently used depending on the envisaged applications with the
34!! chemistry-climate model.
35!!
36!! This software is governed by the CeCILL  license under French law and
37!! abiding by the rules of distribution of free software.  You can  use,
38!! modify and/ or redistribute the software under the terms of the CeCILL
39!! license as circulated by CEA, CNRS and INRIA at the following URL
40!! "http://www.cecill.info".
41!!
42!! As a counterpart to the access to the source code and  rights to copy,
43!! modify and redistribute granted by the license, users are provided only
44!! with a limited warranty  and the software's author,  the holder of the
45!! economic rights,  and the successive licensors  have only  limited
46!! liability.
47!!
48!! In this respect, the user's attention is drawn to the risks associated
49!! with loading,  using,  modifying and/or developing or reproducing the
50!! software by the user in light of its specific status of free software,
51!! that may mean  that it is complicated to manipulate,  and  that  also
52!! therefore means  that it is reserved for developers  and  experienced
53!! professionals having in-depth computer knowledge. Users are therefore
54!! encouraged to load and test the software's suitability as regards their
55!! requirements in conditions enabling the security of their systems and/or
56!! data to be ensured and,  more generally, to use and operate it in the
57!! same conditions as regards security.
58!!
59!! The fact that you are presently reading this means that you have had
60!! knowledge of the CeCILL license and that you accept its terms.
61!! =========================================================================
62
63
64SUBROUTINE DIURNAL_GEOM(  &
65   lat                      ,&
66   time_of_year             ,&
67   polar_night              ,&
68   polar_day                ,&
69   sunon                    ,&
70   sunoff                   ,&
71   loc_angle                ,& 
72   zen_angle                ,&
73   zangtz )
74 
75  !    Stacy Walters, NCAR, 1998.
76
77  USE CONST_MOD
78  USE CHEM_CONS, ONLY : dayspy, d2r, phi, lambda
79  USE INCA_DIM
80 
81  IMPLICIT NONE
82
83  !------------------------------------------------------------------
84  !     ... Dummy arguments
85  !------------------------------------------------------------------
86  INTEGER, INTENT(in)  ::     lat                ! latitude index
87  REAL, INTENT(in)     ::     time_of_year       ! time of year
88  REAL, INTENT(out)    ::     sunon(PLON)        ! sunrise angle in radians
89  REAL, INTENT(out)    ::     sunoff(PLON)       ! sunset angle in radians
90  REAL, INTENT(out)    ::     zen_angle(PLON)    ! solar zenith angle
91  REAL, INTENT(out)    ::     loc_angle(PLON)    ! "local" time angle
92  LOGICAL, INTENT(out) ::     polar_day(PLON)    ! continuous daylight flag
93  LOGICAL, INTENT(out) ::     polar_night(PLON)  ! continuous night flag
94  LOGICAL, INTENT(out) ::     zangtz(PLON)
95 
96  !------------------------------------------------------------------
97  !        ... Local variables
98  !------------------------------------------------------------------
99  INTEGER :: i, k
100  REAL, SAVE    :: twopi, pid2, dec_max
101!$OMP THREADPRIVATE(twopi, pid2, dec_max)
102  REAL    ::  declination
103  REAL    ::  latitude, longitude
104  REAL    ::  doy                ! day of year
105  REAL    ::  tod                ! time of day
106  REAL    ::  sin_dec, cos_dec   ! sin, cos declination
107  REAL    ::  cosphi             ! cos latitude
108  REAL    ::  sinphi             ! sin latitude
109  LOGICAL, SAVE :: entered = .FALSE.
110!$OMP THREADPRIVATE(entered)
111
112
113  sunon(:) = 0. 
114  sunoff(:) = 0. 
115 
116
117  IF( .NOT. entered ) THEN
118      twopi = 2. * pi
119      pid2  = .5 * pi
120      dec_max = 23.45 * d2r
121      entered = .TRUE.
122  END IF
123
124  !------------------------------------------------------------------
125  !        Note: this formula assumes a 365 day year (inconsistent
126  !        with 360 in LMDz!!!)
127  !------------------------------------------------------------------
128  doy = AINT( time_of_year )
129  declination = dec_max * COS((doy - 172.)*twopi/dayspy)
130
131  !------------------------------------------------------------------
132  !     ... Compute base for zenith angle
133  !------------------------------------------------------------------
134  tod = (time_of_year - doy) + .5
135
136  DO k = 1, PLON
137
138    latitude  = phi(k)
139    longitude = lambda(k)
140    sinphi = SIN( latitude )
141    cosphi = COS( latitude )
142    polar_day(k) = .FALSE.
143    polar_night(k) = .FALSE.
144    !------------------------------------------------------------------
145    !        Determine if in polar day or night
146    !        If NOT in polar day or night then
147    !        calculate terminator longitudes
148    !------------------------------------------------------------------
149    IF( ABS(latitude) >= (pid2 - ABS(declination)) ) THEN
150        IF( SIGN(1.,declination) == SIGN(1.,latitude) ) THEN
151            polar_day(k) = .TRUE.
152            sunoff(k) = 2.*twopi
153            sunon(k)  = -twopi
154        ELSE
155            polar_night(k) = .TRUE.
156            zen_angle(k) = -1.0
157            CYCLE
158        END IF
159    ELSE
160        sunoff(k) = ACOS( -TAN(declination)*TAN(latitude) )
161        sunon(k)  = twopi - sunoff(k)
162    END IF
163
164    sin_dec = SIN( declination )
165    cos_dec = COS( declination )
166    !-------------------------------------------------------------------
167    !        Note: Longitude 0 (Greenwich) at 0:00 hrs
168    !              maps to local angle = pi
169    !-------------------------------------------------------------------
170    loc_angle(k) = (tod + longitude/twopi)*twopi
171    loc_angle(k) = MOD( loc_angle(k)+twopi,twopi )
172   
173    IF( polar_day(k) ) THEN
174        zen_angle(k) = ACOS( sinphi*sin_dec + & 
175           cosphi*cos_dec*COS(loc_angle(k)) )
176    ELSE
177        IF ( loc_angle(k) <= sunoff(k) .OR. &
178           loc_angle(k) >= sunon(k) ) THEN
179            zen_angle(k) = ACOS( sinphi*sin_dec + & 
180               cosphi*cos_dec*COS(loc_angle(k)) )
181        ELSE
182            zen_angle(k) = -1.0
183        ENDIF
184    END IF
185   
186  END DO
187
188  zangtz(:) = zen_angle(:) > 0.
189
190END SUBROUTINE DIURNAL_GEOM
191
Note: See TracBrowser for help on using the repository browser.