#!/bin/csh -f
#
if( $#argv == 0 ) then
    echo "sh_coseismic Usage:"
    echo "sh_coseismic -eqf <EQ file> -def <EQ Def name> -type <rapid/final/suppl/best>"
    echo "The EQ file should be the same one used in ../tables/globk_rep.cmd"
    echo "Other options:"
    echo "-cen <CEN>  Processing to used NMT/CWU/COM (default)"
    echo "-ndays <n-days> Number of days before and after (default 2)"
    echo "-rerun      Option will force re-run of globk solution.  Normally, once the" 
    echo "            globk run is completed, only the tables and plot are regenerated"
    echo "            with each run"
    echo "-vs         Vector scale (default 1)"
    echo "-arrow_val  Length of scale vector (default 20 mm)."
    echo " "
    exit
endif
# Suggested File name format
# pbo_120905_1442_EQ13_coseis_final.evt.timestamp
# 
#
# Defaults
set type = 'best'
set cen = 'COM'
set nam = 'PBO'
set ndays = 2
set rerun = 'N'
set vs = '1'
set arrow_value = 20
#

foreach a ($argv)
  switch($a)
    case -e*:
       set eqf = $argv[2]
       breaksw
    case -d*:
       set def = $argv[2]
       breaksw
    case -t*:
       set type = $argv[2]
       breaksw
    case -c*:
       set cen = $argv[2]
       set nam = $cen
       if( $nam == 'COM' ) set nam = "PBO"
       breaksw
    case -n*:
       set ndays = $argv[2]
       breaksw
    case -r*
       set rerun = 'Y'
       breaksw
    case -v*
       set vs = $argv[2]
       breaksw
    case -a*
       set arrow_value = $argv[2]
       breaksw
    default:
#       echo Unknown option $a
#       grep '^## ' $0 | awk '{print substr($0,3)}'
#       exit

  endsw
  shift
end

#
echo "Settting Run for EQ $def from EQ file $eqf"
echo "Type of analysis: $type"
#
# Get letter for type
set tl = 'x'
if( $type == 'rapid' ) then
   set tl = 'a'
else if ( $type == 'final' ) then
   set tl = 'b'
else if ( $type == 'suppl' ) then
   set tl = 'c'
else if ( $type == 'repro' ) then
   set tl = 'e'
else if ( $type == 'best' ) then
   set tl = 'x'
else
   echo "Analysis type $type not found (valid rapid,final,suppl,repro)"
   exit
endif

#
# See of we can find earthquake (requires standard format).
set eqdef = `grep '^ ' $eqf | grep -i " eq_def " | awk -v def=$def '{if ($2 == def ) {print $0}}'` 
if( $#eqdef > 0 ) then   # EQ-found so continue
   # Get the date
   echo $eqdef 
   set eqdoy = `echo $eqdef | awk '{print $7,$8,$9}'`

   # If best type specified; look at available data use best (e/c/b/a in that order).
   if( $type == 'best' ) then
       set dayp = `doy $eqdoy | head -1 | awk -v n=$ndays '{printf("%12.1f",$10+n)}'`
       set tf   = `doy $dayp | head -2 | tail -1 | awk -v cen=$cen '{printf("../%s/%s%4.4d%1d.?.glb\n",cen,tolower(cen),$3,$7)}'`
       set tl   = `ls $tf | tail -1 | awk '{print substr($0,17,1)}'`
       if( $tl == 'a' ) then
           set type = rapid
       else if ( $tl == b ) then
           set type = final
       else if ( $tl == c ) then
           set type = suppl
       else if ( $tl == e ) then
           set type = repro
       else
           echo $tf found unknown type $tl
           exit -2
       endif
       echo "Best Analysis is $type code $tl"
   endif

   # Now set up GDL file.
   set gdlf = `doy $eqdoy | head -1 | awk -v cen=$cen -v nam=$nam -v tl=$tl -v def=$def '{printf("../%s/%s%2.2d%2.2d%2.2d_EQ%s.%s.gdl",cen,tolower(nam),substr($2,3,2),substr($2,6,2),substr($2,9,2),def,tl)}'`
   echo "Creating $gdlf for coseismic offset estimates"
   echo "# GDL File EQ $def from $eqf Analyis type $type Analysic center $cen" >! $gdlf

   # Now generate the names of the files needed
   set n = 0
   while ( $n < $ndays )
        @ n++
       set daym = `doy $eqdoy | head -1 | awk -v n=$n '{printf("%12.1f",$10-n)}'`
       doy $daym | head -2 | tail -1 | awk -v cen=$cen -v tl=$tl '{printf("../%s/%s%4.4d%1d.%s.glb\n",cen,tolower(cen),$3,$7,tl)}' >>  $gdlf
       set dayp = `doy $eqdoy | head -1 | awk -v n=$n '{printf("%12.1f",$10+n)}'`
       doy $dayp | head -2 | tail -1 | awk -v cen=$cen -v tl=$tl '{printf("../%s/%s%4.4d%1d.%s.glb\n",cen,tolower(cen),$3,$7,tl)}' >>  $gdlf
   end 
#  Now run the analysis of the data
   pushd ../$cen
      if( ! -e $gdlf:r.org || $rerun == 'Y' ) then
          globk 6 $gdlf:r.prt $gdlf:r.log $gdlf ../tables/globk_rep.cmd  COSEIS
      else
          echo "$gdlf:r.org exists and using this file"
      endif
      sh_exglk -f $gdlf:r.org -eq $def $gdlf:r.coff

#     Now generate the root for file name of the event
# Suggested File name format
# pbo_120905_1442_EQ13_coseis_final.evt.timestamp
echo $nam $eqdef $def $type
      set cofroot = `echo $nam $eqdef $def $type | awk '{printf("%s_%2.2d%2.2d%2.2d_%2.2d%2.2d_EQ%s_coseis_%s",tolower($1),substr($8,3,2),$9,$10,$11,$12,$14,$15)}'`
echo $cofroot

#Release Date  : 
# Run time               : 2012/10/ 9 13:13 47.00
      set reldate = `grep 'Run time               :' $gdlf:r.org | awk '{gsub("/"," ")} {gsub(":"," ")} {printf("%4d%2.2d%2.2d%2.2d%2.2d%2.2d",$3,$4,$5,$6,$7,$8)}'`
# 
#     Now generate the plot of the results
      set R = `echo $eqdef | awk '{dlat=2*$5/111} {dlng=2*$5/111/cos($3/57.3)} {printf("-R%.4f/%.4f/%.4f/%.4f",$4-dlng,$4+dlng,$3-dlat,$3+dlat)}'`
      set dlat = `echo $eqdef | awk '{print int((2*$5/111)*10)}'`
      if( $dlat < 5 ) then
          set B = '-Ba0.1f0.05/a0.1f0.05WeSn'
          set L = '-linearscale 5'
      else if ( $dlat < 20 ) then
          set B = '-Ba0.5f0.25/a0.5f0.25WeSn'
          set L = '-linearscale 50'
      else
          set B = '-Ba2f0.5/a2f0.5WeSn'
          set L = '-linearscale 200'
      endif
      echo $dlat $B $L
#Release Date  : 
# Run time               : 2012/10/ 9 13:13 47.00
      set reldate = `grep 'Run time               :' $gdlf:r.org | awk '{gsub("/"," ")} {gsub(":"," ")} {printf("%4d%2.2d%2.2d%2.2d%2.2d%2.2d",$3,$4,$5,$6,$7,$8)}'`

      set epc = `echo  $eqdef | awk '{if( $4 >0 ) {printf("-epicentre %.4f %.4f",$3,$4)} else {printf("-epicentre %.4f %.4f",$3,$4+360)}}'`
      set epr = `echo  $eqdef | awk '{printf("-epiradius %.2f",$5)}'`
      sh_plotvel -f $gdlf:r.coff -ps $cofroot $R $epc $epr -vs $vs -arrow_value $arrow_value  $B $L
      sed -i 's/mm\/yr/mm/g' ${cofroot}.ps 
      \mv  ${cofroot}.ps  ${cofroot}.ps.${reldate} 
      echo "Created ${cofroot}.ps.${reldate}"
 
#     Now create the co-seismic offset file
# Get the data analyzed
      set analfs = `tail -n +2 $gdlf | sort -n`
      set cofffl = ${cofroot}.evt.${reldate}
      echo "Creating $cofffl"
      echo $def $eqdef | awk '{printf("PBO Coseismic Offsets for EQ code %s, Date %4d %2.2d %2.2d (ymd) Time %2.2d %2.2d (hr min)\n",$1,$8,$9,$10,$11,$12)}' >! $cofffl
      echo 1.0         | awk '{printf("Format Version         : %4.2f\n",$1)}'  >> $cofffl
      echo $eqdef      | awk '{printf("EQ Location (lat/long) : %8.4f %8.4f deg\n",$3,$4)}' >> $cofffl
      echo $reldate    | awk '{printf("Release Date           : %s\n",$1)}' >> $cofffl
      echo $analfs     | awk '{printf("Data Analyzed          : %s\n",$0)}' >> $cofffl
      echo $nam        | awk '{printf("Analysis Center        : %s\n",$1)}' >> $cofffl
      echo $type       | awk '{printf("Analysis Type          : %s\n",$1)}' >> $cofffl
      echo "Event Type             : Coseismic" >> $cofffl
      echo $def        | awk '{printf("Event Code             : EQ%s\n",$1)}' >>  $cofffl
      echo $def $eqdef | awk '{printf("Event Date             : %4.4d%2.2d%2.2d%2.2d%2.2d\n",$8,$9,$10,$11,$12)}' >>  $cofffl
      echo $gdlf:r.org | awk '{printf("GLOBK Solution file    : %s\n",$1)}' >> $cofffl
      echo ".  Long        Lat        dE       dN       E +-     N +-   RhoEN      dH      H +-  Site"  >> $cofffl
      echo ".   deg        deg        mm       mm        mm       mm     --        mm       mm     --" >>  $cofffl
      grep '^ '  $gdlf:r.coff | awk '{printf(" %9.5f %9.5f %8.2f %8.2f  %8.2f %8.2f   0.000 %8.2f %8.2f  %s\n",$1,$2,$3,$4,$7,$8,$10,$12,$13)}'>> $cofffl

endif
echo 'To Send to Unavco use:'
echo 'ssh -n ldm@everest.mit.edu "cd /net/lehman/raid8/tah/ACC_PBO/anal/'${cen}'; ~/bin/pqinsert *.'${reldate}'"'
