#!/usr/bin/perl
#
# Create an animation showing a standing wave,
#    as the superposition of two travelling waves
#
# MWR 4/7/2019

$mypi = 3.14159;

$datadir = "./stand_data";
$gif_delay = 30;

$num_snapshots = 40;
$speed = 1.0;
$wavelength = 5.0;
$k = (2.0*$mypi)/($wavelength);
$xmin = -10.0;
$xmax = 10.0;
$A = 1.0;
$period = 3.0;
$omega = 2.0*$mypi/($period);
$dx = 0.03;
$t_start = 0;
$t_end = 6.28;
$dt = ($t_end - $t_start) / $num_snapshots;

# we draw the standing wave only in this region
$stand_wave_xmin = -3.75;
$stand_wave_xmax = 3.75;
$stand_wave_other_y = 1000;

# step one: create one datafile per snapshot
for ($i = 0; $i < $num_snapshots; $i++) {
  $datafile[$i] = sprintf "%s/travel_%02d.dat", $datadir, $i;
  open(DATAFILE, ">$datafile[$i]") || die("can't open file $datafile[$i]");

  $t = $t_start + ($i*$dt);

  # constant amplitude
  $amp = $A;

  for ($x = $xmin; $x <= $xmax; $x += $dx) {
    $y = $amp*cos($k*$x - $speed*$t);

    # this is the standing wave version
    if (($x >= $stand_wave_xmin) && ($x <= $stand_wave_xmax)) {
      $y_stand = 2.0*$amp*cos(0.5*$omega*$t)*cos($k*$x);
    }
    else {
      $y_stand = $stand_wave_other_y;
    }
    printf DATAFILE "%lf %lf %lf \n", $x, $y, $y_stand;
  }
  close(DATAFILE);

}


  


# step three: plot each datafile to create an animated GIF image 
#                  which shows
#                      a) wave travelling in positive x dir
#                      b) wave travelling in negative x dir
#                      c) sum of the two waves

  my($output_file, $term_options);

  $sum_gif_delay = 2*$gif_delay;

  $output_file = "travel_c.gif";
  $term_options = "gif animate delay $sum_gif_delay  ";

  $cmdfile = "gnuplot.in";
  
  open (CMDFILE, ">$cmdfile") || die("can't open $cmdfile for writing");
  printf CMDFILE "set output '$output_file' \n";
  printf CMDFILE "set term $term_options \n";

  # revert to some oldish default colors
  printf CMDFILE "set style line 1 lt rgb 'red' lw 3  \n";
  printf CMDFILE "set style line 2 lt rgb 'green' lw 3  \n";
  printf CMDFILE "set style line 3 lt rgb 'blue' lw 3  \n";

  printf CMDFILE "set style line 4 lt rgb 'red' lw 1   \n";
  printf CMDFILE "set style line 5 lt rgb 'blue' lw 1  \n";
  printf CMDFILE "set style line 6 lt rgb 'black' lw 4  \n";


  printf CMDFILE "set grid \n";
  printf CMDFILE "set key top right \n";
  printf CMDFILE "set xrange [$xmin:$xmax]  \n";
  printf CMDFILE "set yrange [0-($A*2.2):0+($A*2.2)]  \n";
  printf CMDFILE "set xlabel 'X position (m)  '  \n";
  printf CMDFILE "set ylabel 'Y position (m) ' \n";
  # printf CMDFILE "set title 'Standing wave ' \n"; 

  for ($i = 0; $i < $num_snapshots; $i++) {
    printf CMDFILE "    plot '$datafile[$i]' using 1:2 with lines ls 4  t '' ";
    printf CMDFILE " , ";
    printf CMDFILE "       '$datafile[$i]' using (0-\$1):2 with lines ls 5  t '' ";
    printf CMDFILE " , ";
    printf CMDFILE "       '$datafile[$i]' using 1:3 with points ls 6 ps 0.2 pt 6  t '' ";
    printf CMDFILE "$cmd \n";
  }

  printf CMDFILE "$cmd \n";


  printf CMDFILE "set output \n";
  printf CMDFILE "quit \n";
  close(CMDFILE) ;
  
  $retval = `gnuplot < $cmdfile`;



  


# step four: plot each datafile to create an animated GIF image 
#                  which shows
#                      c) sum of the two waves only, zoomed in a bit
 
  


  my($output_file, $term_options);

  $sum_gif_delay = 2*$gif_delay;

  # this file will hold a pruned version of the original datafiles
  #    that we use for plotting this figure only
  $xmin = -4.0;
  $xmax = 4.0;

  $output_file = "travel_d.gif";
  $term_options = "gif animate delay $sum_gif_delay  ";

  $cmdfile = "gnuplot.in";
  
  open (CMDFILE, ">$cmdfile") || die("can't open $cmdfile for writing");
  printf CMDFILE "set output '$output_file' \n";
  printf CMDFILE "set term $term_options \n";

  # revert to some oldish default colors
  printf CMDFILE "set style line 1 lt rgb 'red' lw 3  \n";
  printf CMDFILE "set style line 2 lt rgb 'green' lw 3  \n";
  printf CMDFILE "set style line 3 lt rgb 'blue' lw 3  \n";

  printf CMDFILE "set style line 4 lt rgb 'red' lw 1   \n";
  printf CMDFILE "set style line 5 lt rgb 'blue' lw 1  \n";
  printf CMDFILE "set style line 6 lt rgb 'black' lw 4  \n";


  printf CMDFILE "set grid \n";
  printf CMDFILE "set key top right \n";
  printf CMDFILE "set xrange [$xmin:$xmax]  \n";
  printf CMDFILE "set yrange [0-($A*2.2):0+($A*2.2)]  \n";
  printf CMDFILE "set xlabel 'X position (m)  '  \n";
  printf CMDFILE "set ylabel 'Y position (m) ' \n";
  # printf CMDFILE "set title 'Standing wave ' \n"; 


  for ($i = 0; $i < $num_snapshots; $i++) {

    # we need to prune the datafile so that it runs only from
    #   the fixed positions at
    #      $stand_wave_xmin to $stand_wave_xmax
    $pruned_file[$i] = sprintf "$datadir/pruned_%05d.dat", $i;
    open(DATAFILE, "$datafile[$i]") || die("can't open datafile $datafile[$i]");
    open(PRUNED_FILE, ">$pruned_file[$i]") || die("can't create $pruned_file[$i]");
    while (<DATAFILE>) {
      $line = " " . $_;
      @words = split(/\s+/, $line);
      $x = $words[1];
      if (($x >= ($stand_wave_xmin - 0.0)) && ($x <= ($stand_wave_xmax + 0.0))) {
        printf PRUNED_FILE "%s", $line;
      }
    }
    close(PRUNED_FILE);
    close(DATAFILE);



    printf CMDFILE "    plot '$pruned_file[$i]' using 1:3 with lines ls 6 lw 2 t '' ";
    printf CMDFILE "$cmd \n";
  }

  printf CMDFILE "$cmd \n";


  printf CMDFILE "set output \n";
  printf CMDFILE "quit \n";
  close(CMDFILE) ;
  
  $retval = `gnuplot < $cmdfile`;




exit 0;