FMS 2025.01-dev
Flexible Modeling System
Loading...
Searching...
No Matches
astronomy.inc
1!***********************************************************************
2!* GNU Lesser General Public License
3!*
4!* This file is part of the GFDL Flexible Modeling System (FMS).
5!*
6!* FMS is free software: you can redistribute it and/or modify it under
7!* the terms of the GNU Lesser General Public License as published by
8!* the Free Software Foundation, either version 3 of the License, or (at
9!* your option) any later version.
10!*
11!* FMS is distributed in the hope that it will be useful, but WITHOUT
12!* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
13!* FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
14!* for more details.
15!*
16!* You should have received a copy of the GNU Lesser General Public
17!* License along with FMS. If not, see <http://www.gnu.org/licenses/>.
18!***********************************************************************
19!> @defgroup astronomy_mod astronomy_mod
20!> @ingroup astronomy
21!> @brief Provides astronomical variables for use
22!! by other modules within fms. The only currently used interface is
23!! for determination of astronomical values needed by the shortwave
24!! radiation packages.
25!> @author Fei Liu
26
27!> @addtogroup astronomy_mod
28!> @{
29
30!> @brief set_orbital_parameters saves the input values of eccentricity,
31!! obliquity and perihelion time as module variables for use by
32!! astronomy_mod.
33!!
34!! @throw FATAL, "astronomy_mod module has not been initialized"
35!! @throw FATAL, "astronomy_mod ecc must be between 0 and 0.99"
36!! @throw FATAL, "astronomy_mod obliquity must be between -90. and 90. degrees"
37!! @throw FATAL, "astronomy_mod perihelion must be between 0.0 and 360. degrees"
38
39subroutine set_orbital_parameters_(ecc_in, obliq_in, per_in)
40
41real(kind=fms_ast_kind_), intent(in) :: ecc_in !< Eccentricity of orbital ellipse [dimensionless]
42real(kind=fms_ast_kind_), intent(in) :: obliq_in !< Obliquity [degrees]
43real(kind=fms_ast_kind_), intent(in) :: per_in !< Longitude of perihelion with respect to autumnal
44 !! equinox in northern hemisphere [degrees]
45
46integer, parameter :: lkind = fms_ast_kind_
47
48!---------------------------------------------------------------------
49! exit if module has not been initialized.
50!---------------------------------------------------------------------
51 if (.not. module_is_initialized) &
52 call error_mesg('astronomy_mod', 'module has not been initialized', fatal)
53
54!--------------------------------------------------------------------
55! be sure input values are within valid ranges.
56! QUESTION : ARE THESE THE RIGHT LIMITS ???
57!---------------------------------------------------------------------
58
59 if (ecc_in < 0.0_lkind .or. ecc_in > 0.99_lkind) &
60 call error_mesg('astronomy_mod', 'ecc must be between 0 and 0.99', fatal)
61 if (obliq_in < -90.0_lkind .or. real(obliq, fms_ast_kind_) > 90.0_lkind) &
62 call error_mesg('astronomy_mod', 'obliquity must be between -90. and 90. degrees', fatal)
63 if (per_in < 0.0_lkind .or. per_in > 360.0_lkind) &
64 call error_mesg('astronomy_mod', 'perihelion must be between 0.0 and 360. degrees', fatal)
65
66!---------------------------------------------------------------------
67! save input values into module variables.
68!---------------------------------------------------------------------
69
70 ecc = real(ecc_in, r8_kind)
71 obliq = real(obliq_in, r8_kind)
72 per = real(per_in, r8_kind)
73
74!---------------------------------------------------------------------
75! call orbit to define table of orbital angles as function of
76! orbital time using the input values of parameters just supplied.
77!----------------------------------------------------------------------
78
79 call orbit
80
81!----------------------------------------------------------------------
82
83end subroutine set_orbital_parameters_
84
85!> @brief get_orbital_parameters retrieves the orbital parameters for use
86!! by another module.
87!!
88!! @throw FATAL, "astronomy_mod module has not been initialized"
89
90subroutine get_orbital_parameters_(ecc_out, obliq_out, per_out)
91
92!-------------------------------------------------------------------
93! get_orbital_parameters retrieves the orbital parameters for use
94! by another module.
95!--------------------------------------------------------------------
96
97
98real(kind=fms_ast_kind_), intent(out) :: ecc_out !< Eccentricity of orbital ellipse [dimensionless]
99real(kind=fms_ast_kind_), intent(out) :: obliq_out !< Obliquity [degrees]
100real(kind=fms_ast_kind_), intent(out) :: per_out !< Longitude of perihelion with respect to autumnal
101 !! equinox in northern hemisphere [degrees]
102
103!---------------------------------------------------------------------
104! exit if module has not been initialized.
105!---------------------------------------------------------------------
106
107 if (.not. module_is_initialized) &
108 call error_mesg('astronomy_mod', 'module has not been initialized', fatal)
109
110!--------------------------------------------------------------------
111! fill the output arguments with the eccentricity, obliquity and
112! perihelion angle.
113!--------------------------------------------------------------------
114
115 ecc_out = real(ecc, fms_ast_kind_)
116 obliq_out = real(obliq, fms_ast_kind_)
117 per_out = real(per, fms_ast_kind_)
118
119end subroutine get_orbital_parameters_
120
121
122!> @brief diurnal_solar_2d returns 2d fields of cosine of zenith angle,
123!! daylight fraction and earth-sun distance at the specified latitudes,
124!! longitudes and time. These values may be instantaneous or averaged
125!! over a specified time interval.
126!!
127!! @param [in] <lat> Latitudes of model grid points
128!! @param [in] <lon> Longitudes of model grid points
129!! @param [in] <gmt> Time of day at longitude 0.0; midnight = 0.0, one day = 2 * pi
130!! @param [in] <time_since_ae> Time of year; autumnal equinox = 0.0, one year = 2 * pi
131!! @param [out] <cosz> Cosine of solar zenith angle
132!! @param [out] <fracday> Daylight fraction of time interval
133!! @param [out] <rrsun> earth-sun distance (r) relative to semi-major axis of orbital ellipse (a):(a/r)**2
134!! @param [in] <dt> OPTIONAL: time interval after gmt over which the astronomical variables are to be
135!! averaged. this produces averaged output rather than instantaneous.
136!! @param [in] <allow_negative_cosz> Allow negative values for cosz?
137!! @param [out] <half_day_out> half_day_out
138!!
139!! @throw FATAL, "astronomy_mod time_since_ae not between 0 and 2pi"
140!! @throw FATAL, "astronomy_mod gmt not between 0 and 2pi"
141
142subroutine diurnal_solar_2d_(lat, lon, gmt, time_since_ae, cosz, &
143 fracday, rrsun, dt, allow_negative_cosz, &
144 half_day_out)
145
146real(kind=fms_ast_kind_), dimension(:,:), intent(in) :: lat, lon
147real(kind=fms_ast_kind_), intent(in) :: gmt, time_since_ae
148real(kind=fms_ast_kind_), dimension(:,:), intent(out) :: cosz, fracday
149real(kind=fms_ast_kind_), intent(out) :: rrsun
150real(kind=fms_ast_kind_), intent(in), optional :: dt
151logical, intent(in), optional :: allow_negative_cosz
152real(kind=fms_ast_kind_), dimension(:,:), intent(out), optional :: half_day_out
153
154!---------------------------------------------------------------------
155! local variables
156
157
158real(kind=fms_ast_kind_), dimension(size(lat,1),size(lat,2)) :: t, tt, h, aa, bb, st, stt, sh
159real(kind=fms_ast_kind_) :: ang, dec
160logical :: Lallow_negative
161integer, parameter :: lkind = fms_ast_kind_
162
163
164!---------------------------------------------------------------------
165! local variables
166!
167! t time of day with respect to local noon (2 pi = 1 day)
168! [ radians ]
169! tt end of averaging period [ radians ]
170! h half of the daily period of daylight, centered at noon
171! [ radians, -pi --> pi ]
172! aa sin(lat) * sin(declination)
173! bb cos(lat) * cos(declination)
174! st sine of time of day
175! stt sine of time of day at end of averaging period
176! sh sine of half-day period
177! ang position of earth in its orbit wrt autumnal equinox
178! [ radians ]
179! dec earth's declination [ radians ]
180!
181!--------------------------------------------------------------------
182
183!--------------------------------------------------------------------
184! be sure the time in the annual cycle is legitimate.
185!---------------------------------------------------------------------
186
187 if (time_since_ae < 0.0_lkind .or. time_since_ae > real(twopi, fms_ast_kind_)) &
188 call error_mesg('astronomy_mod', 'time_since_ae not between 0 and 2pi', fatal)
189
190!--------------------------------------------------------------------
191! be sure the time at longitude = 0.0 is legitimate.
192!---------------------------------------------------------------------
193
194 if (gmt < 0.0_lkind .or. gmt > real(twopi, fms_ast_kind_)) &
195 call error_mesg('astronomy_mod', 'gmt not between 0 and 2pi', fatal)
196
197!---------------------------------------------------------------------
198!> define the orbital angle (location in year), solar declination and
199!! earth sun distance factor. use functions contained in this module.
200!---------------------------------------------------------------------
201
202 ang = angle(time_since_ae)
203 dec = declination(ang)
204 rrsun = r_inv_squared(ang)
205
206!---------------------------------------------------------------------
207!> define terms needed in the cosine zenith angle equation.
208!--------------------------------------------------------------------
209
210 aa = sin(lat)*sin(dec)
211 bb = cos(lat)*cos(dec)
212
213!---------------------------------------------------------------------
214!> define local time. force it to be between -pi and pi.
215!--------------------------------------------------------------------
216
217 t = gmt + lon - real(pi, fms_ast_kind_)
218 where(t >= real(pi, fms_ast_kind_)) t = t - real(twopi, fms_ast_kind_)
219 where(t < real(-pi,fms_ast_kind_)) t = t + real(twopi, fms_ast_kind_)
220
221 lallow_negative = .false.
222 if (present(allow_negative_cosz)) then
223 if (allow_negative_cosz) lallow_negative = .true.
224 end if
225
226
227!---------------------------------------------------------------------
228!> perform a time integration to obtain cosz and fracday if desired.
229!! output is valid over the period from t to t + dt.
230!--------------------------------------------------------------------
231
232 h = half_day(lat,dec)
233
234 if ( present(half_day_out) ) then
235 half_day_out = h
236 end if
237
238 if ( present(dt) ) then ! (perform time averaging)
239 tt = t + dt
240 st = sin(t)
241 stt = sin(tt)
242 sh = sin(h)
243 cosz = 0.0_lkind
244
245 if (.not. lallow_negative) then
246!-------------------------------------------------------------------
247! case 1: entire averaging period is before sunrise.
248!-------------------------------------------------------------------
249
250 where (t < -h .and. tt < -h) cosz = 0.0_lkind
251
252!-------------------------------------------------------------------
253! case 2: averaging period begins before sunrise, ends after sunrise
254! but before sunset
255!-------------------------------------------------------------------
256
257 where ( (tt+h) /= 0.0_lkind .and. t < -h .and. abs(tt) <= h) &
258 cosz = aa + bb*(stt + sh)/ (tt + h)
259
260!-------------------------------------------------------------------
261! case 3: averaging period begins before sunrise, ends after sunset,
262! but before the next sunrise. modify if averaging period extends
263! past the next day's sunrise, but if averaging period is less than
264! a half- day (pi) that circumstance will never occur.
265!-------------------------------------------------------------------
266
267 where (t < -h .and. h /= 0.0_lkind .and. h < tt) &
268 cosz = aa + bb*( sh + sh)/(h+h)
269
270!-------------------------------------------------------------------
271! case 4: averaging period begins after sunrise, ends before sunset.
272!-------------------------------------------------------------------
273 where ( abs(t) <= h .and. abs(tt) <= h) &
274
275 cosz = aa + bb*(stt - st)/ (tt - t)
276
277!-------------------------------------------------------------------
278! case 5: averaging period begins after sunrise, ends after sunset.
279! modify when averaging period extends past the next day's sunrise.
280!-------------------------------------------------------------------
281
282 where ((h-t) /= 0.0_lkind .and. abs(t) <= h .and. h < tt) &
283 cosz = aa + bb*(sh - st)/(h-t)
284
285!-------------------------------------------------------------------
286! case 6: averaging period begins after sunrise , ends after the
287! next day's sunrise. note that this includes the case when the
288! day length is one day (h = pi).
289!-------------------------------------------------------------------
290
291 where (real(twopi, fms_ast_kind_) - h < tt .and. &
292 (tt+h-real(twopi, fms_ast_kind_)) /= 0.0_lkind .and. t <= h ) &
293 cosz = (cosz*(h - t) + (aa*(tt + h - real(twopi, fms_ast_kind_)) + bb &
294 * (stt + sh))) / ((h - t) + (tt + h - real(twopi, fms_ast_kind_)))
295
296!-------------------------------------------------------------------
297! case 7: averaging period begins after sunset and ends before the
298! next day's sunrise
299!-------------------------------------------------------------------
300
301 where(h < t .and. real(twopi, fms_ast_kind_) - h >= tt) &
302 cosz = 0.0_lkind
303
304!-------------------------------------------------------------------
305! case 8: averaging period begins after sunset and ends after the
306! next day's sunrise but before the next day's sunset. if the
307! averaging period is less than a half-day (pi) the latter
308! circumstance will never occur.
309!-----------------------------------------------------------------
310
311 where (h < t .and. real(twopi, fms_ast_kind_) - h < tt)
312 cosz = aa + bb*(stt + sh) / (tt + h - real(twopi, fms_ast_kind_))
313 end where
314
315 else
316 cosz = aa + bb*(stt - st)/ (tt - t)
317 end if
318
319
320
321!-------------------------------------------------------------------
322! day fraction is the fraction of the averaging period contained
323! within the (-h,h) period.
324!-------------------------------------------------------------------
325
326 where (t < -h .and. tt < -h) fracday = 0.0_lkind
327 where (t < -h .and. abs(tt) <= h) fracday = (tt + h )/dt
328 where (t < -h .and. h < tt) fracday = ( h + h )/dt
329 where (abs(t) <= h .and. abs(tt) <= h) fracday = (tt - t )/dt
330 where (abs(t) <= h .and. h < tt) fracday = ( h - t )/dt
331 where (h < t) fracday = 0.0_lkind
332 where (real(twopi, fms_ast_kind_) - h < tt) &
333 fracday = fracday + (tt + h - real(twopi, fms_ast_kind_))/dt
334!----------------------------------------------------------------------
335!> if instantaneous values are desired, define cosz at time t.
336!----------------------------------------------------------------------
337 else ! (no time averaging)
338 if (.not. lallow_negative) then
339 where (abs(t) < h)
340 cosz = aa + bb*cos(t)
341 fracday = 1.0_lkind
342 elsewhere
343 cosz = 0.0_lkind
344 fracday = 0.0_lkind
345 end where
346 else
347 cosz = aa + bb*cos(t)
348 where (abs(t) < h)
349 fracday = 1.0_lkind
350 elsewhere
351 fracday = 0.0_lkind
352 end where
353 end if
354 end if
355
356
357!----------------------------------------------------------------------
358!> Check that cosz is not negative, if desired.
359!----------------------------------------------------------------------
360
361 if (.not. lallow_negative) then
362 cosz = max(0.0_lkind, cosz)
363 end if
364
365end subroutine diurnal_solar_2d_
366
367
368!> @brief diurnal_solar_1d takes 1-d input fields, adds a second dimension
369!! and calls diurnal_solar_2d. on return, the 2d fields are returned
370!! to the original 1d fields.
371!!
372!! @param [in] <lat> Latitudes of model grid points
373!! @param [in] <lon> Longitudes of model grid points
374!! @param [in] <gmt> Time of day at longitude 0.0; midnight = 0.0, one day = 2 * pi
375!! @param [in] <time_since_ae> Time of year; autumnal equinox = 0.0, one year = 2 * pi
376!! @param [out] <cosz> Cosine of solar zenith angle
377!! @param [out] <fracday> Daylight fraction of time interval
378!! @param [out] <rrsun> earth-sun distance (r) relative to semi-major axis of orbital ellipse (a):(a/r)**2
379!! @param [in] <dt> OPTIONAL: time interval after gmt over which the astronomical variables are to be
380!! averaged. this produces averaged output rather than instantaneous.
381!! @param [in] <allow_negative_cosz> Allow negative values for cosz?
382!! @param [out] <half_day_out> half_day_out
383
384subroutine diurnal_solar_1d_(lat, lon, gmt, time_since_ae, cosz, &
385 fracday, rrsun, dt, allow_negative_cosz, &
386 half_day_out)
387
388!---------------------------------------------------------------------
389
390real(kind=fms_ast_kind_), dimension(:), intent(in) :: lat, lon
391real(kind=fms_ast_kind_), intent(in) :: gmt, time_since_ae
392real(kind=fms_ast_kind_), dimension(:), intent(out) :: cosz, fracday
393real(kind=fms_ast_kind_), intent(out) :: rrsun
394real(kind=fms_ast_kind_), intent(in), optional :: dt
395logical, intent(in), optional :: allow_negative_cosz
396real(kind=fms_ast_kind_), dimension(:), intent(out), optional :: half_day_out
397
398!---------------------------------------------------------------------
399! local variables
400!---------------------------------------------------------------------
401
402real(kind=fms_ast_kind_), dimension(size(lat),1) :: lat_2d, lon_2d, cosz_2d, fracday_2d, halfday_2d
403
404!--------------------------------------------------------------------
405!> define 2-d versions of input data arrays.
406!--------------------------------------------------------------------
407
408 lat_2d(:,1) = lat
409 lon_2d(:,1) = lon
410
411
412!--------------------------------------------------------------------
413!> call diurnal_solar_2d to calculate astronomy fields.
414!--------------------------------------------------------------------
415! if (present(dt)) then
416
417 call diurnal_solar(lat_2d, lon_2d, gmt, time_since_ae, &
418 cosz_2d, fracday_2d, rrsun, dt=dt, &
419 allow_negative_cosz=allow_negative_cosz, half_day_out=halfday_2d)
420! else
421! call diurnal_solar_2d (lat_2d, lon_2d, gmt, time_since_ae, &
422! cosz_2d, fracday_2d, rrsun)
423! endif
424
425!-------------------------------------------------------------------
426!> place output fields into 1-d arguments for return to
427!! calling routine.
428!-------------------------------------------------------------------
429
430 fracday = fracday_2d(:,1)
431 cosz = cosz_2d(:,1)
432 if (present(half_day_out)) then
433 half_day_out = halfday_2d(:,1)
434 end if
435
436end subroutine diurnal_solar_1d_
437
438
439
440!> @brief diurnal_solar_0d takes scalar input fields, makes them into 2d
441!! arrays dimensioned (1,1), and calls diurnal_solar_2d. on return,
442!! the 2d fields are converted back to the desired scalar output.
443!!
444!! @param [in] <lat> Latitudes of model grid points
445!! @param [in] <lon> Longitudes of model grid points
446!! @param [in] <gmt> Time of day at longitude 0.0; midnight = 0.0, one day = 2 * pi
447!! @param [in] <time_since_ae> Time of year; autumnal equinox = 0.0, one year = 2 * pi
448!! @param [out] <cosz> Cosine of solar zenith angle
449!! @param [out] <fracday> Daylight fraction of time interval
450!! @param [out] <rrsun> earth-sun distance (r) relative to semi-major axis of orbital ellipse (a):(a/r)**2
451!! @param [in] <dt> OPTIONAL: time interval after gmt over which the astronomical variables are to be
452!! averaged. this produces averaged output rather than instantaneous.
453!! @param [in] <allow_negative_cosz> Allow negative values for cosz?
454!! @param [out] <half_day_out> half_day_out
455
456subroutine diurnal_solar_0d_(lat, lon, gmt, time_since_ae, cosz, &
457 fracday, rrsun, dt, allow_negative_cosz, &
458 half_day_out)
459
460real(kind=fms_ast_kind_), intent(in) :: lat, lon, gmt, time_since_ae
461real(kind=fms_ast_kind_), intent(out) :: cosz, fracday, rrsun
462real(kind=fms_ast_kind_), intent(in), optional :: dt
463logical, intent(in), optional :: allow_negative_cosz
464real(kind=fms_ast_kind_), intent(out), optional :: half_day_out
465
466!--------------------------------------------------------------------
467! local variables:
468!--------------------------------------------------------------------
469
470real(kind=fms_ast_kind_), dimension(1,1) :: lat_2d, lon_2d, cosz_2d, fracday_2d, halfday_2d
471
472!---------------------------------------------------------------------
473!> create 2d arrays from the scalar input fields.
474!---------------------------------------------------------------------
475
476 lat_2d = lat
477 lon_2d = lon
478
479!--------------------------------------------------------------------
480!> call diurnal_solar_2d to calculate astronomy fields.
481!--------------------------------------------------------------------
482! if (present(dt)) then
483
484 call diurnal_solar(lat_2d, lon_2d, gmt, time_since_ae, &
485 cosz_2d, fracday_2d, rrsun, dt=dt, &
486 allow_negative_cosz=allow_negative_cosz, half_day_out=halfday_2d)
487! else
488! call diurnal_solar_2d (lat_2d, lon_2d, gmt, time_since_ae, &
489! cosz_2d, fracday_2d, rrsun)
490! end if
491
492!-------------------------------------------------------------------
493!> place output fields into scalars for return to calling routine.
494!-------------------------------------------------------------------
495
496 fracday = fracday_2d(1,1)
497 cosz = cosz_2d(1,1)
498 if (present(half_day_out)) then
499 half_day_out = halfday_2d(1,1)
500 end if
501
502end subroutine diurnal_solar_0d_
503
504
505!> @brief diurnal_solar_cal_2d receives time_type inputs, converts
506!! them to real variables and then calls diurnal_solar_2d to
507!! compute desired astronomical variables.
508!!
509!! @param [in] <lat> Latitudes of model grid points
510!! @param [in] <lon> Longitudes of model grid points
511!! @param [in] <gmt> Time of day at longitude 0.0; midnight = 0.0, one day = 2 * pi
512!! @param [in] <time> Time of year (time_type)
513!! @param [out] <cosz> Cosine of solar zenith angle
514!! @param [out] <fracday> Daylight fraction of time interval
515!! @param [out] <rrsun> earth-sun distance (r) relative to semi-major axis of orbital ellipse (a):(a/r)**2
516!! @param [in] <dt> OPTIONAL: time interval after gmt over which the astronomical variables are to be
517!! averaged. this produces averaged output rather than instantaneous.
518!! @param [in] <allow_negative_cosz> Allow negative values for cosz?
519!! @param [out] <half_day_out> half_day_out
520!!
521!! @throw FATAL, "astronomy_mod radiation time step must be no longer than 12 hrs"
522!! @throw FATAL, "astronomy_mod radiation time step must not be an integral number of days"
523
524subroutine diurnal_solar_cal_2d_(lat, lon, time, cosz, fracday, &
525 rrsun, dt_time, allow_negative_cosz, &
526 half_day_out)
527
528!-------------------------------------------------------------------
529
530real(kind=fms_ast_kind_), dimension(:,:), intent(in) :: lat, lon
531type(time_type), intent(in) :: time
532real(kind=fms_ast_kind_), dimension(:,:), intent(out) :: cosz, fracday
533real(kind=fms_ast_kind_), intent(out) :: rrsun
534type(time_type), intent(in), optional :: dt_time
535logical, intent(in), optional :: allow_negative_cosz
536real(kind=fms_ast_kind_), dimension(:,:), intent(out), optional :: half_day_out
537
538!---------------------------------------------------------------------
539
540!---------------------------------------------------------------------
541! local variables
542!---------------------------------------------------------------------
543
544real(kind=fms_ast_kind_) :: dt
545real(kind=fms_ast_kind_) :: gmt, time_since_ae
546integer, parameter :: lkind = fms_ast_kind_
547
548!---------------------------------------------------------------------
549!> Extract time of day (gmt) from time_type variable time with
550!! function universal_time.
551!---------------------------------------------------------------------
552
553 gmt = real(universal_time(time), fms_ast_kind_)
554
555!---------------------------------------------------------------------
556!> Extract the time of year (time_since_ae) from time_type variable
557!! time using the function orbital_time.
558!---------------------------------------------------------------------
559
560 time_since_ae = real(orbital_time(time), fms_ast_kind_)
561
562!---------------------------------------------------------------------
563!> Convert optional time_type variable dt_time (length of averaging
564!! period) to a real variable dt with the function universal_time.
565!---------------------------------------------------------------------
566
567 if (present(dt_time)) then
568 dt = real(universal_time(dt_time), fms_ast_kind_)
569 if (dt > real(pi, fms_ast_kind_)) then
570 call error_mesg('astronomy_mod', 'radiation time step must be no longer than 12 hrs', fatal)
571 endif
572 if (dt == 0.0_lkind) then
573 call error_mesg('astronomy_mod', 'radiation time step must not be an integral number of days', fatal)
574 endif
575
576!--------------------------------------------------------------------
577!> Call diurnal_solar_2d to calculate astronomy fields, with or
578!! without the optional argument dt.
579!--------------------------------------------------------------------
580
581 call diurnal_solar(lat, lon, gmt, time_since_ae, cosz, &
582 fracday, rrsun, dt=dt, &
583 allow_negative_cosz=allow_negative_cosz, &
584 half_day_out=half_day_out)
585 else
586 call diurnal_solar(lat, lon, gmt, time_since_ae, cosz, &
587 fracday, rrsun, allow_negative_cosz=allow_negative_cosz, &
588 half_day_out=half_day_out)
589 end if
590
591end subroutine diurnal_solar_cal_2d_
592
593
594!> @brief diurnal_solar_cal_1d receives time_type inputs, converts
595!! them to real variables and then calls diurnal_solar_2d to
596!! compute desired astronomical variables.
597!!
598!! @param [in] <lat> Latitudes of model grid points
599!! @param [in] <lon> Longitudes of model grid points
600!! @param [in] <gmt> Time of day at longitude 0.0; midnight = 0.0, one day = 2 * pi
601!! @param [in] <time> Time of year (time_type)
602!! @param [out] <cosz> Cosine of solar zenith angle
603!! @param [out] <fracday> Daylight fraction of time interval
604!! @param [out] <rrsun> earth-sun distance (r) relative to semi-major axis of orbital ellipse (a):(a/r)**2
605!! @param [in] <dt> OPTIONAL: time interval after gmt over which the astronomical variables are to be
606!! averaged. this produces averaged output rather than instantaneous.
607!! @param [in] <allow_negative_cosz> Allow negative values for cosz?
608!! @param [out] <half_day_out> half_day_out
609
610subroutine diurnal_solar_cal_1d_(lat, lon, time, cosz, fracday, &
611 rrsun, dt_time, allow_negative_cosz, &
612 half_day_out)
613
614!--------------------------------------------------------------------
615
616real(kind=fms_ast_kind_), dimension(:), intent(in) :: lat, lon
617type(time_type), intent(in) :: time
618real(kind=fms_ast_kind_), dimension(:), intent(out) :: cosz, fracday
619real(kind=fms_ast_kind_), intent(out) :: rrsun
620type(time_type), intent(in), optional :: dt_time
621logical, intent(in), optional :: allow_negative_cosz
622real(kind=fms_ast_kind_), dimension(:), intent(out), optional :: half_day_out
623!--------------------------------------------------------------------
624
625!-------------------------------------------------------------------
626! local variables
627!-------------------------------------------------------------------
628
629real(kind=fms_ast_kind_), dimension(size(lat),1) :: lat_2d, lon_2d, cosz_2d, &
630 fracday_2d, halfday_2d
631
632!--------------------------------------------------------------------
633!> Define 2-d versions of input data arrays.
634!--------------------------------------------------------------------
635
636 lat_2d(:,1) = lat
637 lon_2d(:,1) = lon
638
639!--------------------------------------------------------------------
640!> Call diurnal_solar_cal_2d to convert the time_types to reals and
641!! then calculate the astronomy fields.
642!--------------------------------------------------------------------
643
644 if (present(dt_time)) then
645 call diurnal_solar(lat_2d, lon_2d, time, cosz_2d, &
646 fracday_2d, rrsun, dt_time=dt_time, &
647 allow_negative_cosz=allow_negative_cosz, &
648 half_day_out=halfday_2d)
649 else
650 call diurnal_solar(lat_2d, lon_2d, time, cosz_2d, &
651 fracday_2d, rrsun, allow_negative_cosz=allow_negative_cosz, &
652 half_day_out=halfday_2d)
653 end if
654
655!-------------------------------------------------------------------
656!> Place output fields into 1-d arguments for return to
657!! calling routine.
658!-------------------------------------------------------------------
659
660 fracday = fracday_2d(:,1)
661 cosz = cosz_2d(:,1)
662 if (present(half_day_out)) then
663 half_day_out = halfday_2d(:,1)
664 end if
665
666end subroutine diurnal_solar_cal_1d_
667
668
669!> @brief diurnal_solar_cal_0d receives time_type inputs, converts them to real variables
670!! and then calls diurnal_solar_2d to compute desired astronomical variables.
671!!
672!! @param [in] <lat> Latitudes of model grid points
673!! @param [in] <lon> Longitudes of model grid points
674!! @param [in] <time> Time of year (time_type)
675!! @param [out] <cosz> Cosine of solar zenith angle
676!! @param [out] <fracday> Daylight fraction of time interval
677!! @param [out] <rrsun> Earth-Sun distance (r) relative to semi-major axis of orbital ellipse (a) : (a/r)**2
678!! @param [out] <dt_time> OPTIONAL: time interval after gmt over which the astronomical variables are
679!! to be averaged. this produces averaged output rather than instantaneous.
680!! @param [in] <allow_negative_cosz> allow_negative_cosz
681!! @param [out] <half_day_out> half_day_out
682
683subroutine diurnal_solar_cal_0d_(lat, lon, time, cosz, fracday, &
684 rrsun, dt_time, allow_negative_cosz, &
685 half_day_out)
686
687!---------------------------------------------------------------------
688
689real(kind=fms_ast_kind_), intent(in) :: lat, lon
690type(time_type), intent(in) :: time
691real(kind=fms_ast_kind_), intent(out) :: cosz, fracday, rrsun
692type(time_type), intent(in), optional :: dt_time
693logical, intent(in), optional :: allow_negative_cosz
694real(kind=fms_ast_kind_), intent(out), optional :: half_day_out
695
696!---------------------------------------------------------------------
697
698!---------------------------------------------------------------------
699! local variables
700!---------------------------------------------------------------------
701
702real(kind=fms_ast_kind_), dimension(1,1) :: lat_2d, lon_2d, cosz_2d, fracday_2d, halfday_2d
703
704!--------------------------------------------------------------------
705!> Define 2-d versions of input data arrays.
706!--------------------------------------------------------------------
707
708 lat_2d = lat
709 lon_2d = lon
710
711!--------------------------------------------------------------------
712!> Call diurnal_solar_cal_2d to convert the time_types to reals and
713!! then calculate the astronomy fields.
714!--------------------------------------------------------------------
715
716 if (present(dt_time)) then
717 call diurnal_solar(lat_2d, lon_2d, time, cosz_2d, &
718 fracday_2d, rrsun, dt_time=dt_time, &
719 allow_negative_cosz=allow_negative_cosz, &
720 half_day_out=halfday_2d)
721 else
722 call diurnal_solar(lat_2d, lon_2d, time, cosz_2d, &
723 fracday_2d, rrsun, allow_negative_cosz=allow_negative_cosz, &
724 half_day_out=halfday_2d)
725 end if
726
727!-------------------------------------------------------------------
728!> Place output fields into 1-d arguments for return to
729!! calling routine.
730!-------------------------------------------------------------------
731
732 fracday = fracday_2d(1,1)
733 cosz = cosz_2d(1,1)
734 if (present(half_day_out)) then
735 half_day_out = halfday_2d(1,1)
736 end if
737
738end subroutine diurnal_solar_cal_0d_
739
740
741!> @brief daily_mean_solar_2d computes the daily mean astronomical parameters for
742!! the input points at latitude lat and time of year time_since_ae.
743!!
744!! @param [in] <lat> Latitudes of model grid points
745!! @param [in] <time_since_ae> Time of year; autumnal equinox = 0.0, one year = 2 * pi
746!! @param [out] <cosz> Cosine of solar zenith angle
747!! @param [out] <h_out> 2-d array of half-day lengths at the latitudes
748!! @param [out] <rr_out> the inverse of the square of the earth-sun distance relative
749!! to the mean distance at angle ang in the earth's orbit.
750!!
751!! @throw FATAL, "astronomy_mod time_since_ae not between 0 and 2pi"
752
753subroutine daily_mean_solar_2d_(lat, time_since_ae, cosz, h_out, rr_out)
754
755!----------------------------------------------------------------------
756real(kind=fms_ast_kind_), dimension(:,:), intent(in) :: lat
757real(kind=fms_ast_kind_), intent(in) :: time_since_ae
758real(kind=fms_ast_kind_), dimension(:,:), intent(out) :: cosz, h_out
759real(kind=fms_ast_kind_), intent(out) :: rr_out
760!----------------------------------------------------------------------
761
762!--------------------------------------------------------------------
763! local variables
764!--------------------------------------------------------------------
765
766real(kind=fms_ast_kind_), dimension(size(lat,1),size(lat,2)) :: h
767real(kind=fms_ast_kind_) :: ang, dec, rr
768integer, parameter :: lkind = fms_ast_kind_
769
770!--------------------------------------------------------------------
771! be sure the time in the annual cycle is legitimate.
772!---------------------------------------------------------------------
773
774 if (time_since_ae < 0.0_lkind .or. time_since_ae > real(twopi, fms_ast_kind_)) &
775 call error_mesg('astronomy_mod', 'time_since_ae not between 0 and 2pi', fatal)
776
777!---------------------------------------------------------------------
778!> Define the orbital angle (location in year), solar declination,
779!! half-day length and earth sun distance factor. Use functions
780!! contained in this module.
781!---------------------------------------------------------------------
782
783 ang = angle(time_since_ae)
784 dec = declination(ang)
785 h = half_day(lat, dec)
786 rr = r_inv_squared(ang)
787
788!---------------------------------------------------------------------
789!> Where the entire day is dark, define cosz to be zero. otherwise
790!! use the standard formula. Define the daylight fraction and earth-
791!! sun distance.
792!---------------------------------------------------------------------
793
794 where (h == 0.0_lkind)
795 cosz = 0.0_lkind
796 elsewhere
797 cosz = sin(lat)*sin(dec) + cos(lat)*cos(dec)*sin(h)/h
798 end where
799 h_out = h/real(pi, fms_ast_kind_)
800 rr_out = rr
801
802end subroutine daily_mean_solar_2d_
803
804
805!> @brief daily_mean_solar_1d takes 1-d input fields, adds a second dimension
806!! and calls daily_mean_solar_2d. on return, the 2d fields are
807!! returned to the original 1d fields.
808!!
809!! @param [in] <lat> Latitudes of model grid points
810!! @param [in] <time_since_ae> Time of year; autumnal equinox = 0.0, one year = 2 * pi
811!! @param [out] <cosz> Cosine of solar zenith angle
812!! @param [out] <h_out> 2-d array of half-day lengths at the latitudes
813!! @param [out] <rr_out> the inverse of the square of the earth-sun distance relative
814!! to the mean distance at angle ang in the earth's orbit.
815
816subroutine daily_mean_solar_1d_(lat, time_since_ae, cosz, h_out, rr_out)
817
818!----------------------------------------------------------------------
819real(kind=fms_ast_kind_), intent(in), dimension(:) :: lat
820real(kind=fms_ast_kind_), intent(in) :: time_since_ae
821real(kind=fms_ast_kind_), intent(out), dimension(size(lat(:))) :: cosz
822real(kind=fms_ast_kind_), intent(out), dimension(size(lat(:))) :: h_out
823real(kind=fms_ast_kind_), intent(out) :: rr_out
824
825!----------------------------------------------------------------------
826
827!----------------------------------------------------------------------
828! local variables
829!----------------------------------------------------------------------
830
831real(kind=fms_ast_kind_), dimension(size(lat),1) :: lat_2d, cosz_2d, hout_2d
832
833!--------------------------------------------------------------------
834!> define 2-d versions of input data array.
835!--------------------------------------------------------------------
836
837 lat_2d(:,1) = lat
838
839!--------------------------------------------------------------------
840!> call daily_mean_solar_2d to calculate astronomy fields.
841!--------------------------------------------------------------------
842
843 call daily_mean_solar(lat_2d, time_since_ae, cosz_2d, &
844 hout_2d, rr_out)
845
846!-------------------------------------------------------------------
847!> place output fields into 1-d arguments for return to
848!! calling routine.
849!-------------------------------------------------------------------
850
851 h_out = hout_2d(:,1)
852 cosz = cosz_2d(:,1)
853
854end subroutine daily_mean_solar_1d_
855
856
857!> @brief daily_mean_solar_2level takes 1-d input fields, adds a second
858!! dimension and calls daily_mean_solar_2d. on return, the 2d fields
859!! are returned to the original 1d fields.
860!!
861!! @param [in] <lat> Latitudes of model grid points
862!! @param [in] <time_since_ae> Time of year; autumnal equinox = 0.0, one year = 2 * pi
863!! @param [out] <cosz> Cosine of solar zenith angle
864!! @param [out] <solar> Shortwave flux factor: cosine of zenith angle * daylight fraction /
865!! (earth-sun distance squared)
866
867subroutine daily_mean_solar_2level_(lat, time_since_ae, cosz, solar)
868
869!----------------------------------------------------------------------
870real(kind=fms_ast_kind_), intent(in), dimension(:) :: lat
871real(kind=fms_ast_kind_), intent(in) :: time_since_ae
872real(kind=fms_ast_kind_), intent(out), dimension(size(lat(:))) :: cosz, solar
873
874!----------------------------------------------------------------------
875
876!----------------------------------------------------------------------
877! local variables
878!----------------------------------------------------------------------
879
880real(kind=fms_ast_kind_), dimension(size(lat),1) :: lat_2d, cosz_2d, hout_2d
881real(kind=fms_ast_kind_) :: rr_out
882
883!--------------------------------------------------------------------
884!> define 2-d versions of input data array.
885!--------------------------------------------------------------------
886
887 lat_2d(:,1) = lat
888
889!--------------------------------------------------------------------
890!> call daily_mean_solar_2d to calculate astronomy fields.
891!--------------------------------------------------------------------
892
893 call daily_mean_solar(lat_2d, time_since_ae, cosz_2d, &
894 hout_2d, rr_out)
895
896!-------------------------------------------------------------------
897!> place output fields into 1-d arguments for return to
898!! calling routine.
899!-------------------------------------------------------------------
900
901 solar = cosz_2d(:,1)*hout_2d(:,1)*rr_out
902 cosz = cosz_2d(:,1)
903
904end subroutine daily_mean_solar_2level_
905
906
907!> @brief daily_mean_solar_1d takes 1-d input fields, adds a second dimension
908!! and calls daily_mean_solar_2d. on return, the 2d fields are
909!! returned to the original 1d fields.
910!!
911!! @param [in] <lat> Latitudes of model grid points
912!! @param [in] <time_since_ae> Time of year; autumnal equinox = 0.0, one year = 2 * pi
913!! @param [out] <cosz> Cosine of solar zenith angle
914!! @param [out] <h_out> 2-d array of half-day lengths at the latitudes
915!! @param [out] <rr_out> the inverse of the square of the earth-sun distance relative to
916!! the mean distance at angle ang in the earth's orbit.
917
918subroutine daily_mean_solar_0d_(lat, time_since_ae, cosz, h_out, rr_out)
919
920real(kind=fms_ast_kind_), intent(in) :: lat, time_since_ae
921real(kind=fms_ast_kind_), intent(out) :: cosz, h_out, rr_out
922
923!--------------------------------------------------------------------
924! local variables
925!--------------------------------------------------------------------
926
927real(kind=fms_ast_kind_), dimension(1,1) :: lat_2d, cosz_2d, hout_2d
928
929!--------------------------------------------------------------------
930!> define 2-d versions of input data array.
931!--------------------------------------------------------------------
932
933 lat_2d = lat
934
935!--------------------------------------------------------------------
936!> call daily_mean_solar_2d to calculate astronomy fields.
937!--------------------------------------------------------------------
938
939 call daily_mean_solar(lat_2d, time_since_ae, cosz_2d, &
940 hout_2d, rr_out)
941
942!-------------------------------------------------------------------
943!> return output fields to scalars for return to calling routine.
944!-------------------------------------------------------------------
945
946 h_out = hout_2d(1,1)
947 cosz = cosz_2d(1,1)
948
949end subroutine daily_mean_solar_0d_
950
951
952!> @brief daily_mean_solar_cal_2d receives time_type inputs, converts
953!! them to real variables and then calls daily_mean_solar_2d to
954!! compute desired astronomical variables.
955!!
956!! @param [in] <lat> Latitudes of model grid points
957!! @param [in] <time> Time of year (time_type)
958!! @param [out] <cosz> Cosine of solar zenith angle
959!! @param [out] <fracday> Daylight fraction of time interval
960!! @param [out] <rrsun> Earth-Sun distance (r) relative to semi-major axis of orbital ellipse (a):(a/r)**2
961!!
962!! @throw FATAL, "astronomy_mod time_since_ae not between 0 and 2pi"
963
964subroutine daily_mean_solar_cal_2d_(lat, time, cosz, fracday, rrsun)
965
966!-------------------------------------------------------------------
967real(kind=fms_ast_kind_), dimension(:,:), intent(in) :: lat
968type(time_type), intent(in) :: time
969real(kind=fms_ast_kind_), dimension(:,:), intent(out) :: cosz, fracday
970real(kind=fms_ast_kind_), intent(out) :: rrsun
971!-------------------------------------------------------------------
972
973!-------------------------------------------------------------------
974! local variables
975!-------------------------------------------------------------------
976
977real(kind=fms_ast_kind_) :: time_since_ae
978integer, parameter :: lkind = fms_ast_kind_
979
980!--------------------------------------------------------------------
981! be sure the time in the annual cycle is legitimate.
982!---------------------------------------------------------------------
983
984 time_since_ae = real(orbital_time(time), fms_ast_kind_)
985 if (time_since_ae < 0.0_lkind .or. time_since_ae > real(twopi, fms_ast_kind_)) &
986 call error_mesg('astronomy_mod', 'time_since_ae not between 0 and 2pi', fatal)
987
988!--------------------------------------------------------------------
989! call daily_mean_solar_2d to calculate astronomy fields.
990!--------------------------------------------------------------------
991
992 call daily_mean_solar(lat, time_since_ae, cosz, fracday, rrsun)
993
994end subroutine daily_mean_solar_cal_2d_
995
996
997!> @brief daily_mean_solar_cal_1d receives time_type inputs, converts
998!! them to real, 2d variables and then calls daily_mean_solar_2d to
999!! compute desired astronomical variables.
1000!!
1001!! @param [in] <lat> Latitudes of model grid points
1002!! @param [in] <time> Time of year (time_type)
1003!! @param [out] <cosz> Cosine of solar zenith angle
1004!! @param [out] <fracday> Daylight fraction of time interval
1005!! @param [out] <rrsun> Earth-Sun distance (r) relative to semi-major axis of orbital ellipse (a):(a/r)**2
1006
1007subroutine daily_mean_solar_cal_1d_(lat, time, cosz, fracday, rrsun)
1008
1009real(kind=fms_ast_kind_), dimension(:), intent(in) :: lat
1010type(time_type), intent(in) :: time
1011real(kind=fms_ast_kind_), dimension(:), intent(out) :: cosz, fracday
1012real(kind=fms_ast_kind_), intent(out) :: rrsun
1013
1014!---------------------------------------------------------------------
1015! local variables
1016!---------------------------------------------------------------------
1017
1018real(kind=fms_ast_kind_), dimension(size(lat),1) :: lat_2d, cosz_2d, fracday_2d
1019
1020!--------------------------------------------------------------------
1021!> define 2-d versions of input data array.
1022!--------------------------------------------------------------------
1023
1024 lat_2d(:,1) = lat
1025
1026!--------------------------------------------------------------------
1027!> call daily_mean_solar_cal_2d to convert the time_types to reals and
1028!! then calculate the astronomy fields.
1029!--------------------------------------------------------------------
1030
1031 call daily_mean_solar(lat_2d, time, cosz_2d, fracday_2d, rrsun)
1032
1033!-------------------------------------------------------------------
1034!> place output fields into 1-d arguments for return to
1035!! calling routine.
1036!-------------------------------------------------------------------
1037
1038 fracday = fracday_2d(:,1)
1039 cosz = cosz_2d(:,1)
1040
1041end subroutine daily_mean_solar_cal_1d_
1042
1043
1044!> @brief daily_mean_solar_cal_2level receives 1d arrays and time_type input,
1045!! converts them to real, 2d variables and then calls
1046!! daily_mean_solar_2d to compute desired astronomical variables.
1047!!
1048!! @param [in] <lat> Latitudes of model grid points
1049!! @param [in] <time> Time of year (time_type)
1050!! @param [out] <cosz> Cosine of solar zenith angle
1051!! @param [out] <solar> Shortwave flux factor: cosine of zenith angle * daylight fraction /
1052!! (earth-sun distance squared)
1053
1054subroutine daily_mean_solar_cal_2level_(lat, time, cosz, solar)
1055
1056real(kind=fms_ast_kind_), dimension(:), intent(in) :: lat
1057type(time_type), intent(in) :: time
1058real(kind=fms_ast_kind_), dimension(:), intent(out) :: cosz, solar
1059!---------------------------------------------------------------------
1060! local variables
1061!---------------------------------------------------------------------
1062
1063real(kind=fms_ast_kind_), dimension(size(lat),1) :: lat_2d, cosz_2d, fracday_2d
1064real(kind=fms_ast_kind_) :: rrsun
1065
1066
1067!--------------------------------------------------------------------
1068!> define 2-d versions of input data array.
1069!--------------------------------------------------------------------
1070
1071 lat_2d(:,1) = lat
1072
1073!--------------------------------------------------------------------
1074!> call daily_mean_solar_cal_2d to convert the time_types to reals and
1075!! then calculate the astronomy fields.
1076!--------------------------------------------------------------------
1077
1078 call daily_mean_solar(lat_2d, time, cosz_2d, fracday_2d, rrsun)
1079
1080!-------------------------------------------------------------------
1081!> place output fields into 1-d arguments for return to
1082!! calling routine.
1083!-------------------------------------------------------------------
1084
1085 solar = cosz_2d(:,1)*fracday_2d(:,1)*rrsun
1086 cosz = cosz_2d(:,1)
1087
1088end subroutine daily_mean_solar_cal_2level_
1089
1090
1091!> @brief daily_mean_solar_cal_0d converts scalar input fields to real, 2d variables and
1092!! then calls daily_mean_solar_2d to compute desired astronomical variables.
1093!!
1094!! @param [in] <lat> Latitudes of model grid points
1095!! @param [in] <time> Time of year (time_type)
1096!! @param [out] <cosz> Cosine of solar zenith angle
1097!! @param [out] <fracday> Daylight fraction of time interval
1098!! @param [out] <rrsun> Earth-Sun distance (r) relative to semi-major axis of orbital ellipse (a):(a/r)**2
1099
1100subroutine daily_mean_solar_cal_0d_(lat, time, cosz, fracday, rrsun)
1101
1102!--------------------------------------------------------------------
1103real(kind=fms_ast_kind_), intent(in) :: lat
1104type(time_type), intent(in) :: time
1105real(kind=fms_ast_kind_), intent(out) :: cosz, fracday, rrsun
1106!--------------------------------------------------------------------
1107
1108!--------------------------------------------------------------------
1109! local variables
1110!--------------------------------------------------------------------
1111
1112real(kind=fms_ast_kind_), dimension(1,1) :: lat_2d, cosz_2d, fracday_2d
1113
1114!--------------------------------------------------------------------
1115!> define 2-d versions of input data array.
1116!--------------------------------------------------------------------
1117
1118 lat_2d = lat
1119
1120!--------------------------------------------------------------------
1121!> call daily_mean_solar_cal_2d to convert the time_types to reals and
1122!! then calculate the astronomy fields.
1123!--------------------------------------------------------------------
1124
1125 call daily_mean_solar(lat_2d, time, cosz_2d, fracday_2d, rrsun)
1126
1127!-------------------------------------------------------------------
1128!> place output fields into scalar arguments for return to
1129!! calling routine.
1130!-------------------------------------------------------------------
1131
1132 fracday = fracday_2d(1,1)
1133 cosz = cosz_2d(1,1)
1134
1135end subroutine daily_mean_solar_cal_0d_
1136
1137
1138!> @brief annual_mean_solar_2d returns 2d fields of annual mean values of the cosine of
1139!! zenith angle, daylight fraction and earth-sun distance at the specified latitude.
1140!!
1141!! @param [in] <jst> Starting index of latitude window
1142!! @param [in] <jnd> Ending index of latitude window
1143!! @param [in] <lat> Latitudes of model grid points
1144!! @param [out] <cosz> Cosine of solar zenith angle
1145!! @param [out] <solar> Shortwave flux factor: cosine of zenith angle * daylight fraction /
1146!! (earth-sun distance squared)
1147!! @param [out] <fracday> Daylight fraction of time interval
1148!! @param [out] <rrsun> Earth-Sun distance (r) relative to semi-major axis of orbital ellipse (a):(a/r)**2
1149
1150subroutine annual_mean_solar_2d_(js, je, lat, cosz, solar, fracday, &
1151 rrsun)
1152
1153!--------------------------------------------------------------------
1154integer, intent(in) :: js, je
1155real(kind=fms_ast_kind_), dimension(:,:), intent(in) :: lat
1156real(kind=fms_ast_kind_), dimension(:,:), intent(out) :: solar, cosz, fracday
1157real(kind=fms_ast_kind_), intent(out) :: rrsun
1158!--------------------------------------------------------------------
1159
1160!--------------------------------------------------------------------
1161! local variables
1162!--------------------------------------------------------------------
1163
1164real(kind=fms_ast_kind_), dimension(size(lat,1),size(lat,2)) :: s,z
1165real(kind=fms_ast_kind_) :: t
1166integer :: n
1167integer, parameter :: lkind = fms_ast_kind_
1168
1169!--------------------------------------------------------------------
1170! if the calculation has not yet been done, do it here.
1171!--------------------------------------------------------------------
1172
1173 if (.not. annual_mean_calculated) then
1174
1175!----------------------------------------------------------------------
1176!> determine annual mean values of solar flux and product of cosz
1177!! and solar flux by integrating the annual cycle in num_angles
1178!! orbital increments.
1179!----------------------------------------------------------------------
1180
1181 solar = 0.0_lkind
1182 cosz = 0.0_lkind
1183 do n =1, num_angles
1184 t = real((n-1), fms_ast_kind_)*real(twopi, fms_ast_kind_)/real(num_angles, fms_ast_kind_)
1185 call daily_mean_solar(lat,t, z, fracday, rrsun)
1186 s = z*rrsun*fracday
1187 solar = solar + s
1188 cosz = cosz + z*s
1189 end do
1190 solar = solar/real(num_angles, fms_ast_kind_)
1191 cosz = cosz/real(num_angles, fms_ast_kind_)
1192
1193!--------------------------------------------------------------------
1194!> define the flux-weighted annual mean cosine of the zenith angle.
1195!--------------------------------------------------------------------
1196
1197 where(solar .eq. 0.0_lkind)
1198 cosz = 0.0_lkind
1199 elsewhere
1200 cosz = cosz/solar
1201 end where
1202
1203!-------------------------------------------------------------------
1204!> define avg fracday such as to make the avg flux (solar) equal to
1205!! the product of the avg cosz * avg fracday * assumed mean avg
1206!! radius of 1.0. it is unlikely that these avg fracday and avg rr
1207!! will ever be used.
1208!--------------------------------------------------------------------
1209
1210 where(solar .eq. 0.0_lkind)
1211 fracday = 0.0_lkind
1212 elsewhere
1213 fracday = solar/cosz
1214 end where
1215 rrsun = 1.00_lkind
1216
1217!---------------------------------------------------------------------
1218!> save the values that have been calculated as module variables, if
1219!! those variables are present; i.e., not the spectral 2-layer model.
1220!---------------------------------------------------------------------
1221 if (allocated (cosz_ann)) then
1222
1223 cosz_ann = real(cosz, r8_kind)
1224 solar_ann = real(solar, r8_kind)
1225 fracday_ann = real(fracday, r8_kind)
1226 rrsun_ann = real(rrsun, r8_kind)
1227
1228!--------------------------------------------------------------------
1229!> increment the points computed counter. set flag to end execution
1230!! once values have been calculated for all points owned by the
1231!! processor.
1232!--------------------------------------------------------------------
1233 num_pts = num_pts + size(lat,1)*size(lat,2)
1234 if ( num_pts == total_pts) annual_mean_calculated = .true.
1235 endif
1236
1237!--------------------------------------------------------------------
1238!> if the calculation has been done, return the appropriate module
1239!! variables.
1240!--------------------------------------------------------------------
1241 else
1242 if (allocated (cosz_ann)) then
1243
1244 cosz = real(cosz_ann, fms_ast_kind_)
1245 solar = real(solar_ann, fms_ast_kind_)
1246 fracday = real(fracday_ann, fms_ast_kind_)
1247 rrsun = real(rrsun_ann, fms_ast_kind_)
1248 endif
1249 endif
1250
1251end subroutine annual_mean_solar_2d_
1252
1253
1254!> @brief annual_mean_solar_1d creates 2-d input fields from 1-d input fields and then calls
1255!! annual_mean_solar_2d to obtain 2-d output fields which are then stored into 1-d
1256!! fields for return to the calling subroutine.
1257!!
1258!! @param [in] <jst> Starting index of latitude window
1259!! @param [in] <jnd> Ending index of latitude window
1260!! @param [in] <lat> Latitudes of model grid points
1261!! @param [out] <cosz> Cosine of solar zenith angle
1262!! @param [out] <solar> Shortwave flux factor: cosine of zenith angle * daylight fraction /
1263!! (earth-sun distance squared)
1264!! @param [out] <fracday> Daylight fraction of time interval
1265!! @param [out] <rrsun_out> Earth-Sun distance (r) relative to semi-major axis of orbital ellipse (a):(a/r)**2
1266
1267subroutine annual_mean_solar_1d_(jst, jnd, lat, cosz, solar, &
1268 fracday, rrsun_out)
1269
1270!---------------------------------------------------------------------
1271integer, intent(in) :: jst, jnd
1272real(kind=fms_ast_kind_), dimension(:), intent(in) :: lat(:)
1273real(kind=fms_ast_kind_), dimension(:), intent(out) :: cosz, solar, fracday
1274real(kind=fms_ast_kind_), intent(out) :: rrsun_out
1275!---------------------------------------------------------------------
1276
1277!---------------------------------------------------------------------
1278! local variables
1279
1280
1281real(kind=fms_ast_kind_), dimension(size(lat),1) :: lat_2d, solar_2d, cosz_2d, fracday_2d
1282real(kind=fms_ast_kind_) :: rrsun
1283
1284!--------------------------------------------------------------------
1285! if the calculation has not been done, do it here.
1286!--------------------------------------------------------------------
1287
1288 if ( .not. annual_mean_calculated) then
1289
1290!--------------------------------------------------------------------
1291!> define 2-d versions of input data array.
1292!--------------------------------------------------------------------
1293
1294 lat_2d(:,1) = lat
1295
1296!--------------------------------------------------------------------
1297!> call annual_mean_solar_2d to calculate the astronomy fields.
1298!--------------------------------------------------------------------
1299
1300 call annual_mean_solar(jst, jnd, lat_2d, cosz_2d, solar_2d, fracday_2d, rrsun)
1301
1302!-------------------------------------------------------------------
1303!> place output fields into 1-D arrays for return to calling routine.
1304!-------------------------------------------------------------------
1305
1306 fracday = fracday_2d(:,1)
1307 rrsun_out = rrsun
1308 solar = solar_2d(:,1)
1309 cosz = cosz_2d(:,1)
1310
1311!--------------------------------------------------------------------
1312!> if the calculation has been done, simply return the module
1313!! variables contain the results at the desired latitudes.
1314!--------------------------------------------------------------------
1315
1316 else
1317 cosz(:) = real(cosz_ann(1,jst:jnd), fms_ast_kind_)
1318 solar(:) = real(solar_ann(1,jst:jnd), fms_ast_kind_)
1319 fracday(:) = real(fracday_ann(1,jst:jnd), fms_ast_kind_)
1320 rrsun = real(rrsun_ann, fms_ast_kind_)
1321 endif
1322
1323end subroutine annual_mean_solar_1d_
1324
1325
1326!> @brief annual_mean_solar_2level creates 2-d input fields from 1-d input fields
1327!! and then calls annual_mean_solar_2d to obtain 2-d output fields which are
1328!! then stored into 1-d fields for return to the calling subroutine. This
1329!! subroutine will be called during model initialization.
1330!!
1331!! @throw FATAL, "astronomy_mod annual_mean_solar_2level should be called only once"
1332
1333subroutine annual_mean_solar_2level_(lat, cosz, solar)
1334
1335!---------------------------------------------------------------------
1336real(kind=fms_ast_kind_), dimension(:), intent(in) :: lat !< Latitudes of model grid points
1337real(kind=fms_ast_kind_), dimension(:), intent(out) :: cosz !< Cosine of solar zenith angle
1338real(kind=fms_ast_kind_), dimension(:), intent(out) :: solar !< shortwave flux factor: cosine of zenith angle *
1339 !! daylight fraction / (earth-sun distance squared)
1340!---------------------------------------------------------------------
1341
1342!---------------------------------------------------------------------
1343! local variables
1344
1345
1346real(kind=fms_ast_kind_), dimension(size(lat),1) :: lat_2d, solar_2d, cosz_2d, fracday_2d
1347integer :: jst, jnd
1348real(kind=fms_ast_kind_) :: rrsun
1349
1350!--------------------------------------------------------------------
1351! if the calculation has not been done, do it here.
1352!--------------------------------------------------------------------
1353
1354 if ( .not. annual_mean_calculated) then
1355
1356!--------------------------------------------------------------------
1357!> define 2-d versions of input data array.
1358!--------------------------------------------------------------------
1359
1360 lat_2d(:,1) = lat
1361 jst = 1
1362 jnd = size(lat(:))
1363
1364!--------------------------------------------------------------------
1365!> call annual_mean_solar_2d to calculate the astronomy fields.
1366!--------------------------------------------------------------------
1367
1368 call annual_mean_solar(jst, jnd, lat_2d, cosz_2d, solar_2d, fracday_2d, rrsun)
1369
1370!-------------------------------------------------------------------
1371!> place output fields into 1-D arrays for return to calling routine.
1372!-------------------------------------------------------------------
1373
1374 solar = solar_2d(:,1)
1375 cosz = cosz_2d(:,1)
1376
1377!--------------------------------------------------------------------
1378!> if the calculation has been done, print an error message since
1379!! this subroutine should be called only once.
1380!--------------------------------------------------------------------
1381
1382 else
1383 call error_mesg('astronomy_mod', 'annual_mean_solar_2level should be called only once', fatal)
1384 endif
1385 annual_mean_calculated = .true.
1386
1387end subroutine annual_mean_solar_2level_
1388
1389
1390!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1391!
1392! PRIVATE SUBROUTINES
1393!
1394!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1395
1396
1397
1398!> @brief r_inv_squared returns the inverse of the square of the earth-sun
1399!! distance relative to the mean distance at angle ang in the Earth's orbit.
1400function r_inv_squared_(ang)
1401
1402!--------------------------------------------------------------------
1403real(kind=fms_ast_kind_), intent(in) :: ang !< angular position of earth in its orbit, relative to a
1404 !! value of 0.0 at the NH autumnal equinox, value between
1405 !! 0.0 and 2 * pi [radians]
1406!--------------------------------------------------------------------
1407
1408!---------------------------------------------------------------------
1409! local variables
1410
1411
1412real(kind=fms_ast_kind_) :: r_inv_squared_ !< The inverse of the square of the earth-sun distance relative
1413 !! to the mean distance [dimensionless]
1414real(kind=fms_ast_kind_) :: r !< Earth-Sun distance relative to mean distance [dimensionless]
1415real(kind=fms_ast_kind_) :: rad_per !< Angular position of perihelion [radians]
1416integer, parameter :: lkind = fms_ast_kind_
1417
1418!--------------------------------------------------------------------
1419!> define the earth-sun distance (r) and then return the inverse of
1420!! its square (r_inv_squared) to the calling routine.
1421!--------------------------------------------------------------------
1422
1423 rad_per = real(per, fms_ast_kind_)*real(deg_to_rad, fms_ast_kind_)
1424 r = (1.0_lkind - real(ecc, fms_ast_kind_)**2)/(1.0_lkind &
1425 + real(ecc, fms_ast_kind_)*cos(ang - rad_per))
1426 r_inv_squared_ = r**(-2)
1427
1428
1429end function r_inv_squared_
1430
1431
1432!> @brief angle determines the position within the earth's orbit at time t
1433!! in the year (t = 0 at NH autumnal equinox) by interpolating
1434!! into the orbital position table.
1435
1436function angle_(t)
1437
1438!--------------------------------------------------------------------
1439real(kind=fms_ast_kind_), intent(in) :: t !< time of year (between 0 and 2*pi; t=0 at NH autumnal equinox
1440!--------------------------------------------------------------------
1441
1442!--------------------------------------------------------------------
1443! local variables
1444!--------------------------------------------------------------------
1445
1446real(kind=fms_ast_kind_) :: angle_ !< Orbital position relative to NH autumnal equinox [radians]
1447real(kind=fms_ast_kind_) :: norm_time !< Index into orbital table corresponding to input time [dimensionless]
1448real(kind=fms_ast_kind_) :: x !< Fractional distance between the orbital table entries bracketing
1449 !! the input time [dimensionless]
1450integer :: int !< Table index which is lower than actual position, but closest to it [dimensionless]
1451integer :: int_1 !< Next table index just larger than actual orbital position [dimensionless]
1452integer, parameter :: lkind = fms_ast_kind_
1453
1454!--------------------------------------------------------------------
1455!> Define orbital tables indices bracketing current orbital time
1456!! (int and int_1). Define table index distance between the lower
1457!! table value (int) and the actual orbital time (x). Define orbital
1458!! position as being x of the way between int and int_1. Renormalize
1459!! angle to be within the range 0 to 2*pi.
1460!--------------------------------------------------------------------
1461
1462 norm_time = t * real(num_angles, fms_ast_kind_)/real(twopi, fms_ast_kind_)
1463 int = floor(norm_time)
1464 int = modulo(int,num_angles)
1465 int_1 = int+1
1466 x = norm_time - real(floor(norm_time), fms_ast_kind_)
1467 angle_ = (1.0_lkind - x) * real(orb_angle(int), fms_ast_kind_) &
1468 + x * real(orb_angle(int_1), fms_ast_kind_)
1469 angle_ = modulo(angle_, real(twopi, fms_ast_kind_))
1470
1471end function angle_
1472
1473
1474!> @brief Declination returns the solar declination angle at orbital
1475!! position ang in earth's orbit.
1476
1477function declination_(ang)
1478
1479!--------------------------------------------------------------------
1480real(kind=fms_ast_kind_), intent(in) :: ang !< solar orbital position ang in earth's orbit
1481
1482!--------------------------------------------------------------------
1483
1484!--------------------------------------------------------------------
1485! local variables
1486
1487
1488real(kind=fms_ast_kind_) :: declination_ !< Solar declination angle [radians]
1489real(kind=fms_ast_kind_) :: rad_obliq !< Obliquity of the ecliptic [radians]
1490real(kind=fms_ast_kind_) :: sin_dec !< Sine of the solar declination [dimensionless]
1491
1492!---------------------------------------------------------------------
1493! compute the solar declination.
1494!---------------------------------------------------------------------
1495
1496 rad_obliq = real(obliq, fms_ast_kind_)*real(deg_to_rad, fms_ast_kind_)
1497 sin_dec = - sin(rad_obliq)*sin(ang)
1498 declination_ = asin(sin_dec)
1499
1500end function declination_
1501
1502
1503!> @brief half_day_2d returns a 2-d array of half-day lengths at the
1504!! latitudes and declination provided.
1505!!
1506
1507function half_day_2d_(latitude, dec) result(h)
1508
1509!---------------------------------------------------------------------
1510real(kind=fms_ast_kind_), dimension(:,:), intent(in) :: latitude !< Latitutde of view point
1511real(kind=fms_ast_kind_), intent(in) :: dec !< Solar declination angle at view point
1512real(kind=fms_ast_kind_), dimension(size(latitude,1),size(latitude,2)) :: h
1513!---------------------------------------------------------------------
1514
1515!---------------------------------------------------------------------
1516! local variables
1517!---------------------------------------------------------------------
1518
1519real(kind=fms_ast_kind_), dimension (size(latitude,1),size(latitude,2)):: &
1520 cos_half_day, & !< Cosine of half-day length [dimensionless]
1521 lat !< Model latitude, adjusted so that it is never 0.5*pi or -0.5*pi
1522real(kind=fms_ast_kind_) :: tan_dec !< tangent of solar declination [dimensionless]
1523integer, parameter :: lkind = fms_ast_kind_
1524real(kind=fms_ast_kind_) :: eps = 1.0e-05_lkind !< small increment
1525
1526
1527!--------------------------------------------------------------------
1528!> define tangent of the declination.
1529!--------------------------------------------------------------------
1530
1531 tan_dec = tan(dec)
1532
1533!--------------------------------------------------------------------
1534!> adjust latitude so that its tangent will be defined.
1535!--------------------------------------------------------------------
1536
1537 lat = latitude
1538 where (latitude == 0.5_lkind*real(pi, fms_ast_kind_)) lat = latitude - eps
1539 where (latitude == -0.5_lkind*real(pi, fms_ast_kind_)) lat = latitude + eps
1540
1541!--------------------------------------------------------------------
1542!> define the cosine of the half-day length. adjust for cases of
1543!! all daylight or all night.
1544!--------------------------------------------------------------------
1545
1546 cos_half_day = -tan(lat)*tan_dec
1547 where (cos_half_day <= -1.0_lkind) h = real(pi, fms_ast_kind_)
1548 where (cos_half_day >= 1.0_lkind) h = 0.0_lkind
1549 where (cos_half_day > -1.0_lkind .and. &
1550 cos_half_day < 1.0_lkind) h = acos(cos_half_day)
1551
1552end function half_day_2d_
1553
1554
1555!> @brief half_day_0d takes scalar input fields, makes them into 2-d fields
1556!! dimensioned (1,1), and calls half_day_2d. On return, the 2-d
1557!! fields are converted to the desired scalar output.
1558!!
1559!! @param [in] <latitude> Latitutde of view point
1560!! @param [in] <dec> Solar declination angle at view point
1561
1562function half_day_0d_(latitude, dec) result(h)
1563
1564real(kind=fms_ast_kind_), intent(in) :: latitude, dec
1565real(kind=fms_ast_kind_) :: h
1566
1567!----------------------------------------------------------------------
1568! local variables
1569!----------------------------------------------------------------------
1570
1571real(kind=fms_ast_kind_), dimension(1,1) :: lat_2d, h_2d
1572
1573!---------------------------------------------------------------------
1574! create 2d array from the input latitude field.
1575!---------------------------------------------------------------------
1576
1577 lat_2d = latitude
1578
1579!---------------------------------------------------------------------
1580! call half_day with the 2d arguments to calculate half-day length.
1581!---------------------------------------------------------------------
1582
1583 h_2d = half_day(lat_2d, dec)
1584
1585!---------------------------------------------------------------------
1586! create scalar from 2d array.
1587!---------------------------------------------------------------------
1588
1589 h = h_2d(1,1)
1590
1591end function half_day_0d_
1592
1593
1594!> @}
1595! close documentation grouping
daily_mean_solar_2d_
subroutine daily_mean_solar_2d_(lat, time_since_ae, cosz, h_out, rr_out)
daily_mean_solar_2d computes the daily mean astronomical parameters for the input points at latitude ...
Definition astronomy.inc:754
annual_mean_solar_2d_
subroutine annual_mean_solar_2d_(js, je, lat, cosz, solar, fracday, rrsun)
annual_mean_solar_2d returns 2d fields of annual mean values of the cosine of zenith angle,...
Definition astronomy.inc:1152
diurnal_solar_cal_0d_
subroutine diurnal_solar_cal_0d_(lat, lon, time, cosz, fracday, rrsun, dt_time, allow_negative_cosz, half_day_out)
diurnal_solar_cal_0d receives time_type inputs, converts them to real variables and then calls diurna...
Definition astronomy.inc:686
get_orbital_parameters_
subroutine get_orbital_parameters_(ecc_out, obliq_out, per_out)
get_orbital_parameters retrieves the orbital parameters for use by another module.
Definition astronomy.inc:91
diurnal_solar_0d_
subroutine diurnal_solar_0d_(lat, lon, gmt, time_since_ae, cosz, fracday, rrsun, dt, allow_negative_cosz, half_day_out)
diurnal_solar_0d takes scalar input fields, makes them into 2d arrays dimensioned (1,...
Definition astronomy.inc:459
diurnal_solar_1d_
subroutine diurnal_solar_1d_(lat, lon, gmt, time_since_ae, cosz, fracday, rrsun, dt, allow_negative_cosz, half_day_out)
diurnal_solar_1d takes 1-d input fields, adds a second dimension and calls diurnal_solar_2d....
Definition astronomy.inc:387
declination_
real(kind=fms_ast_kind_) function declination_(ang)
Declination returns the solar declination angle at orbital position ang in earth's orbit.
Definition astronomy.inc:1478
daily_mean_solar_cal_2d_
subroutine daily_mean_solar_cal_2d_(lat, time, cosz, fracday, rrsun)
daily_mean_solar_cal_2d receives time_type inputs, converts them to real variables and then calls dai...
Definition astronomy.inc:965
diurnal_solar_2d_
subroutine diurnal_solar_2d_(lat, lon, gmt, time_since_ae, cosz, fracday, rrsun, dt, allow_negative_cosz, half_day_out)
diurnal_solar_2d returns 2d fields of cosine of zenith angle, daylight fraction and earth-sun distanc...
Definition astronomy.inc:145
daily_mean_solar_0d_
subroutine daily_mean_solar_0d_(lat, time_since_ae, cosz, h_out, rr_out)
daily_mean_solar_1d takes 1-d input fields, adds a second dimension and calls daily_mean_solar_2d....
Definition astronomy.inc:919
half_day_2d_
real(kind=fms_ast_kind_) function, dimension(size(latitude, 1), size(latitude, 2)) half_day_2d_(latitude, dec)
half_day_2d returns a 2-d array of half-day lengths at the latitudes and declination provided.
Definition astronomy.inc:1508
r_inv_squared_
real(kind=fms_ast_kind_) function r_inv_squared_(ang)
r_inv_squared returns the inverse of the square of the earth-sun distance relative to the mean distan...
Definition astronomy.inc:1401
diurnal_solar_cal_1d_
subroutine diurnal_solar_cal_1d_(lat, lon, time, cosz, fracday, rrsun, dt_time, allow_negative_cosz, half_day_out)
diurnal_solar_cal_1d receives time_type inputs, converts them to real variables and then calls diurna...
Definition astronomy.inc:613
set_orbital_parameters_
subroutine set_orbital_parameters_(ecc_in, obliq_in, per_in)
set_orbital_parameters saves the input values of eccentricity, obliquity and perihelion time as modul...
Definition astronomy.inc:40
daily_mean_solar_cal_0d_
subroutine daily_mean_solar_cal_0d_(lat, time, cosz, fracday, rrsun)
daily_mean_solar_cal_0d converts scalar input fields to real, 2d variables and then calls daily_mean_...
Definition astronomy.inc:1101
half_day_0d_
real(kind=fms_ast_kind_) function half_day_0d_(latitude, dec)
half_day_0d takes scalar input fields, makes them into 2-d fields dimensioned (1,1),...
Definition astronomy.inc:1563
angle_
real(kind=fms_ast_kind_) function angle_(t)
angle determines the position within the earth's orbit at time t in the year (t = 0 at NH autumnal eq...
Definition astronomy.inc:1437
daily_mean_solar_2level_
subroutine daily_mean_solar_2level_(lat, time_since_ae, cosz, solar)
daily_mean_solar_2level takes 1-d input fields, adds a second dimension and calls daily_mean_solar_2d...
Definition astronomy.inc:868
annual_mean_solar_1d_
subroutine annual_mean_solar_1d_(jst, jnd, lat, cosz, solar, fracday, rrsun_out)
annual_mean_solar_1d creates 2-d input fields from 1-d input fields and then calls annual_mean_solar_...
Definition astronomy.inc:1269
daily_mean_solar_cal_2level_
subroutine daily_mean_solar_cal_2level_(lat, time, cosz, solar)
daily_mean_solar_cal_2level receives 1d arrays and time_type input, converts them to real,...
Definition astronomy.inc:1055
annual_mean_solar_2level_
subroutine annual_mean_solar_2level_(lat, cosz, solar)
annual_mean_solar_2level creates 2-d input fields from 1-d input fields and then calls annual_mean_so...
Definition astronomy.inc:1334
daily_mean_solar_1d_
subroutine daily_mean_solar_1d_(lat, time_since_ae, cosz, h_out, rr_out)
daily_mean_solar_1d takes 1-d input fields, adds a second dimension and calls daily_mean_solar_2d....
Definition astronomy.inc:817
diurnal_solar_cal_2d_
subroutine diurnal_solar_cal_2d_(lat, lon, time, cosz, fracday, rrsun, dt_time, allow_negative_cosz, half_day_out)
diurnal_solar_cal_2d receives time_type inputs, converts them to real variables and then calls diurna...
Definition astronomy.inc:527
daily_mean_solar_cal_1d_
subroutine daily_mean_solar_cal_1d_(lat, time, cosz, fracday, rrsun)
daily_mean_solar_cal_1d receives time_type inputs, converts them to real, 2d variables and then calls...
Definition astronomy.inc:1008