//----高精曆算數據（許劍偉先生提供的數據代碼，我有改動，因我摻合了天文臺的數據）-----//
function MQC (h) { return  0.0002967/Math.tan( h +0.003138/( h+0.08919) ); } //大氣折射,h是真高度
function MQC2(ho){ return -0.0002909/Math.tan( ho+0.002227/(ho+0.07679) ); } //大氣折射,ho是視高度
function llr2xyz(JW){ //球面轉直角坐標
  var r=new Array(),J=JW[0],W=JW[1],R=JW[2];
  r[0]=R*cos(W)*cos(J);
  r[1]=R*cos(W)*sin(J);
  r[2]=R*sin(W);
  return r;
}
function xyz2llr(xyz){ //直角坐標轉球
  var r=new Array(), x=xyz[0],y=xyz[1],z=xyz[2];
  r[2] = sqrt(x*x+y*y+z*z);
  r[1] = asin(z/r[2]);
  r[0] = rad2mrad( atan2(y,x) );
  return r;
}
function parallax(z,H,fa,high){ //視差修正z赤道坐標,fa地理緯度,H時角,high海拔(千米)
  var dw=1; if(z[2]<500) dw=cs_AU;
  z[2]*=dw;
  var r0,x0,y0,z0, f=cs_ba, u=Math.atan(f*Math.tan(fa)), g=z[0]+H;
  r0 = cs_rEar*Math.cos(u)  + high*Math.cos(fa); //站點與地心向徑的赤道投影長度
  z0 = cs_rEar*Math.sin(u)*f+ high*Math.sin(fa); //站點與地心向徑的軸向投影長度
  x0 = r0*Math.cos(g);
  y0 = r0*Math.sin(g);

  var s = llr2xyz(z);
  s[0]-=x0, s[1]-=y0, s[2]-=z0;
  s = xyz2llr(s);
  z[0]=s[0], z[1]=s[1], z[2]=s[2]/dw;
}
function nutation2(t){ //中精度章動計算,t世紀數
 var i,c,a, t2=t*t, B=this.nutB, dL=0,dE=0;
 for(i=0;i<B.length;i+=5){
  c = B[i]+B[i+1]*t+B[i+2]*t2;
  if(i==0) a=-1.742*t; else a=0;
  dL+=(B[i+3]+a)*Math.sin(c);
  dE+= B[i+4]   *Math.cos(c);
 }
 return [dL/100/rad, dE/100/rad]; //黃經章動,交角章動
}
function mod2(a,b){ //臨界余數(a與最近的整倍數b相差的距離)
  var c=a/b;
  c -= Math.floor(c);
  if(c>0.5) c-=1;
  return c*b;
}
function gxc_sunLat(t) { return 0;       }
function gxc_moonLat(t){ return 0.063*Math.sin(0.057+8433.4662*t+0.000064*t*t)/rad; } //月緯光行差,誤差0.006"
function pGST(T,dt){ //T是2000年首起算的日數(UT),dt是deltatT(日),精度不高時dt可取值為0
  //返回平恒星時,即格林尼治子午圈的平春分點起算的赤經
  var t=(T+dt)/36525,t2=t*t, t3=t2*t, t4=t3*t;
  return pi2*(0.7790572732640 + 1.00273781191135448*T) //T是UT,下行的t是力學時(世紀數)
      + (0.014506 + 4612.15739966*t + 1.39667721*t2 - 0.00009344*t3 + 0.00001882*t4)/rad;
}
function m_coord(t,n1,n2,n3){ //月球坐標,n1,n2,n3為所取的項數
  var re = new Array();
  re[0] = XL1_calc( 0, t, n1 );
  re[1] = XL1_calc( 1, t, n2 );
  re[2] = XL1_calc( 2, t, n3 );
  return re;
}
function e_coord(t,n1,n2,n3){ //地球坐標,t為世紀數
  var re = new Array();
  re[0]= XL0_calc( 0,0, t, n1);
  re[1]= XL0_calc( 0,1, t, n2);
  re[2]= XL0_calc( 0,2, t, n3);
  return re;
}
function rad2strE(d,tim,ext){ //將弧度轉為字串,ext為小數保留位數
 //tim=0例: -23°59' 48.23"tim=1例:  18h 29m 44.52s
 var s=" ";
 var w1="°",w2="′",w3="″";
 if(d<0)  d=-d,s='-';
 if(tim){ d*=12/Math.PI; w1="h ",w2="m ",w3="s " }
 else     d*=180/Math.PI;
 var a=Math.floor(d); d=(d-a)*60;
 var b=Math.floor(d); d=(d-b)*60;
 var c=Math.floor(d);
 var Q = Math.pow(10,ext);
 d=Math.floor( (d-c)*Q +0.5 );
 if(d>=Q)   d-=Q,  c++;
 if(c>=60)  c-=60, b++;
 if(b>=60)  b-=60, a++;
 a=""+a, b="0"+b, c="0"+c, d="00000"+d;
 s+=a.substr(a.length-3,3)+w1;
 s+=b.substr(b.length-2,2)+w2;
 s+=c.substr(c.length-2,2);
 if(ext) s+="."+d.substr(d.length-ext,ext)+w3;
 return s;
}
function rad2str(d,tim){ return rad2strE(d,tim,2); }
var dt_at=new Array( // TD - UT1 計算表
-4000,108371.7,-13036.80,392.000, 0.0000,
 -500, 17201.0,  -627.82, 16.170,-0.3413,
 -150, 12200.6,  -346.41,  5.403,-0.1593,
  150,  9113.8,  -328.13, -1.647, 0.0377,
  500,  5707.5,  -391.41,  0.915, 0.3145,
  900,  2203.4,  -283.45, 13.034,-0.1778,
 1300,   490.1,   -57.35,  2.085,-0.0072,
 1600,   120.0,    -9.81, -1.532, 0.1403,
 1700,    10.2,    -0.91,  0.510,-0.0370,
 1800,    13.4,    -0.72,  0.202,-0.0193,
 1830,     7.8,    -1.81,  0.416,-0.0247,
 1860,     8.3,    -0.13, -0.406, 0.0292,
 1880,    -5.4,     0.32, -0.183, 0.0173,
 1900,    -2.3,     2.06,  0.169,-0.0135,
 1920,    21.2,     1.69, -0.304, 0.0167,
 1940,    24.2,     1.22, -0.064, 0.0031,
 1960,    33.2,     0.51,  0.231,-0.0109,
 1980,    51.0,     1.29, -0.026, 0.0032,
 2000,    63.87,    0.1,   0,     0,
 2005,    64.7,     0.4,   0,     0,   //1次項記為x,則 10x=0.4秒/年*(2015-2005),x=0.4
 2015,    69
);