Astronomy: Effemeridi (Ephemeris)

[salvomic]

here there is an ancillary program to Effemeridi (Astrolabio) and Astro Lab 4, with some tools to get interpolation (3, 5 points and Lagrange).

The tools are:
interpolation of 3 or 5 values
finding the "n" interpolation factor
finding the extremum (minimum, maximum of a list to interpolate)
finding the zero (3 or 5 values)
"halves" (4 values)
get a value with Lagrange interpolation.

Typical example: interpolation({1.9556788, 1.97345409, 1.965049586}, 0.18), nzero3({1.9556788, 1.97345409, 1.965049586}) ...
Very useful for RA (right ascension) and Declination (or lambda, beta: longitude, latitude) for an aster or planet or to find conjunction of planets, least angular distance and so on...

Salvo

Code:


// Ancillary program for Astro Lab 4 (Marcel Pelletier) and Effemeridi - Astrolabio (Salvo Micciché)
// lista = a list like {data1, data2, data3...}
// n = interpolation factor

EXPORT interpolation3(lista,n)
BEGIN
  LOCAL a,b,c,R,dif1,dif2;
  dif1:=ΔLIST(lista);
  dif2:=ΔLIST(dif1);
  a:=dif1(1);
  b:=dif1(2);
  c:=dif2(1);
  R:=lista(2)+(n/2)*(a+b+n*c);
  RETURN R;
END;

EXPORT interpolation5(lista,n)
BEGIN
  LOCAL a,b,c,R,dif1,dif2, dif3, dif4;
  LOCAL d,e,f,g,h,j,k;
  dif1:=ΔLIST(lista);
  dif2:=ΔLIST(dif1);
  dif3:=ΔLIST(dif2);
  dif4:=ΔLIST(dif3);
  a:=dif1(1);
  b:=dif1(2);
  c:=dif1(3);
  d:=dif1(4);
  e:=dif2(1);
  f:=dif2(2);
  g:=dif2(3);
  h:=dif3(1);
  j:=dif3(2);
  k:=dif4(1);
  R:=lista(3)+(n/2)*(b+c)+(n^2/2)*f+n*((n^2-1)/12)*(h+j)+(n^2*(n^2-1)/24)*k;
  RETURN R;
END;

EXPORT extremum3(lista)
BEGIN
  LOCAL a,b,c,R,dif1,dif2, ym, nm;
  dif1:=ΔLIST(lista);
  dif2:=ΔLIST(dif1);
  a:=dif1(1);
  b:=dif1(2);
  c:=dif2(1);
  ym:=lista(2)-((a+b)^2)/(8*c);
  nm:= -(a+b)/(2*c);
  RETURN ({ym, nm});
END;

EXPORT extremum5(lista)
BEGIN
  LOCAL a,b,c,R,dif1,dif2, dif3, dif4;
  LOCAL d,e,f,g,h,j,k, nm;
  LOCAL jj, temp;
  dif1:=ΔLIST(lista);
  dif2:=ΔLIST(dif1);
  dif3:=ΔLIST(dif2);
  dif4:=ΔLIST(dif3);
  a:=dif1(1);
  b:=dif1(2);
  c:=dif1(3);
  d:=dif1(4);
  e:=dif2(1);
  f:=dif2(2);
  g:=dif2(3);
  h:=dif3(1);
  j:=dif3(2);
  k:=dif4(1);
  jj:= 0; nm:=0; temp:= 0;
  WHILE 1 DO
  temp:= nm;
  nm:= (6*b+6*c-h-j+3*jj^2*(h+j)+2*jj^3*k)/(k-12*f);
    IF nm==temp THEN BREAK; END;
  jj:= nm;
  END;
  RETURN nm;
END;

EXPORT nzero3(lista)
BEGIN
  LOCAL a,b,c,R,dif1,dif2, n0, j;
  LOCAL temp;
  dif1:=ΔLIST(lista);
  dif2:=ΔLIST(dif1);
  a:=dif1(1);
  b:=dif1(2);
  c:=dif2(1);
  j:= 0; n0:=0; temp:= 0;
  WHILE 1 DO
  temp:= n0;
  n0:= -2*lista(2)/(a+b+c*j);
    IF n0==temp THEN BREAK; END;
  j:= n0;
  END;
  RETURN n0;
END;

EXPORT nzero5(lista)
BEGIN
  LOCAL a,b,c,R,dif1,dif2, dif3, dif4;
  LOCAL d,e,f,g,h,j,k, n0;
  LOCAL jj, temp;
  dif1:=ΔLIST(lista);
  dif2:=ΔLIST(dif1);
  dif3:=ΔLIST(dif2);
  dif4:=ΔLIST(dif3);
  a:=dif1(1);
  b:=dif1(2);
  c:=dif1(3);
  d:=dif1(4);
  e:=dif2(1);
  f:=dif2(2);
  g:=dif2(3);
  h:=dif3(1);
  j:=dif3(2);
  k:=dif4(1);
  jj:= 0; n0:=0; temp:= 0;
  WHILE 1 DO
  temp:= n0;
  n0:= (-24*lista(3)+jj^2*(k-12*f)-2*jj^3*(h+j)-jj^4*k)/(2*(6*b+6*c-h-j));
    IF n0==temp THEN BREAK; END;
  jj:= n0;
  END;
  RETURN n0;
END;

EXPORT halves4(lista)
BEGIN
  LOCAL y;
  y:= (9*(lista(2)+lista(3))-lista(1)-lista(4))/16;
  RETURN y;
END;

EXPORT lagrangeValue(lista_x, lista_y, xValue)
BEGIN
  LOCAL y, temp;
  temp:= lagrange(lista_x, lista_y);
  y:= EVAL(EXPR(temp+"|X="+xValue));
  RETURN y;
END;