Post Reply 
Astronomy: interpolation
07-14-2015, 08:56 PM (This post was last modified: 07-14-2015 10:35 PM by salvomic.)
Post: #1
Astronomy: interpolation
hi Marcel and every friends,
about Astronomy it's very important to have some tools to get interpolation (in part used in our Astro Lab and Effemeridi for Prime).
I think this is the base for an ancillary program for almost every task.
Please, help to improve and doing the necessary controls to avoid logic errors and crashes too...

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...
Formulae by Jean Meeus, Astronomical Algorithms, very precise, even more with the power of the Prime!

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;

∫aL√0mic (IT9CLU) :: HP Prime 50g 41CX 71b 42s 39s 35s 12C 15C - DM42, DM41X - WP34s Prime Soft. Lib
Visit this user's website Find all posts by this user
Quote this message in a reply
Post Reply 




User(s) browsing this thread: 1 Guest(s)