time   calls  line
                  1 function jd = date2jd(varargin)
                  2 %DATE2JD Julian day number from Gregorian date.
                  3 %
                  4 %   JD = DATE2JD(YEAR, MONTH, DAY, HOUR, MINUTE, SECOND) returns the Julian
                  5 %   day number of the given date (Gregorian calendar) plus a fractional part
                  6 %   depending on the time of day.
                  7 %
                  8 %   Any missing MONTH or DAY will be replaced by ones.  Any missing HOUR,
                  9 %   MINUTE or SECOND will be replaced by zeros.
                 10 %
                 11 %   If no date is specified, the current date and time is used.
                 12 %
                 13 %   Start of the JD (Julian day) count is from 0 at 12 noon 1 January -4712
                 14 %   (4713 BC), Julian proleptic calendar.  Note that this day count conforms
                 15 %   with the astronomical convention starting the day at noon, in contrast
                 16 %   with the civil practice where the day starts with midnight.
                 17 %
                 18 %   Astronomers have used the Julian period to assign a unique number to
                 19 %   every day since 1 January 4713 BC.  This is the so-called Julian Day
                 20 %   (JD).  JD 0 designates the 24 hours from noon UTC on 1 January 4713 BC
                 21 %   (Julian proleptic calendar) to noon UTC on 2 January 4713 BC.
                 22 
                 23 %   Sources:  - http://tycho.usno.navy.mil/mjd.html
                 24 %             - The Calendar FAQ (http://www.faqs.org)
                 25 
                 26 %   Author:      Peter J. Acklam
                 27 %   Time-stamp:  2002-05-24 13:30:06 +0200
                 28 %   E-mail:      pjacklam@online.no
                 29 %   URL:         http://home.online.no/~pjacklam
                 30 
            13   31    nargsin = nargin; 
            13   32    error(nargchk(0, 6, nargsin)); 
            13   33    if nargsin 
            13   34       argv = {1 1 1 0 0 0}; 
            13   35       argv(1:nargsin) = varargin; 
                 36    else
                 37       argv = num2cell(clock);
                 38    end
            13   39    [year, month, day, hour, minute, second] = deal(argv{:}); 
                 40 
                 41    % The following algorithm is a modified version of the one found in the
                 42    % Calendar FAQ.
                 43 
            13   44    a = floor((14 - month)/12); 
            13   45    y = year + 4800 - a; 
            13   46    m = month + 12*a - 3; 
                 47 
                 48    % For a date in the Gregorian calendar:
            13   49    jd = day + floor((153*m + 2)/5) ... 
                 50         + y*365 + floor(y/4) - floor(y/100) + floor(y/400) - 32045 ...
            13   51         + ( second + 60*minute + 3600*(hour - 12) )/86400;