dtmath.lisp 6.11 KB
Newer Older
1
2
;; Mathematics on timepoints.
;; Liam Healy Fri Oct 12 2001 - 16:32
3
;; Time-stamp: <2014-10-06 23:10:23EDT dtmath.lisp>
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23

;; Copyright 2011 Liam M. Healy
;; Distributed under the terms of the GNU General Public License
;;
;; This program is free software: you can redistribute it and/or modify
;; it under the terms of the GNU General Public License as published by
;; the Free Software Foundation, either version 3 of the License, or
;; (at your option) any later version.
;;
;; This program is distributed in the hope that it will be useful,
;; but WITHOUT ANY WARRANTY; without even the implied warranty of
;; MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
;; GNU General Public License for more details.
;;
;; You should have received a copy of the GNU General Public License
;; along with this program.  If not, see <http://www.gnu.org/licenses/>.

(in-package :antik)

(export '(timepoint-j2000sec make-timepoint-J2000sec
Liam M. Healy's avatar
Liam M. Healy committed
24
25
	  datime+ datime- days+ days-
	  start-of-day seconds-per-day datime-J2000day))
healy's avatar
healy committed
26
27
28

;;; Mathematics can only be done on timepoints, not on dtspecs,
;;; because the scale must be known first.
29
30
31
32
;;; Keep in mind that CL's *-universal-time are handy for
;;; linearizing and delinearizing a specification to do math,
;;; but they carry no scale; they are simply a dtspec.

33
34
35
(defvar *j2000-epoch* (read-timepoint "2000-01-01T12:00:00"))

(defvar *clut-j2000-tai*
36
    (timeparse-to-clut
37
     (timeparse (convert-time-scale *j2000-epoch* :tai)))
38
  "J2000 epoch in CL's `universal-time' time linearization, in the TAI scale.")
healy's avatar
healy committed
39
40
41
42

(defun timepoint-j2000sec (tp)
  "Seconds from the J2000 epoch, 2000-01-01T12:00:00UTC."
  (flet ((tp-to-clut (tp)
43
	   (timeparse-to-clut (timeparse tp))))
healy's avatar
healy committed
44
    (check-type tp timepoint)
45
    (cl:- (tp-to-clut (convert-time-scale tp :tai)) *clut-j2000-tai*)))
healy's avatar
healy committed
46

47
48
49
50
51
(defun make-timepoint-J2000sec (J2000sec &optional day-only)
  "Make a timepoint from J2000sec."
  (convert-time-scale
   (make-timepoint
    :timeparse
52
    (clut-to-timeparse (cl:+ *clut-j2000-tai* J2000sec))
53
54
55
    :scale :tai
    :day-only day-only)
   :utc))
healy's avatar
healy committed
56
57

;;;****************************************************************************
58
59
60
;;; Whole days
;;;****************************************************************************

61
62
(defun start-of-day
    (datime &optional (relative-day 0) (day-only (day-only datime)))
63
64
  "From the given datime, return the start of the given day,
   relative-day = 0: this day, -1: day before, +1 next day, etc."
65
66
  (let ((answer
	 (dtspec-from-julian-day-number
67
	  (cl:+ (julian-day-number datime t) relative-day)
68
69
70
71
	  t)))
    (change-class answer (class-of datime))
    (if (typep datime 'timepoint)
	(setf (slot-value answer 'scale) (scale datime)))
72
73
    (if (or day-only (day-only datime))
	(setf (slot-value answer 'day-only) day-only))
74
75
76
77
78
    answer))

;;; (seconds-per-day #d1998-12-31) -> 86401
(defun seconds-per-day (datime)
  "Compute the number of seconds in the given day."
79
80
  (pq-magnitude
   (datime- (start-of-day datime 1) (start-of-day datime))))
81
82

(defun datime-J2000day (datime)
83
84
  "Julian Day from J2000 epoch, and seconds from the start of the day."
  (values
85
    ;; Days from J2000 epoch to the start of the UTC day; always integer + 0.5.
86
87
    (cl:- (days- datime (MAKE-TIMEPOINT :TIMEPARSE '(0 0 12 1 1 2000) :DAY-ONLY NIL :SCALE :UTC))
	  1/2)
88
89
90
91
92
    ;; Seconds from the start of the UTC day to the current timepoint.
    ;; No timescale is required because this is
    ;; only a difference from the start of the day.
    (cl:- (timeparse-to-clut (timeparse datime)) ;(timepoint-J2000sec datime)
	  (timeparse-to-clut (timeparse (start-of-day datime))))))
93

94
95
96
97
(defun days- (datime1 datime2)
  "Integer number of days separating the days.  Differs from datime-
   in that the time of day doesn't matter, just the day on which the two 
   datimes occur."
98
  (cl:- (julian-day-number datime1) (julian-day-number datime2)))
99

100
101
102
(defun days+ (datime &optional (increment 1))
  "Add the number of days to the given datime."
  (start-of-day datime increment))
103

104
105
;;;****************************************************************************
;;; Addition, subtraction
healy's avatar
healy committed
106
107
108
109
110
111
;;;****************************************************************************

(defun datime+ (tp time)
  "Add the time to the timepoint.  If time is a number,
   it is assumed to represent seconds."
  (check-type tp timepoint)
112
  (with-units
113
114
115
116
117
      (with-pq ((time time))
	(make-timepoint-j2000sec
	 (cl:+ (timepoint-J2000sec tp)
	       (if (numberp time)
		   time
118
		   (pq-magnitude time)))
119
	 (day-only tp)))))
healy's avatar
healy committed
120

121
122
123
;;; N.B.:
;;; (datime- #d1999-01-01T12:00:00 #d1998-12-31T12:00:00)
;;; 86401.0 second
healy's avatar
healy committed
124
125
(defun datime- (datime1 b)
  "The difference of the datimes, or subtract a time from a datime."
126
  (if (typep b 'timepoint)
127
128
129
      ;; add conditional to check if either are day-only
      ;; if so, just count the number of days.
      (if (or (day-only datime1) (day-only b))
130
	  (make-pq (days- datime1 b) 'day)
131
132
133
134
135
	  (make-pq		 ; time difference between two datimes
	   (cl:- (timepoint-J2000sec datime1)
		 (timepoint-J2000sec b))
	   'second))
      (datime+ datime1 (antik:- b))))	; subtract a time
healy's avatar
healy committed
136
	  
137
(defmethod +i ((a dtspec) b)
healy's avatar
healy committed
138
139
  (datime+ a b))

140
(defmethod +i (a (b dtspec))
healy's avatar
healy committed
141
142
  (datime+ b a))

143
(defmethod -i ((a dtspec) b)
healy's avatar
healy committed
144
145
146
147
148
149
  (datime- a b))

;;;****************************************************************************
;;; Relational
;;;****************************************************************************

150
151
152
153
(defun timepoint-relation (relation x y)
  "Do the two timepoints obey the relation?"
  (datime-relation
   relation
154
155
156
157
158
159
   (multiple-value-bind (converted error)
       (ignore-errors (convert-time-scale y (scale x)))
     (if (typep error 'error)		; can't convert y to x type?
	 ;; then convert x to y type
	 (list (convert-time-scale x (scale y)) y)
       (list x (if (eq (scale x) (scale y)) y converted))))))
healy's avatar
healy committed
160

161
;;; define <, etc.
162
163
164
165
(defmethod <i ((x timepoint) (y timepoint)) (timepoint-relation '< x y))
(defmethod >i ((x timepoint) (y timepoint)) (timepoint-relation '> x y))
(defmethod <=i ((x timepoint) (y timepoint)) (timepoint-relation '<= x y))
(defmethod >=i ((x timepoint) (y timepoint)) (timepoint-relation '>= x y))
166

167
(defmethod =i ((x timepoint) (y timepoint))
168
  (if (timepoint-relation '= x y) x))