;
; Copywrite 2004 Peter Mills. All rights reserved.
;
; A class definition for working with dates and times.
;
@time_lib
pro time_class__define
dum={time_class, yr:0, dy:0.D}
end
; Most primitive operators: nleap, isleap, time_class_add_days, time_class_diff
; definition of reference day, time_class_ref_day, for converting to double
; utility for equalizing the number of elements in binary operators, equalize_elements
; The reference day:
function time_class_ref_date
return, {time_class, 1900, 0}
end
; Returns the number of leap years between y1 and y2:
function nleap, y1, y2
y1m1=long(y1)-1
y2m1=long(y2)-1
return, y2m1/4 - y1m1/4 + y1m1/100 - y2m1/100 + y2m1/400 - y1m1/400
end
;returns a 1 if year is leap, 0 if not:
function isleap, year
return, (year mod 400 eq 0 or year mod 4 eq 0) and (year mod 100 ne 0)
end
pro equalize_elements, t1, t2, t1a, t2a
n1=n_elements(t1)
n2=n_elements(t2)
if n1 eq n2 then begin
t1a=t1
t2a=t2
endif else if n1 eq 1 then begin
t1a=replicate(t1, n2)
t2a=t2
endif else if n2 eq 1 then begin
t1a=t1
t2a=replicate(t2, n1)
endif else begin
n=min([n1, n2])
t1a=t1[0:n-1]
t2a=t2[0:n-1]
endelse
end
; Adds d days to the object instance:
pro time_class_add_days, t, d1
;if the days have only one element, add that to all the times
;if the days have fewer elements than the times, keep the time
;array the same size and add only those elements
;if the days have more elements, add only up to the number of times...
n=n_elements(t)
if n_elements(d1) eq 1 then d=replicate(double(d1), n) else d=double(d1)
n1=min([n_elements(d), n])
t[0:n1-1].dy=t[0:n1-1].dy+d[0:n1-1]
ind1= where(t[0:n1-1].yr ne 0 and t[0:n1-1].dy ge 365+isleap(t[0:n1-1].yr), cnt1)
if cnt1 gt 0 then begin
oldyear=t[ind1].yr
t[ind1].yr=t[ind1].yr+long(t[ind1].dy)/365
t[ind1].dy=t[ind1].dy-365.*long(t[ind1].dy/365.)
t[ind1].dy=t[ind1].dy-nleap(oldyear, t[ind1].yr)
endif
ind2=where(t.dy lt 0, cnt2)
while cnt2 gt 0 do begin
t[ind2].yr=t[ind2].yr-1
t[ind2].dy=365+t[ind2].dy+isleap(t[ind2].yr)
ind2=where(t.dy lt 0, cnt2)
endwhile
;stop
end
;returns the number of days elapsed between this and the argument
;as a floating point value:
function time_class_diff, t1a, t2a
equalize_elements, t1a, t2a, t1, t2
return, 365.*(t1.yr-t2.yr)+t1.dy-t2.dy+nleap(t2.yr, t1.yr)
end
;Initialization operators:
;time_class_init initializes from individual fields
;time_class_convert converts from a string or the old time structures
;Given year, month, day, hour, minute, second, initialized an
;instance of the time class:
function time_class_init, y1, mon1, d1, h1, min1, s1, err=err
ny=n_elements(y1)
nmon=n_elements(mon1)
nd=n_elements(d1)
nh=n_elements(h1)
nmin=n_elements(min1)
ns=n_elements(s1)
n=max([ny, nmon, nd, nh, nmin, ns])
;if n ne 1 then begin
if ny eq 1 then y=replicate(y1, n) else y=long(y1)
if nmon eq 1 then mon=replicate(mon1, n) else mon=long(mon1)
if nd eq 1 then d=replicate(d1, n) else d=long(d1)
if nh eq 1 then h=replicate(h1, n) else h=long(h1)
if nmin eq 1 then min=replicate(min1, n) else min=long(min1)
if ns eq 1 then s=replicate(double(s1), n) else s=double(s1)
;endif
err=0
if ny ne n or nmon ne n or nd ne n or nh ne n or nmin ne n or ns ne n then begin
message, "Year, month, day, hour, minute and seconds fields should have same number of elements", /continue
n=min([ny, nmon, nd, nh, nmin, ns])
y=y[0:n-1]
mon=mon[0:n-1]
d=d[0:n-1]
h=h[0:n-1]
min=min[0:n-1]
s=s[0:n-1]
err=1
;return, null
endif
ind=where(mon lt 0 or mon gt 12, cnt)
if cnt gt 0 then begin
message, "Month field out of bounds", /continue
err=2
endif
;ind=where(day lt 1 or day gt 31, cnt)
;if cnt gt 0 then begin
; message, "Day field out of bounds", /continue
; err=3
;endif
ind=where(h lt 0 or h gt 59, cnt)
if cnt gt 0 then begin
message, "Hour field out of bounds", /continue
err=4
endif
ind=where(min lt 0 or min gt 59, cnt)
if cnt gt 0 then begin
message, "Minute field out of bounds", /continue
err=5
endif
ind=where(s lt 0 or s gt 60, cnt)
if cnt gt 0 then begin
message, "Second field out of bounds", /continue
err=6
endif
month_days=[0, 31, 59, 90, 120, 151, 181, 212, 243, 273, 304, 334]
ind1=where(d ne 0, cnt1)
if cnt1 gt 0 then d[ind1]=d[ind1]-1
t=replicate({time_class}, n)
t.yr=y+(mon-1) / 12
mon=(mon-1) mod 12 + 1
d=d+month_days[mon-1]
dd=d+((s/60.D0+min)/60.D0+h)/24.D0
t.dy=0
time_class_add_days, t, dd
;stop
ind2=where(isleap(t.yr) and mon gt 2, cnt2)
if cnt2 gt 0 then t[ind2].dy=t[ind2].dy+1
;stop
return, t
end
function time_class_convert, t1, ref_date=ref_date
type=size(t1, /type)
if type eq 8 then begin
tag_name=tag_names(t1, /struct)
if tag_name eq "TIME_STR" then begin
return, time_class_init(t1.year, t1.month, t1.day, t1.hour, t1.minute, t1.second)
endif else if tag_name eq "TIME_CLASS" then begin
return, t1
endif else begin
message, "Unknown structure type"
return, null
endelse
endif else if type eq 7 then begin
t2=convert_date(t1)
return, time_class_init(t2.year, t2.month, t2.day, t2.hour, t2.minute, t2.second)
endif else begin
if n_elements(ref_date) eq 0 then ref_date=time_class_ref_date()
if n_elements(ref_date) eq 1 then ref_date=replicate(ref_date, n_elements(t1))
time_class_add_days, ref_date, t1
return, ref_date
endelse
end
;Operators for extracting fields:
;time_class_get_fields, time_class_doy, time_class_year, time_class_month, time_class_day
;time_class_hour, time_class_minute, time_class_second and, best of all, time_class_dow
; Returns all the fields of the date. Two versions--one for longword and one for
; short.
pro time_class_get_fields, t, y, mon, d, h, min, s
month_days=[0, 31, 59, 90, 120, 151, 181, 212, 243, 273, 304, 334]
y=t.yr;
dayfrac=t.dy-long(t.dy)
dayfrac=dayfrac*24.
h=long(dayfrac)
dayfrac=(dayfrac-h)*60
min=long(dayfrac)
s=(dayfrac-min)*60
ind1=where(y eq 0, cnt1, compl=ind2, ncompl=cnt2)
;n=n_elements(t)
;d=lonarr(n)
;mon=lonarr(n)
dim=size(t, /dim)
d=make_array(dim=size(t, /dim), /long)
mon=make_array(dim=size(t, /dim), /long)
if cnt1 gt 0 then begin
d[ind1]=long(t[ind1].dy)
mon[ind1]=0
endif
if cnt2 gt 0 then begin
dy_int=long(t[ind2].dy)
mon[ind2]=value_locate(month_days, dy_int)+1
d[ind2]=1
ind3=where(isleap(t[ind2].yr), cnt3)
if cnt3 gt 0 then begin
ind4=where(mon[ind2[ind3]] gt 2, cnt4)
if cnt4 gt 0 then begin
ind5=where(dy_int[ind3[ind4]] eq month_days[mon[ind2[ind3[ind4]]]-1], cnt5)
if cnt5 gt 0 then begin
mon[ind2[ind3[ind4[ind5]]]]=mon[ind2[ind3[ind4[ind5]]]]-1
endif
ind6=where(dy_int[ind3[ind4]] ne 59, cnt6)
if cnt6 gt 0 then d[ind2[ind3[ind4[ind6]]]]=0 ;fuck that's a lot of indirection...
endif
endif
d[ind2]=d[ind2]+dy_int-month_days[mon[ind2]-1]
endif
;stop
end
;returns the day of the year:
function time_class_doy, t1
return, t1.dy;
end
;returns the year field:
function time_class_year, t1
return, t1.yr;
end
; Returns the month field:
function time_class_month, t
month_days=[0, 31, 59, 90, 120, 151, 181, 212, 243, 273, 304, 334]
;mon=lonarr(n_elements(t))
mon=make_array(dim=size(t, /dim), /long)
ind1=where(t.yr eq 0, cnt1, compl=ind2, ncompl=cnt2)
if cnt1 gt 0 then begin
mon[ind1]=0
endif
if cnt2 gt 0 then begin
dy_int=long(t[ind2].dy)
mon[ind2]=value_locate(month_days, long(t[ind2].dy))+1
ind3=where(isleap(t[ind2].yr), cnt3)
if cnt3 gt 0 then begin
ind4=where(mon[ind2[ind3]] gt 2 and dy_int[ind3] eq month_days[mon[ind2[ind3]]-1], cnt4)
if cnt4 gt 0 then mon[ind2[ind3[ind4]]]=mon[ind2[ind3[ind4]]]-1
endif
endif
return, mon
end
; Returns the day-of-the-month field:
function time_class_day, t
month_days=[0, 31, 59, 90, 120, 151, 181, 212, 243, 273, 304, 334]
ind1=where(t.yr eq 0, cnt1, compl=ind2, ncompl=cnt2)
;d=lonarr(n_elements(t))
d=make_array(dim=size(t, /dim), /long)
if cnt1 gt 0 then begin
d[ind1]=t[ind1].dy
endif
if cnt2 gt 0 then begin
dy_int=long(t[ind2].dy)
mon=value_locate(month_days, dy_int)+1
d[ind2]=1
ind3=where(isleap(t[ind2].yr), cnt3)
if cnt3 gt 0 then begin
ind4=where(mon[ind3] gt 2, cnt4)
if cnt4 gt 0 then begin
ind5=where(dy_int[ind3[ind4]] eq month_days[mon[ind3[ind4]]-1], cnt5)
if cnt5 gt 0 then mon[ind3[ind4[ind5]]]=mon[ind3[ind4[ind5]]]-1
ind6=where(dy_int[ind3[ind4]] ne 59, cnt6)
if cnt6 gt 0 then d[ind2[ind3[ind4[ind6]]]]=0
endif
endif
d[ind2]=d[ind2]+dy_int-month_days[mon-1]
endif
;stop
return, d
end
; Returns the hour field:
function time_class_hour, t
dayfrac=t.dy-long(t.dy)
dayfrac=dayfrac*24.;
h=long(dayfrac)
return, h
end
; Returns the minutes field:
function time_class_minute, t
dayfrac=t.dy-long(t.dy)
dayfrac=dayfrac*24.
h=long(dayfrac)
dayfrac=(dayfrac-h)*60
min=long(dayfrac)
return, min
end
; Returns the seconds field:
function time_class_second, t
min=(t.dy-long(t.dy))*24*60;
s=(min-long(min))*60;
return, s
end
;returns the day of the week:
function time_class_dow, t
dd=time_class_double(t)
ind1=where(dd ge 0, cnt1, compl=ind2, ncompl=cnt2)
dow=lonarr(n_elements(t))
if cnt1 gt 0 then begin
dow[ind1]=long(dd[ind1]) mod 7 + 1
endif
if cnt2 gt 0 then begin
dow[ind2]=7+long(floor(dd[ind2])) mod 7
endif
;stop
return, dow
end
;type conversion: time_class_double, time_class_time_str, time_class_string
;converts to a floating point. Gives the number of days elapsed since ref_date.
function time_class_double, t
ind1=where(t.yr ne 0, cnt1, compl=ind2, ncompl=cnt2)
n=n_elements(t)
result=dblarr(n)
if cnt1 gt 0 then begin
result[ind1]=time_class_diff(t[ind1], time_class_ref_date())
endif
if cnt2 gt 0 then begin
result[ind2]=t[ind2].dy
endif
;stop
return, result
end
;converts to the old style time structures:
function time_class_time_str, t1
time_class_get_fields, t1, year, month, day, hour, minute, second
;t2=replicate({time_str}, size(t1, /dim))
t2=make_array(value={time_str}, dim=size(t1, /dim))
t2.year=year
t2.month=month
t2.day=day
t2.hour=hour
t2.minute=minute
t2.second=second
return, t2
end
; Writes an object instance of type time to a string:
function time_class_string, t1
;do this the easy way by extracting the fields, converting to the old time structures
;and then writing them to a string
t2=time_class_time_str(t1)
return, time_string(t2);
end
;binary operators: time_class_compare, time_class_add, time_class_subtract, time_class_multiply,
; time_class_divide
; compare two time classes:
function time_class_compare, t1, t2
if n_elements(t1) eq 1 or n_elements(t2) eq 1 then begin
if t1[0].yr gt t2[0].yr then return, -1
if t2[0].yr gt t1[0].yr then return, 1
if t1[0].dy gt t2[0].dy then return, -1
if t2[0].dy gt t1[0].dy then return, 1
return, 0
endif else begin
return, (t2.yr gt t1.yr or t2.dy gt t2.dy)-(t1.yr gt t2.yr or t1.dy gt t2.dy)
endelse
end
; Addition operator. May not be well defined if the second argument has
; a year field.
function time_class_add, t1a, t2a, err=ind
equalize_elements, t1a, t2a, t1, t2
;n=min([n_elements(t1), n_elements(t2)])
n=n_elements(t1)
ind=where(t1.dy ne 0 and t1.yr ne 0 and t2.yr ne 0 and t2.dy ne 0, cnt)
if cnt gt 0 then begin
t1str=time_class_write_string(t1[ind])
t2str=time_class_write_string(t2[ind])
message, "Addition of time classes:", /continue
message, strjoin(t1str, ", ")+" with:", /continue
message, strjoin(t2str, ", "), /continue
message, "is ambiguous", /continue
endif
newt=replicate({time_class}, n)
newt.yr=t1.yr+t2.yr
newt.dy=t1.dy
time_class_add_days, newt, t2.dy
return, newt;
end
; Subtraction operator. The returned time will not have a year field.
function time_class_subtract, t1a, t2a
equalize_elements, t1a, t2a, t1, t2
;n=min([n_elements(t1), n_elements(t2)])
n=n_elements(t1)
newt=replicate({time_class}, n)
ind1=where(t2.yr ne 0, cnt1, compl=ind2, ncompl=cnt2)
if cnt1 gt 0 then begin
newt[ind1].dy=0
newt[ind1].dy=time_class_diff(t1[ind1], t2[ind1])
endif
if cnt2 gt 0 then begin
newt[ind2].yr=t1[ind2].yr
newt[ind2].dy=t1[ind2].dy
temp=newt[ind2]
time_class_add_days, temp, -t2[ind2].dy
newt[ind2]=temp
endif
return, newt
end
; Multiplication operator. May not be well defined if the first argument
; has a year field.
function time_class_multiply, t1a, k1
equalize_elements, t1a, double(k1), t1, k
;n=min([n_elements(t1), n_elements(t2)])
n=n_elements(k)
newt=replicate({time_class}, n)
ind1=where(t1.yr ne 0 and t1.dy ne 0, cnt1)
ind2=where(t1.yr eq 0, cnt2)
ind3=where(t1.dy eq 0, cnt3)
if cnt1 gt 0 then begin
t1str=strjoin(time_class_write_string(t1[ind1]), ", ")
str2=strjoin(string(k[ind1]), ", ")
message, "Multiplication of time classes:", /continue
message, t1str, " with:", /continue
message, str2, /continue
message, "is ambiguous", /continue
newt[ind1].dy=time_class_double(t1[ind1])*k[ind1]
endif
if cnt2 gt 0 then begin
newt[ind2].dy=t1[ind2].dy*k[ind2]
endif
if cnt3 gt 0 then begin
newt[ind3].yr=t1[ind3].yr*k[ind3]
endif
;stop
return, newt
end
; Division operator. May not be well defined if the first argument has a
; year field.
function time_class_divide, t1a, t2a
equalize_elements, t1a, t2a, t1, t2
type=size(t2, /type)
;n=min([n_elements(t1), n_elements(t2)])
n=n_elements(t1)
if type ne 8 then begin
t2=double(t2)
newt=replicate({time_class}, n)
ind1=where(t1.yr ne 0 and t1.dy ne 0, cnt1)
ind2=where(t1.yr eq 0, cnt2)
ind3=where(t1.dy eq 0, cnt3)
if cnt1 gt 0 then begin
t1str=strjoin(time_class_write_string(t1[ind1]), ", ")
str2=strjoin(string(t2[ind1]), ", ")
message, "Division of time classes:", /continue
message, t1str, " with:", /continue
message, str2, /continue
message, "is ambiguous", /continue
newt[ind1].dy=time_class_double(t1[ind1])/t2[ind1]
endif
if cnt2 gt 0 then begin
newt[ind2].dy=t1[ind2].dy/t2[ind2]
endif
if cnt3 gt 0 then begin
newt[ind3].yr=t1[ind3].yr/t2[ind3]
endif
return, newt
endif else begin
; Division operator. May not be well defined if the second field has a
; year field.
result=dblarr(n)
flag1=t1.yr eq 0 and t2.yr eq 0
flag2=t1.dy eq 0 and t2.dy eq 0
ind1=where(flag1 eq 0 and flag2 eq 0, cnt1)
ind2=where(flag1, cnt2)
ind3=where(flag2, cnt3)
if cnt1 gt 0 then begin
t1str=strjoin(time_class_write_string(t1[ind1]), ", ")
t2str=strjoin(time_class_write_string(t2[ind1]), ", ")
message, "Division of time classes:", /continue
message, t1str, " with:", /continue
message, t2str, /continue
message, "is ambiguous", /continue
result[ind1].dy=time_class_double(t1[ind1])/time_class_double(t2[ind1])
endif
if cnt2 gt 0 then begin
result[ind2]=t1[ind2].dy/t2[ind2].dy
endif
if cnt3 gt 0 then begin
result[ind3]=t1[ind3].yr/t2[ind2].yr
endif
return, result
endelse
end