#!/bin/tcsh -f 

# # ##################################################
# # Free fall: projectle
# ./rk2.csh -n -f 2 -- 10 0 0 0 1 1 20000 0 0.2
# gnuplot <<EOF
# set term post enhanced color eps 24
# set out "rk2_g_x.eps"
# plot "rk2.dat" using 1:2 with lines title "x(t)"
# set out "rk2_g_y.eps"
# plot "rk2.dat" using 1:3 with lines title "y(t)"
# set out "rk2_g_vx.eps"
# plot "rk2.dat" using 1:4 with lines title "v_x(t)"
# set out "rk2_g_vy.eps"
# plot "rk2.dat" using 1:5 with lines title "v_y(t)"
# set out "rk2_g_E.eps"
# plot "rk2.dat" using 1:(\$6-1.0) with lines title "E(t)-E(0)"
# set size square
# set xlabel "x"
# set ylabel "y"
# set out "rk2_g_xy.eps"
# plot "rk2.dat" using 2:3 with lines notit
# EOF
# # ##################################################
# # Air resistance ~ v + fall
# rk2.csh -n -f 3 -- 10.0 30  0 0 1 1 20000 0 0.2; mv rk2.dat rk2.k30.dat
# rk2.csh -n -f 3 -- 10.0 0   0 0 1 1 20000 0 0.2; mv rk2.dat rk2.k0.dat
# rk2.csh -n -f 3 -- 10.0 0.2 0 0 1 1 20000 0 0.2; mv rk2.dat rk2.k0.2.dat
# rk2.csh -n -f 3 -- 10.0 1.0 0 0 1 1 20000 0 0.2; mv rk2.dat rk2.k1.dat
# rk2.csh -n -f 3 -- 10.0 5.0 0 0 1 1 20000 0 0.2; mv rk2.dat rk2.k5.dat
# rk2.csh -n -f 3 -- 10.0 10  0 0 1 1 20000 0 0.2; mv rk2.dat rk2.k10.dat
# rk2.csh -n -f 3 -- 10.0 20  0 0 1 1 20000 0 0.2; mv rk2.dat rk2.k20.dat
# gnuplot <<EOF
# set term post enhanced color eps 24
# set out "rk2_vg.eps"
# set xlabel "x"
# set ylabel "y"
# plot [0:0.23][0.0:0.055] \
# "rk2.k0.dat"   u 2:3 w l t " 0.0", \
# "rk2.k0.2.dat" u 2:3 w l t " 0.2", \
# "rk2.k1.dat"   u 2:3 w l t " 1.0", \
# "rk2.k5.dat"   u 2:3 w l t " 5.0", \
# "rk2.k10.dat"  u 2:3 w l t "10.0", \
# "rk2.k20.dat"  u 2:3 w l t "20.0", \
# "rk2.k30.dat"  u 2:3 w l t "30.0"
# set out
# EOF
# # Air resistance ~ v^2 + fall
# rk2.csh -n -f 4 -- 10.0 30  0 0 1 1 20000 0 0.2; mv rk2.dat rk2.k30.dat
# rk2.csh -n -f 4 -- 10.0 0   0 0 1 1 20000 0 0.2; mv rk2.dat rk2.k0.dat
# rk2.csh -n -f 4 -- 10.0 0.2 0 0 1 1 20000 0 0.2; mv rk2.dat rk2.k0.2.dat
# rk2.csh -n -f 4 -- 10.0 1.0 0 0 1 1 20000 0 0.2; mv rk2.dat rk2.k1.dat
# rk2.csh -n -f 4 -- 10.0 5.0 0 0 1 1 20000 0 0.2; mv rk2.dat rk2.k5.dat
# rk2.csh -n -f 4 -- 10.0 10  0 0 1 1 20000 0 0.2; mv rk2.dat rk2.k10.dat
# rk2.csh -n -f 4 -- 10.0 20  0 0 1 1 20000 0 0.2; mv rk2.dat rk2.k20.dat
# gnuplot <<EOF
# set term post enhanced color eps 24
# set out "rk2_v2g.eps"
# set xlabel "x"
# set ylabel "y"
# plot [0:0.23][0.0:0.055] \
# "rk2.k0.dat"   u 2:3 w l t " 0.0", \
# "rk2.k0.2.dat" u 2:3 w l t " 0.2", \
# "rk2.k1.dat"   u 2:3 w l t " 1.0", \
# "rk2.k5.dat"   u 2:3 w l t " 5.0", \
# "rk2.k10.dat"  u 2:3 w l t "10.0", \
# "rk2.k20.dat"  u 2:3 w l t "20.0", \
# "rk2.k30.dat"  u 2:3 w l t "30.0"
# set out
# EOF
# \rm *.dat
# rk2.csh -n -f 5 -- -10 0 1 0 0 3.6  40000 0 10; mv rk2.dat rk2.k3.6.dat
# rk2.csh -n -f 5 -- -10 0 1 0 0 3.8  40000 0 10; mv rk2.dat rk2.k3.8.dat
# rk2.csh -n -f 5 -- -10 0 1 0 0 4    20000 0  8; mv rk2.dat rk2.k4.0.dat
# rk2.csh -n -f 5 -- -10 0 1 0 0 4.1  50000 0 20; mv rk2.dat rk2.k4.1.dat
# rk2.csh -n -f 5 -- -10 0 1 0 0 4.3 100000 0 60; mv rk2.dat rk2.k4.3.dat
# gnuplot <<EOF
# set size square
# set term post enhanced color eps 24
# set out "rk2_planet1.eps"
# set xlabel "x"
# set ylabel "y"
# plot [-14:2][-8:8] \
# "rk2.k3.6.dat"   u 2:3 w l t " 1.68",\
# "rk2.k3.8.dat"   u 2:3 w l t " 2.40",\
# "rk2.k4.0.dat"   u 2:3 w l t " 3.93",\
# "rk2.k4.1.dat"   u 2:3 w l t " 5.52",\
# "rk2.k4.3.dat"   u 2:3 w l t "16.95"
# set out
# EOF
# gnuplot <<EOF
# set log
# set xlabel "a^3"
# set ylabel "T^2"
# set term post enhanced color eps 24
# set out "rk2_planetT.eps"
# plot [:12000] \
#  "planet_T" u ((\$3/2)**3):((2*\$2)**2) notit , ((4*pi**2)/10)*x notit
# set out
# EOF
# # #############################################
# # Other force laws:
# \rm *.dat
# rk2.csh -n -f 6 -- -10 3 1 0 0 3. 80000 0 1.  # 1/r^3 - spiral
# gnuplot <<EOF
# set xlabel "x"
# set ylabel "y"
# set size square
# set term post enhanced color eps 24
# set out "rk2_r3_spiral.eps"
# plot "rk2.dat"  u 2:3 w l notit
# set out
# EOF
# # #############################################################
# # #############################################################
# # #############################################################
#
# Plot of scattering trajectories in Coulomb potential:
#
# \rm rk2b*.dat
# ./rk2.csh -f 5 -n -- 1 0 -50 0.008 3 0 1000 0 30;cat rk2.dat >> rk2bb.dat;mv rk2.dat rk2b0.008.dat
# ./rk2.csh -f 5 -n -- 1 0 -50 0.015 3 0 1000 0 30;cat rk2.dat >> rk2bb.dat;mv rk2.dat rk2b0.015.dat
# ./rk2.csh -f 5 -n -- 1 0 -50 0.020 3 0 1000 0 30;cat rk2.dat >> rk2bb.dat;mv rk2.dat rk2b0.020.dat
# ./rk2.csh -f 5 -n -- 1 0 -50 0.035 3 0 1000 0 30;cat rk2.dat >> rk2bb.dat;mv rk2.dat rk2b0.035.dat
# ./rk2.csh -f 5 -n -- 1 0 -50 0.080 3 0 1000 0 30;cat rk2.dat >> rk2bb.dat;mv rk2.dat rk2b0.080.dat
# ./rk2.csh -f 5 -n -- 1 0 -50 0.120 3 0 1000 0 30;cat rk2.dat >> rk2bb.dat;mv rk2.dat rk2b0.120.dat
# ./rk2.csh -f 5 -n -- 1 0 -50 0.200 3 0 1000 0 30;cat rk2.dat >> rk2bb.dat;mv rk2.dat rk2b0.200.dat
# ./rk2.csh -f 5 -n -- 1 0 -50 0.240 3 0 1000 0 30;cat rk2.dat >> rk2bb.dat;mv rk2.dat rk2b0.240.dat
# ./rk2.csh -f 5 -n -- 1 0 -50 0.320 3 0 1000 0 30;cat rk2.dat >> rk2bb.dat;mv rk2.dat rk2b0.320.dat
# ./rk2.csh -f 5 -n -- 1 0 -50 0.450 3 0 1000 0 30;cat rk2.dat >> rk2bb.dat;mv rk2.dat rk2b0.450.dat
# ./rk2.csh -f 5 -n -- 1 0 -50 0.600 3 0 1000 0 30;cat rk2.dat >> rk2bb.dat;mv rk2.dat rk2b0.600.dat
# ./rk2.csh -f 5 -n -- 1 0 -50 1.500 3 0 1000 0 30;cat rk2.dat >> rk2bb.dat;mv rk2.dat rk2b1.500.dat
# gnuplot <<EOF
# set term post enhanced color eps 24
# set out "scatter_traj.eps"
# plot [-20:20]"rk2bb.dat" u 2:3 w dots notit
# set out 
# EOF
# # #############################################################
# # #############################################################
# #                Scattering Cross section
# # #############################################################
#./scatter.csh -f 1 -- 1 0 -50 0.02 3 0 5000 0 30 1 0.0002
# grep ds= scatter.dat | awk '{print log(sin($2/2)),log($4)}' | awk '$1>-1.5' | fit-straight-line -n
# gnuplot <<EOF
# set term post enhanced color eps 24
# set out "scatter_sigma.02.eps"
# set log
# set xlabel "sin^{-4}({/Symbol q}/2)"
# set ylabel "{/Symbol s({/Symbol q})}"
# plot [:1000] "<grep ds= scatter.dat" u ((sin(\$2/2))**(-4)):4 notit,(1./(4.*3.**4))*x notit
# set out
# unset log
# set out "scatter_sigma.01.eps"
# set log y
# set xlabel "{/Symbol q}"
# plot [:] "<grep ds= scatter.dat" u 2:4 notit,(1./(4.*3.**4))*(sin(x/2))**(-4) notit
# set out
# unset log
# EOF
# #         scattering of positronium by a hydrogen atom:
# # Determine this way the sigma_tot = pi a^2
# foreach a (3 2 1 0.5)
#  ./scatter.csh -f 3 -- $a 0 -5 0.01 0.5 0 2000 0 40 0.99999 0.00002
# end
# # Make qualitative plots of b(theta),sigma(theta) as functions of v:
# # gives crude sigma_tot = 3.1422 instead of 3.1459
# foreach v (2.5 2.0 1.5 1.0 0.5 0.25 0.125)
#  ./scatter.csh -f 3 -- 1 0 -5 0.01 $v 0 4000 0 40 0.99999 0.0002
#  mv scatter.dat scatter_v$v
# end
# gnuplot <<EOF
#  set term post enhanced color eps 24
#  set out "scatter_hy_b_v.eps"
#  set xlabel "b"
#  set ylabel "{/Symbol q}"
#  plot [:] "<grep @ scatter_v2.0"   u 2:3 t "2.0"   w l,\
#           "<grep @ scatter_v1.5"   u 2:3 t "1.5"   w l,\
#           "<grep @ scatter_v1.0"   u 2:3 t "1.0"   w l,\
#           "<grep @ scatter_v0.5"   u 2:3 t "0.5"   w l,\
#           "<grep @ scatter_v0.25"  u 2:3 t "0.25"  w l,\
#           "<grep @ scatter_v0.125" u 2:3 t "0.125" w l
#  set out
#  set out "scatter_hy_s_v.eps"
#  set xlabel "{/Symbol q}"
#  set ylabel "{/Symbol s(q)}"
#  set log y
#  plot [:][0.008:]\
#           "<grep ds= scatter_v2.0"   u 2:4 t "2.0"   w l,\
#           "<grep ds= scatter_v1.5"   u 2:4 t "1.5"   w l,\
#           "<grep ds= scatter_v1.0"   u 2:4 t "1.0"   w l,\
#           "<grep ds= scatter_v0.5"   u 2:4 t "0.5"   w l,\
#           "<grep ds= scatter_v0.25"  u 2:4 t "0.25"  w l,\
#           "<grep ds= scatter_v0.125" u 2:4 t "0.125" w l
#  set out
# EOF
# ##############################################################
# #                       Yukawa scattering
# ##############################################################
# foreach v (10.0 15.0) # (5.0 4.0 3.0 2.0)
#   ./scatter.csh -f 4 -n -- 1 1 -50 0.0001 $v 0 5000 0 30 3 0.001
#   mv scatter.dat scatter_yu_v$v
# end
# gnuplot <<EOF
#  set term post enhanced color eps 24
#  set out "scatter_yu_b_v.eps"
#  set xlabel "b"
#  set ylabel "{/Symbol q}"
#  set log y
#   plot [0.001:3] \
#            "<grep @ scatter_yu_v4.0"    u 2:3 t "yu v= 4.0"   w l,\
#            "<grep @ scatter_yu_v4.0"    u 5:3 t "cb v= 4.0"    w l,\
#            "<grep @ scatter_yu_v15.0"   u 2:3 t "yu v=15.0"   w l,\
#            "<grep @ scatter_yu_v15.0"   u 5:3 t "cb v=15.0"   w l
#  set out
# EOF
# foreach a (0.25 0.5 1.0 1.5 2.0)
#   ./scatter.csh -f 4 -n -- 1 $a -50 0.0001 4.0 0 5000 0 30 4 0.001
#   mv scatter.dat scatter_yu_a$a
# end
gnuplot <<EOF
  set term post enhanced color eps 24
  set out "scatter_yu_b_a.eps"
  set xlabel "b"
  set ylabel "{/Symbol q}"
  set log y
   plot [0.001:3] \
            "<grep @ scatter_yu_a0.25"  u (\$2):3 t "0.25"   w l,\
            "<grep @ scatter_yu_a0.5"   u (\$2):3 t "0.50"   w l,\
            "<grep @ scatter_yu_a1.0"   u (\$2):3 t "1.00"   w l,\
            "<grep @ scatter_yu_a1.5"   u (\$2):3 t "1.50"   w l,\
            "<grep @ scatter_yu_a2.0"   u (\$2):3 t "2.00"   w l
  set out
EOF
#  ---------------------------------------------------------------------
#  Copyright by Konstantinos N. Anagnostopoulos (2004-2014)
#  Physics Dept., National Technical University,
#  konstant@mail.ntua.gr, www.physics.ntua.gr/~konstant
#  
#  This program is free software: you can redistribute it and/or modify
#  it under the terms of the GNU General Public License as published by
#  the Free Software Foundation, version 3 of the License.
#  
#  This program is distributed in the hope that it will be useful, but
#  WITHOUT ANY WARRANTY; without even the implied warranty of
#  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
#  General Public License for more details.
#  
#  You should have received a copy of the GNU General Public Liense along
#  with this program.  If not, see <http://www.gnu.org/licenses/>.
#  -----------------------------------------------------------------------
