umoon.pas

支持公历、农历及公历转农历使用公式法（不是查表法）的日历控件

  d2:=EncodeDate(1500,2,28);
  julian_calendar_before_1582:=((d1-d2)=2);
  end;
(*@\\\0000001107*)
(*@/// function EncodedateCorrect(year,month,day: word):TDateTime; *)
function EncodedateCorrect(year,month,day: word):TDateTime;
begin
  result:=delphi_date(Calc_Julian_date(year,month,day));
  end;
(*@\\\*)
(*@/// procedure DecodedateCorrect(date:TDateTime; var year,month,day: word); *)
procedure DecodedateCorrect(date:TDateTime; var year,month,day: word);
begin
  Calc_Calendar_date(julian_date(date),year,month,day);
  end;
(*@\\\*)
(*@/// procedure DecodetimeCorrect(date:TDateTime; var hour,min,sec,msec: word); *)
procedure DecodetimeCorrect(date:TDateTime; var hour,min,sec,msec: word);
begin
  Decodetime(1+frac(date),hour,min,sec,msec);
  end;
(*@\\\*)
(*@/// function FalsifyTdateTime(date:TDateTime):TdateTime; *)
function FalsifyTdateTime(date:TDateTime):TdateTime;
var
  d,m,y: word;
begin
  DecodedateCorrect(date,d,m,y);
  result:=Encodedate(d,m,y);
  result:=result+frac(date);
  if negative_dates_broken and (result<0) and (frac(result)<>0) then
    result:=int(result)-(1-abs(frac(result)));
  end;
(*@\\\*)

{ Calendar functions }
(*@/// function WeekNumber(date:TDateTime):integer; *)
function WeekNumber(date:TDateTime):integer;
var
  y,m,d: word;
  h: integer;
  FirstofJanuary,
  FirstThursday,
  FirstWeekStart: TDateTime;
begin
  DecodedateCorrect(date,y,m,d);
  FirstofJanuary:=EncodedateCorrect(y,1,1);
  h:=dayOfWeek(FirstofJanuary);
  FirstThursday:=FirstofJanuary+((12-h) mod 7);
  FirstWeekStart:=FirstThursday-3;
  if trunc(date)<FirstWeekStart then
    result:=WeekNumber(FirstofJanuary-1) (* 12-31 of previous year *)
  else
    result:=(round(trunc(date)-FirstWeekStart) div 7)+1;
  end;
(*@\\\*)
(*@/// function EasterDateGregorian(year:integer):TDateTime; *)
function EasterDateGregorian(year:integer):TDateTime;
var
  a,b,c,d,e,m,n,day,month: integer;
begin
  case year of
    1583..1699:  begin  m:=22; n:=2;  end;
    1700..1799:  begin  m:=23; n:=3;  end;
    1800..1899:  begin  m:=23; n:=4;  end;
    1900..2099:  begin  m:=24; n:=5;  end;
    2100..2199:  begin  m:=24; n:=6;  end;
    2200..2399:  begin  m:=25; n:=0;  end;
    else raise E_OutOfAlgorithmRange.Create('Out of range of the algorithm');
    end;
  a:=year mod 19;
  b:=year mod 4;
  c:=year mod 7;
  d:=(19*a+m) mod 30;
  e:=(2*b+4*c+6*d+n) mod 7;
  day:=(22+d+e);
  if day<=31 then
    month:=3
  else begin
    day:=(d+e-9);
    month:=4;
    end;
  if (day=26) and (month=4) then  day:=19;
  if (day=25) and (month=4) and (d=28) and (e=6) and (a>10) then  day:=18;
  result:=EncodedateCorrect(year,month,day);
  end;
(*@\\\*)
(*@/// function EasterDate(year:integer):TDateTime; *)
function EasterDate(year:integer):TDateTime;
begin
  if year<1583 then
    result:=EasterDateJulian(year)
  else
    result:=EasterDateGregorian(year);
  end;
(*@\\\*)
(*@/// function EasterDateJulian(year:integer):TDateTime; *)
function EasterDateJulian(year:integer):TDateTime;
var
  a,b,c,d,e,f,g: integer;
begin
  a:=year mod 4;
  b:=year mod 7;
  c:=year mod 19;
  d:=(19*c+15) mod 30;
  e:=(2*a+4*b-d+34) mod 7;
  f:=(d+e+114) div 31;
  g:=(d+e+114) mod 31;
  result:=EncodedateCorrect(year,f,g+1);
  end;
(*@\\\*)
(*@/// function PesachDate(year:integer):TDateTime; *)
function PesachDate(year:integer):TDateTime;
var
  a,b,c,d,j,s: integer;
  q,r: extended;
begin
  if year<359 then
    raise E_OutOfAlgorithmRange.Create('Out of range of the algorithm');
  c:=year div 100;
  if year<1583 then
    s:=0
  else
    s:=(3*c-5) div 4;
  a:=(12*year+12) mod 19;
  b:=year mod 4;
  q:=-1.904412361576+1.554241796621*a+0.25*b-0.003177794022*year+s;
  j:=(trunc(q)+3*year+5*b+2-s) mod 7;
  r:=frac(q);
  if false then
  else if j in [2,4,6] then
    d:=trunc(q)+23
  else if (j=1) and (a>6) and (r>=0.632870370) then
    d:=trunc(q)+24
  else if (j=0) and (a>11) and (r>=0.897723765) then
    d:=trunc(q)+23
  else
    d:=trunc(q)+22;

  if d>31 then
    result:=EncodedateCorrect(year,4,d-31)
  else
    result:=EncodedateCorrect(year,3,d);
  end;
(*@\\\*)
(*@/// function JewishYearStyle(year:word):TJewishYearStyle; *)
function JewishYearStyle(year:word):TJewishYearStyle;
var
  i: TJewishYearStyle;
  yearlength: integer;
begin
  yearlength:=round(pesachdate(year-3760)-pesachdate(year-3761));
  result:=low(TJewishYearStyle);
  for i:=low(TJewishYearStyle) to high(TJewishYearStyle) do
    if yearlength=Jewish_year_length[i] then
      result:=i;
  end;
(*@\\\*)
(*@/// function EncodeDateJewish(year,month,day: word):TDateTime; *)
function EncodeDateJewish(year,month,day: word):TDateTime;
var
  yearstyle: TJewishYearStyle;
  offset,i: integer;
begin
  yearstyle:=JewishYearStyle(year);
  if (month<1) or (month>13) then
    raise EConvertError.Create('Invalid month');
  if (month=13) and
     (yearstyle in [ys_common_deficient,ys_common_regular,ys_common_complete]) then
    raise EConvertError.Create('Invalid month');
  if (day<1) or (day>Jewish_Month_length[month,yearstyle]) then
    raise EConvertError.Create('Invalid day');
  offset:=day-1;
  (* count months from tishri *)
  month:=(month+6) mod 13 +1;
  for i:=1 to month-1 do
    offset:=offset+Jewish_Month_length[(i+5) mod 13 +1,yearstyle];
  result:=pesachdate(year-3761)+163+offset;
  end;
(*@\\\*)
(*@/// procedure DecodeDateJewish(date: TDateTime; var year,month,day: word); *)
procedure DecodeDateJewish(date: TDateTime; var year,month,day: word);
var
  year_g,month_g,day_g: word;
  yearstyle: TJewishYearStyle;
  tishri1: TDateTime;
begin
  DecodedateCorrect(date,year_g,month_g,day_g);
  tishri1:=pesachdate(year_g)+163;
  if tishri1>date then begin
    tishri1:=pesachdate(year_g-1)+163;
    year:=year_g+3760;
    end
  else
    year:=year_g+3761;
  yearstyle:=JewishYearStyle(year);
  month:=7;
  day:=round(date-tishri1+1);
  while day>Jewish_Month_length[month,yearstyle] do begin
    dec(day,Jewish_Month_length[month,yearstyle]);
    month:=(month mod 13) +1;
    end;
  end;
(*@\\\*)

{ Misc }
(*@/// procedure calc_epsilon_phi(date:TDateTime; var delta_phi,epsilon:extended); *)
procedure calc_epsilon_phi(date:TDateTime; var delta_phi,epsilon:extended);
(*$ifndef low_accuracy *)
const
  (*@/// arg_mul:array[0..30,0..4] of shortint = (..); *)
  arg_mul:array[0..30,0..4] of shortint = (
     ( 0, 0, 0, 0, 1),
     (-2, 0, 0, 2, 2),
     ( 0, 0, 0, 2, 2),
     ( 0, 0, 0, 0, 2),
     ( 0, 1, 0, 0, 0),
     ( 0, 0, 1, 0, 0),
     (-2, 1, 0, 2, 2),
     ( 0, 0, 0, 2, 1),
     ( 0, 0, 1, 2, 2),
     (-2,-1, 0, 2, 2),
     (-2, 0, 1, 0, 0),
     (-2, 0, 0, 2, 1),
     ( 0, 0,-1, 2, 2),
     ( 2, 0, 0, 0, 0),
     ( 0, 0, 1, 0, 1),
     ( 2, 0,-1, 2, 2),
     ( 0, 0,-1, 0, 1),
     ( 0, 0, 1, 2, 1),
     (-2, 0, 2, 0, 0),
     ( 0, 0,-2, 2, 1),
     ( 2, 0, 0, 2, 2),
     ( 0, 0, 2, 2, 2),
     ( 0, 0, 2, 0, 0),
     (-2, 0, 1, 2, 2),
     ( 0, 0, 0, 2, 0),
     (-2, 0, 0, 2, 0),
     ( 0, 0,-1, 2, 1),
     ( 0, 2, 0, 0, 0),
     ( 2, 0,-1, 0, 1),
     (-2, 2, 0, 2, 2),
     ( 0, 1, 0, 0, 1)
                   );
  (*@\\\*)
  (*@/// arg_phi:array[0..30,0..1] of longint = (); *)
  arg_phi:array[0..30,0..1] of longint = (
     (-171996,-1742),
     ( -13187,  -16),
     (  -2274,   -2),
     (   2062,    2),
     (   1426,  -34),
     (    712,    1),
     (   -517,   12),
     (   -386,   -4),
     (   -301,    0),
     (    217,   -5),
     (   -158,    0),
     (    129,    1),
     (    123,    0),
     (     63,    0),
     (     63,    1),
     (    -59,    0),
     (    -58,   -1),
     (    -51,    0),
     (     48,    0),
     (     46,    0),
     (    -38,    0),
     (    -31,    0),
     (     29,    0),
     (     29,    0),
     (     26,    0),
     (    -22,    0),
     (     21,    0),
     (     17,   -1),
     (     16,    0),
     (    -16,    1),
     (    -15,    0)
    );
  (*@\\\*)
  (*@/// arg_eps:array[0..30,0..1] of longint = (); *)
  arg_eps:array[0..30,0..1] of longint = (
     ( 92025,   89),
     (  5736,  -31),
     (   977,   -5),
     (  -895,    5),
     (    54,   -1),
     (    -7,    0),
     (   224,   -6),
     (   200,    0),
     (   129,   -1),
     (   -95,    3),
     (     0,    0),
     (   -70,    0),
     (   -53,    0),
     (     0,    0),
     (   -33,    0),
     (    26,    0),
     (    32,    0),
     (    27,    0),
     (     0,    0),
     (   -24,    0),
     (    16,    0),
     (    13,    0),
     (     0,    0),
     (   -12,    0),
     (     0,    0),
     (     0,    0),
     (   -10,    0),
     (     0,    0),
     (    -8,    0),
     (     7,    0),
     (     9,    0)
    );
  (*@\\\*)
(*$endif *)
var
  t,omega: extended;
(*$ifdef low_accuracy *)
  l,ls: extended;
(*$else *)
  d,m,ms,f,s: extended;
  i: integer;
(*$endif *)
  epsilon_0,delta_epsilon: extended;
begin
  t:=(julian_date(date)-2451545.0)/36525;

  (* longitude of rising knot *)
  omega:=put_in_360(125.04452+(-1934.136261+(0.0020708+1/450000*t)*t)*t);

(*$ifdef low_accuracy *)
  (*@/// delta_phi and delta_epsilon - low accuracy *)
  (* mean longitude of sun (l) and moon (ls) *)
  l:=280.4665+36000.7698*t;
  ls:=218.3165+481267.8813*t;

  (* correction due to nutation *)
  delta_epsilon:=9.20*cos_d(omega)+0.57*cos_d(2*l)+0.10*cos_d(2*ls)-0.09*cos_d(2*omega);

  (* longitude correction due to nutation *)
  delta_phi:=(-17.20*sin_d(omega)-1.32*sin_d(2*l)-0.23*sin_d(2*ls)+0.21*sin_d(2*omega))/3600;
  (*@\\\*)
(*$else *)
  (*@/// delta_phi and delta_epsilon - higher accuracy *)
  (* mean elongation of moon to sun *)
  d:=put_in_360(297.85036+(445267.111480+(-0.0019142+t/189474)*t)*t);

  (* mean anomaly of the sun *)
  m:=put_in_360(357.52772+(35999.050340+(-0.0001603-t/300000)*t)*t);

  (* mean anomly of the moon *)
  ms:=put_in_360(134.96298+(477198.867398+(0.0086972+t/56250)*t)*t);

  (* argument of the latitude of the moon *)
  f:=put_in_360(93.27191+(483202.017538+(-0.0036825+t/327270)*t)*t);

  delta_phi:=0;
  delta_epsilon:=0;

  for i:=0 to 30 do begin
    s:= arg_mul[i,0]*d
       +arg_mul[i,1]*m
       +arg_mul[i,2]*ms
       +arg_mul[i,3]*f
       +arg_mul[i,4]*omega;
    delta_phi:=delta_phi+(arg_phi[i,0]+arg_phi[i,1]*t*0.1)*sin_d(s);
    delta_epsilon:=delta_epsilon+(arg_eps[i,0]+arg_eps[i,1]*t*0.1)*cos_d(s);
    end;

  delta_phi:=delta_phi*0.0001/3600;
  delta_epsilon:=delta_epsilon*0.0001/3600;
  (*@\\\*)
(*$endif *)

  (* angle of ecliptic *)
  epsilon_0:=84381.448+(-46.8150+(-0.00059+0.001813*t)*t)*t;

  epsilon:=(epsilon_0+delta_epsilon)/3600;
  end;
(*@\\\0000000A0A*)
(*@/// function star_time(date:TDateTime):extended;            // degrees *)
function star_time(date:TDateTime):extended;
var
  jd, t: extended;
  delta_phi, epsilon: extended;
begin
  jd:=julian_date(date);
  t:=(jd-2451545.0)/36525;
  calc_epsilon_phi(date,delta_phi,epsilon);
  result:=put_in_360(280.46061837+360.98564736629*(jd-2451545.0)+
                     t*t*(0.000387933-t/38710000)+
                     delta_phi*cos_d(epsilon) );
  end;
(*@\\\*)

{ Coordinate functions }
(*@/// procedure calc_geocentric(var coord:t_coord; date:TDateTime); *)
{ Based upon Chapter 13 (12) and 22 (21) of Meeus }

procedure calc_geocentric(var coord:t_coord; date:TDateTime);
var
  epsilon: extended;
  delta_phi: extended;
  alpha,delta: extended;
begin
  calc_epsilon_phi(date,delta_phi,epsilon);
  coord.longitude:=put_in_360(coord.longitude+delta_phi);

  (* geocentric coordinates *)
{   alpha:=arctan2_d(cos_d(epsilon)*sin_d(o),cos_d(o)); }
{   delta:=arcsin_d(sin_d(epsilon)*sin_d(o)); }
  alpha:=arctan2_d( sin_d(coord.longitude)*cos_d(epsilon)
                   -tan_d(coord.latitude)*sin_d(epsilon)
                  ,cos_d(coord.longitude));
  delta:=arcsin_d( sin_d(coord.latitude)*cos_d(epsilon)
                  +cos_d(coord.latitude)*sin_d(epsilon)*sin_d(coord.longitude));

  coord.rektaszension:=alpha;
  coord.declination:=delta;
  end;
(*@\\\0000000129*)
(*@/// procedure calc_horizontal(var coord:t_coord; date:TDateTime longitude,latitude: extended); *)
procedure calc_horizontal(var coord:t_coord; date:TDateTime; longitude,latitude: extended);
var
  h: extended;
begin
  h:=put_in_360(star_time(date)-coord.rektaszension-longitude);
  coord.azimuth:=arctan2_d(sin_d(h),
                           cos_d(h)*sin_d(latitude)-
                           tan_d(coord.declination)*cos_d(latitude) );
  coord.elevation:=arcsin_d(sin_d(latitude)*sin_d(coord.declination)+
                            cos_d(latitude)*cos_d(coord.declination)*cos_d(h));
  end;