Return the daylight duration in seconds (i.e, the time between sunrise and sunset) for its MAXIMUM at the summer solstice. The routine is accurate enough for most purposes, but does not solve for precise sunrise/sunset times influenced by changes in the orbital parameters between those times. The time input does NOT need to be noon — it is used simply to evaluate the obliquity needed for the maximum daylight duration calculation. In order of preference, time is taken from currTime, if present, or else the currTime of CLOCK, if present, or else the currTime of the ORBIT’s associated clock. Note: Unless ORBIT_ANAL2B, the obliquity is fixed and the time is irrelevant.
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
type(MAPL_SunOrbit), | intent(in) | :: | ORBIT | |||
real, | intent(in), | dimension(:) | :: | LATS | ||
real, | intent(out), | dimension(:) | :: | DAYL | ||
type(ESMF_Time), | intent(in), | optional | :: | currTime | ||
type(ESMF_Clock), | intent(in), | optional | :: | CLOCK | ||
integer, | intent(out), | optional | :: | RC |
subroutine MAPL_SunGetDaylightDurationMax(ORBIT,LATS,DAYL,currTime,CLOCK,RC) ! !ARGUMENTS: type(MAPL_SunOrbit), intent(IN ) :: ORBIT real, dimension(:) , intent(IN ) :: LATS real, dimension(:) , intent(OUT) :: DAYL type(ESMF_Time) , optional, intent(IN ) :: currTime type(ESMF_Clock) , optional, intent(IN ) :: CLOCK integer, optional, intent(OUT) :: RC !EOPI character(len=ESMF_MAXSTR), parameter :: IAm = "MAPL_SunGetDaylightDurationMax" integer :: STATUS type(ESMF_Time) :: CURRENTTIME real(ESMF_KIND_R8) :: days real :: OBQ _ASSERT(MAPL_SunOrbitCreated(ORBIT),'MAPL_SunOrbit not yet created!') ! Which time? if (present(currTime)) then CURRENTTIME = CURRTIME else if (present(CLOCK)) then call ESMF_ClockGet( CLOCK, currTime=CURRENTTIME, RC=STATUS) else call ESMF_ClockGet(ORBIT%CLOCK, currTime=CURRENTTIME, RC=STATUS) end if _VERIFY(STATUS) end if if (ORBIT%ANAL2B) then ! time variation in obliquity from ref time call ESMF_TimeIntervalGet( & CURRENTTIME - ORBIT%ORB2B_TIME_REF, & d_r8=days, rc=STATUS) _VERIFY(STATUS) OBQ = ORBIT%ORB2B_OBQ_REF + days * ORBIT%ORB2B_OBQ_RATE else ! obliquity fixed in this case OBQ = ORBIT%OB * (MAPL_PI/180.) endif _ASSERT(OBQ >= 0, 'Obliquity less than 0 detected!') _ASSERT(OBQ < MAPL_PI, 'Obliquity greater than pi detected!') ! Maximum daylight duration at summer solstice [secs] ! (an even function of latitude) DAYL = (86400./MAPL_PI)*acos(min(1.,max(-1., & -tan(ABS(LATS))*tan(OBQ)))) _RETURN(ESMF_SUCCESS) end subroutine MAPL_SunGetDaylightDurationMax