Документ взят из кэша поисковой машины. Адрес оригинального документа : http://www.stsci.edu/~mperrin/software/sources/suntimes.pro Дата изменения: Sat Aug 11 03:00:40 2007 Дата индексирования: Sat Mar 1 14:34:59 2014 Кодировка: Поисковые слова: horizon

;+
; NAME: suntimes
;
; PURPOSE:
; compute APPROXIMATE sunrise and sunset times
;
; NOTES:
; Uses algorithms from the
; Buie Lowell IDL library.
;
; This is not likely to be the most accurate thing in the world.
;
; Particularly since it doesn't account for the motion of the sun
; during the day in question.
;
;
; USAGE:
; suntimes,jd,rise,set, [/today, /lst, obsname="lick", degrees=12]
;
; INPUTS:
; KEYWORDS:
; degrees= how many degrees below the horizon?
; Default is 0.25, so the sun is just tangent to the horizon
; /lst return results in LST, not CT
;
; OUTPUTS:
;
; HISTORY:
; Began 2005-07-14 16:07:24 by Marshall Perrin
; see 'obsplan.pro' from the Buie library for algorithm source.
;-

PRO suntimes,jd,rise,set,today=today,date=date,obsname=obsname,$
degrees=degrees,lst=uselst,print=print

gmt_offset_days = gmt_offsec()/84600.
gmt_offset_hours = gmt_offsec()/3600.

if n_elements(degrees) eq 0 then degrees=0.25
if ~(keyword_set(jd)) and ~(keyword_set(date)) then today=1
if keyword_set(today) then begin
; julday by itself doesn't know about time zones.
; So you have to add in an offset to take into account the
; fact that the computer's clock may not be in GMT.
jd = julday()+gmt_offset_days
daycnv,jd-gmt_offset_days,yr,mon,day,hr
date = [yr,mon,day,hr]
endif else $
if keyword_set(date) then begin
if n_elements(date) eq 3 then hr = 0 else hr = date[3]
jdcnv,date[0],date[1],date[2],hr,jd
jd += gmt_offset_days
endif

if ~(keyword_set(obsname)) then obsname='lick'

observatory,obsname,obs
lat=obs.latitude*!dtor
lon=obs.longitude*!dtor
alt=obs.altitude

; Sun position at input JD
sunpos,jd,sunra,sundec,/RADIAN

; Define night, Sun set to sun rise.

; find jd of midnight
; daycnv,jd,y,m,d,h
; jdcnv,y,m,d,gmt_offset_hours,jdlclmid

am = airmass(jd,sunra,sundec,lat,lon,alt=alt,lha=lha,lst=lst)
hatojd,!dpi,sunra,lst,jd,jdlclmid ; jd of nearest local midnight
lsidtim,jdlclmid,lon,midlst ; LST at local midnight
jdofmid = float(long(jdlclmid+0.5d0))-0.5d0
jdstr,jdofmid,100,thisdate

; alt needs to be NEGATIVE to be below the horizon!
altoha,(-1)*degrees/!radeg,sundec,lat,sunhorzha,sunhorztype

; JD of sunset/sunrise, AT, NT, CT
jdsset = jdlclmid - (!dpi-sunhorzha)/2.0d0/!dpi
jdsrise = jdlclmid + (!dpi-sunhorzha)/2.0d0/!dpi

daycnv,jdsset,set_y,set_m,set_d,set_h
daycnv,jdsrise,rise_y,rise_m,rise_d,rise_h

set = set_h-gmt_offset_hours
if set lt 0 then set += 24.0
rise = rise_h-gmt_offset_hours

; print,"set: ", set
; print,"rise: ", rise

if keyword_set(uselst) then begin
ct2lst,rise,obs.longitude*(-1),tz,jdsrise
ct2lst,set,obs.longitude*(-1),tz,jdsset
endif


if keyword_set(print) then begin
if keyword_set(uselst) then print,"Sun times in LST:" else print, "Sun times in PST:"
print,"Set:", set
print,"Rise:", rise
endif


end