#read the satellite ephemerids by AB
 #containing hh mm ss ra dec (decimal deg) and ...
 #compute velocity in ra and dec in arcsec/sec
#warning: if data begins after RA=359 deg or ends  after RA=359deg
#you may have problems


go

 data "../uploads/sat26867.txt"
 #data "uploads/sat26867.txt"
 lines 1 1
 read {yy 1 mm 2 dd 3 name 5.s}
 
 lines 2 0
 read {hh 1 mm 2 ss 3 radeg 4 decdeg 5}
 erase
 ctype black
 expand 1.2
 lweight 1.2
 ltype 0
 ptype 10 0

 
 set t=hh+mm/60+ss/3600

 vecminmax radeg radegmin radegmax
 vecminmax radeg radegmin1 radegmax1
 
 define n dimen(radeg)

     if($radegmax>359&&$radegmin<1){write standard "ra through 0"     
                          set xm= radeg-360 if(radeg<360&&radeg>180)
                          set ym= decdeg if(radeg<360&&radeg>180)
                          set tm=t if(radeg<360&&radeg>180)
                          set xp= radeg if(radeg>=0&&radeg<180) 
                          set yp= decdeg if(radeg>=0&&radeg<180) 
                          set tp=t if(radeg>=0&&radeg<180) 
                          }
     if($radegmax<359&&$radegmin>=0){write standard "ra through 0"     
                          set xm= radeg 
                          set ym= decdeg
                          set tm=t
                          set xp= radeg
                          set yp= decdeg 
                          set tp=t 
                          }


 
 
 if($radegmax>359&&$radegmin<1){
     set xx=xm concat xp
     set yy=ym concat yp
     set tt=tm concat tp
                   }

 if($radegmax<359&&$radegmin>=0){
     set xx=xm 
     set yy=ym 
     set tt=tm 
                 }

 define n dimen(xx)
 define xmin (xx[$0])
 define xmax (xx[$n-1])
 
 LOCATION 4500 31000 23000 31000
 limits $xmin $xmax decdeg
 box
 xlabel RA [deg]
 ylabel DEC [deg]
 points  radeg decdeg
 ptype 10 3
 ctype red
 points xm ym
 ctype blue
 points xp yp
 ptype 10 0
 ctype black
 points xx yy


 vecminmax xp xpmin xpmax
 vecminmax xm xmmin xmmax
 vecminmax tp tpmin tpmax
 vecminmax tm tmmin tmmax

 set xpspline =$xpmin,$xpmax,0.1
 set xmspline =$xmmin,$xmmax,0.1
 define dt (1./3600) #1 second
 set tpspline =$tpmin,$tpmax,$dt
 set tmspline =$tmmin,$tmmax,$dt

 set xspline=xmspline concat xpspline
 set tspline=tmspline concat tpspline
  
 
 spline xx yy xspline yspline
 spline t xx tspline xxspline 
 spline t yy tspline yyspline

 ctype green
 expand 0.7
 points xspline yspline
 
 ctype cyan
 points xxspline yyspline 
 expand 1.2

 define n dimen(tspline)
 set dx=xxspline*0
 set dy=dx
 set deltat=dx
 do j=0,$n-1,1 {
                set dx[$j]=xxspline[$j+1]-xxspline[$j]
                set dy[$j]=yyspline[$j+1]-yyspline[$j]
	        set deltat[$j]=(tspline[$j+1]-tspline[$j])*3600 #(delta t in sec)

                                 }
  set vx=dx/deltat
  set vy=dy/deltat 
  print satvel.dat '%3.5f %3.5f  %3.5f \n' { tspline  vx vy}
  write standard " " 
  write standard ">>>>>>>>>> Data in satvel.dat <<<<<<<<<<<<"
  write standard " " 
  
  LOCATION 4500 31000 13000 20000
  ctype black
  expand 1.2
  limits tspline vx
  box
  xlabel t [hh.hh]
  ylabel \mu RA [''/sec]
  points tspline vx 
  

  LOCATION 4500 31000 3000 10000
  limits tspline vy

  box
  xlabel t [hh.hh]
  ylabel \mu DEC [''/sec]
  points tspline vy 
  
   

  
