/*
Copyright (C) 2008-2009 Mario Rodriguez Riotorto
This library is released under the terms of the GNU General Public License.

Simulacion de traxectorias parabolicas. 

v0: modulo da velocidade (m/s) de lanzamento
y0: altura de disparo (m)
a: angulo de disparo (entre 0 e 90 graos)
m: astro dende cuxa superficie se fai o lanzamento 
   (T, terra; L, Lua; S, Sol; M, Marte; X, Xupiter).

mario ARROBA edu PUNTO xunta PUNTO es
*/



trajectory(v0, y0, a, m):=
  block([vv0,aa,yy0],

    [vv0,aa,yy0]: float([v0,a,y0]),

    if not numberp(vv0) or vv0 < 0.0
      then "O m&oacute;dulo da velocidade debe ser positivo"
    elseif not numberp(aa) or aa >= 90.0 or aa < 0.0
      then "O &aacute;ngulo debe estar no rango [0,90)"
    elseif not numberp(yy0) or yy0 < 0.0
      then "A altura debe ser positiva"
    elseif not member(m, [T, L, S, M, X])
      then "Astro desco&ntilde;ecido"
    else  /* Imos coas contas */
      block([ar: aa*%pi/180, x, y, g,
             vx, vy,
             H, L, tL, eq, anima:[], aux,
             Hmax, Lmax],

        /* Selecciona gravidade */
        if m='T 
          then g: 9.81
        elseif m='L
          then g: 1.61
        elseif m='S
          then g: 273.7
        elseif m='M
          then g: 3.86
          else g: 26.00,

        /* Calcula os parametros do movemento */
        x : vv0*t*cos(ar),
        y : yy0 + vv0*t*sin(ar)-0.5*g*t^2,
        H : float(vv0^2*sin(ar)^2/(2*g)+yy0),
        tL: (vv0*sin(ar) + sqrt(2*g*yy0+vv0^2*sin(ar)^2))/g,
        L : float(vv0*cos(ar)*tL),
        eq: float('y = -g*'x^2/(2*vv0^2*cos(ar)^2) + tan(ar)*'x + yy0),
        Lmax: float(vv0*cos(%pi/4)*(vv0+sqrt(2*g*yy0+vv0^2))/g),
        Hmax: yy0+vv0^2/(2*g),

        /* Componentes do vector velocidade */
        vx: vv0*cos(ar),
        vy: vv0*sin(ar) - g*t,

        /* Grafico */
        for k:0 thru 15 do 
          (aux: [t=k*tL/15],
           anima: endcons(
                     gr2d(grid   = true,
                          xrange = [0,Lmax],
                          yrange = [0,Hmax],
                          dimensions = [400,350],
                          color  = blue,
                          parametric(x,y,t,0,tL),
                          point_size = 2,
                          point_type = filled_circle,
                          color      = red,
                          points([[subst(aux, x), subst(aux, y)]])),
                     anima)  ),
        apply(draw, append([terminal=animated_gif,
                              delay = 20],
                              anima)),

        /* Devolve resultados numericos */
        matrix(
          ['g = g * 'm/'s^2],
          [eq],
          ['x['max] = L * 'm],
          ['y['max] = H * 'm])) )$

