function [decyear,yearday,year] = invjul(juliandays,str)

%
% Calcuate decimal year and/or (yearday, year) from julian day
% only after 1969.
%
% Usage:
% [decyear] = invjul(juliandays, 1)
% [yearday, year] = invjul(juliandays, 2)
% [decyear, yearday, year] = invjul(juliandays, 3)
%
% str == 1, returns decyear only
% str == 2, returns in yearday and year format
% str == 3, returns all decyear, yearday and year
%

% $Log$
%

% $Id$
% klee 09/01/99

decyear = []; yearday = []; year = [];

if nargin == 1 
  % if no str is specified, return decyear only
  decyear = invjul(juliandays, 1);
  return
end

% days since 1969
if juliandays > 2440000
  % I don't know why, but most of the program create one more julian day 
  % than jul.m, so I substract 222 instead of 221.
  days = juliandays - 2440000 - 222;  
end

leapy = fix(days/1461);
dayleft = round((days/1461 - leapy) * 1461);
%if leapy > 0
%  dayleft = dayleft - 1461*leapy;
%end
count = 0;
while dayleft > 0
  if count < 3
    dayleft = dayleft - 365;
    count = count + 1;
  else 
    dayleft = dayleft - 366;
    count = count + 1;
  end
end
isleapyear = 0;
if count < 4
  if dayleft == 0
    dayleft = 365;
  else
    dayleft = dayleft + 365;
  end
else
  isleapyear = 1;
  if dayleft == 0
    dayleft = 366;
  else
    dayleft = dayleft + 366;
  end
end
yearday = dayleft;
year = 1969 + count + leapy*4 - 1;
if isleapyear == 1
  decyear = year + yearday / 366;
else
  decyear = year + yearday / 365;
end

if str == 1
  year = []; yearday = [];
elseif str == 2
  decyear = [];
end

return
  
  