#!/bin/csh -f

source /star/etc/login  >& /dev/null
source /star/etc/cshrc  >& /dev/null
kappa >& /dev/null

if (`alias kappa` == "") then
  echo ""
  echo "Fatal: Initialization error of kappa environment"
  echo ""
  exit
endif

set kap = $KAPPA_DIR

if ($#argv < 2) set argv = "help"

if ($argv[1] =~ *"help"* || $argv[1] == "-h") then
  echo ""
  echo "Use: sp_stats file v1:v2 [[v3:v4] ...]"
  echo "     e.g. 'sp_stats myfile -300:-100 100:300 ...'"
  echo ""
  echo "Note: the boundaries must be in increasing order!"
  echo "There is limited support for Starlink sections: a command"
  echo "like sp_stats myfile'(~1,~1)' -200 -100 100 200 will result"
  echo "in the combined stats of:"
  echo "  myfile'(~1,~1,-200.:-100.)' and myfile'(~1,~1,100.:200.)'"
  echo "i.e. the boundaries get appended to any section specified in"
  echo "the filename."
  echo ""
  exit
endif

# Get infile; avoid '~' character which is causing a lot of problems and
# strip suffix
set in  = $argv[1]
set in  = `echo "$in" | sed s/\~/\@/g` 
set in  = `echo "$in" | cut -d\. -f1`

# Handle section. Re-instate '~', if any, at end
if (`echo "$in" | fgrep -c '('` > 0) then
  set section = `echo "$in" | cut -d\( -f2`
  set section = `echo "$section" | cut -d\) -f1`
else
  set section = ","
endif
set section = "`echo $section | sed s/\@/\~/g`"

# Gte filename
set infile  = `echo "$in" | cut -d\( -f1`
set infile  = "${infile}.sdf"
shift argv
if !(-e "$infile") then
  echo ""
  echo "Fatal: file ${infile} not found"
  echo ""
  exit
endif

echo ""
echo "Stats for file ${infile}(${section}):"


# ----------------------------------------------------------------------
# Statistics
#

STATS:
    # Stats requested, can handle multiple areas.
    @ ii = 0
    set lims = `echo "$argv" | sed s/":"/" "/g`
    set clims = ""
    while ($#lims > 0)
      @ ii ++

      if (`echo $lims[1] | fgrep -c \.` == 0) set lims[1] = "${lims[1]}."
      if (`echo $lims[2] | fgrep -c \.` == 0) set lims[2] = "${lims[2]}."
      set clims = "${lims[1]}:${lims[2]}";
      
      echo "${kap}/stats ${infile}'(${section},${clims})'"
      ${kap}/stats ${infile}"(${section},${clims})"

      set pmean = "`${kap}/parget mean    stats | head -1`"
      set psig  = "`${kap}/parget sigma   stats | head -1`"
      set pmax  = "`${kap}/parget maximum stats | head -1`"
      set pmin  = "`${kap}/parget minimum stats | head -1`"
      set psum  = "`${kap}/parget total   stats | head -1`"
      set pnum  = "`${kap}/parget numgood stats | head -1`"

      if ($ii == 1) then
        set range = `echo "scale=6;((${lims[2]})-(${lims[1]}))" | bc`
        set mean = $pmean
        set sig  = $psig
        set max  = $pmax
        set min  = $pmin
        set sum  = $psum
        set num  = $pnum
      else 
        set range = `echo "scale=6;(($range)+(${lims[2]})-(${lims[1]}))" | bc`
        set mean = `echo "scale=6;((${mean}*${num})+(${pmean}*${pnum}))/(${num}+${pnum})" | bc`
        set sig  = `echo "scale=6;sqrt((${sig}*${sig}*${num}+${psig}*${psig}*${pnum})/(${num}+${pnum}))" | bc`
        if (`echo "1000000*(($pmax)-($max))"| bc | cut -d\. -f1` > 0) set max = "$pmax" 
        if (`echo "1000000*(($pmin)-($min))"| bc | cut -d\. -f1` < 0) set min = "$pmin" 
        set sum = `echo "scale=0;(${sum})+(${psum})" | bc`
        set num = `echo "scale=0;(${num})+(${pnum})" | bc`
      endif

      shift lims
      shift lims

    end

    set intint = `echo "scale=6;((${mean})*(${range}))" | bc`

    echo "      ========== COMBINED RESULTS =========="
    echo "      range = ${argv} (npoint = ${num}):"
    printf "      max = %8.3g,   min = %8.3g\n" ${max} ${min} 
    printf "      sum = %8.3g,  mean = %8.3g\n" ${sum} ${mean} 
    printf "      rms = %8.3g,  intI = %8.3g\n" ${sig} ${intint} 
    echo "      ======================================"
    echo ""