#!/usr/bin/perl
use Getopt::Std;
$blq=0;
@phf  = ( 0, 2, 2, 2, 2, 1, 1, 1, 1, 0,  0,  0 );
@indx = ( 0, 1, 2, 5, 6, 3, 7, 4, 8, 9, 10, 11 ); 
$amf=1.e-6;
$ampl=7;
$phal=8; $phald=10;
$endl=8;
$amfdm=1.e-5;
$show="GRAVITY [m/s^2]";
$tilt=0; $tiltamf=1.e-6;
$sign=1.;
$phaa=0.;
$st="";

getopts('Raehnrs');
if ( $opt_h )
{ 
    print "\nUSAGE:\n";
    print "            olgt.pl [options] blqfile\n\n";
    print "PURPOSE:\n";
    print "          Converts ocean loading gravity parameters\n";
    print "          from Automatic ocean loading provider format\n";
    print "          (a.k.a. BLQ) to TSOFT format\n\n";
    print "          Only one file can be processed at one time,\n";
    print "          and the file must be headered, i.e. the complete\n";
    print "          email that http://froste.oso.chalmers.se/loading\n";
    print "          upon a user request must be input.\n\n";
    print "          This procedure copes with the following BLQ variants:\n";
    print "          One-component    (Gravity only)\n";
    print "          Three-components (Gravity and Tilt)\n";
    print "          um/s**2 or mGal units\n\n";
    print "OPTIONS:\n";
    print " -e  -n   - Extract the tilt information, -e for east, -n for north.\n";
    print "            Default is gravity. -e and -n are mutually exclusive.\n\n";
    print " -r       - Reverse sign of amplitudes\n";
    print " -R       - Reverse sign by adding 180 degrees\n";
    print "            You can specify both -r -R but to no avail\n\n";
    print " -s       - Add the \"s\" before the tide symbol (old tradition)\n\n";
    print "OUTPUT:   Header record 1 contains the component and its units,\n";
    print "          the tide symbols and the name of the ocean model, \n";
    print "          taken from the BLQ input header.\n\n";
    print "          Header record 2 is issued per site segment of the\n";
    print "           BLQ file and contains\n";
    print "          station name, position, and a CompCode\n\n";
    print "          The CompCode gives you information about the ocean\n";
    print "          model and whether the post-processor was used\n";
    print "          CompCode examples:\n";
    print "          HAMTIDE_PP = the Hamburg ocean model, post-processed\n";
    print "          SCHW_1S    = the Schwiderski model in one step\n\n";
    print "          Post processing is applied when stations have been\n";
    print "          located near a coast, else the parameters come from\n";
    print "          the 1-step solution\n\n";
    exit 
}
# PROCESS OPTIONS 
#
if ( $opt_e ){ $ampl=8; $phald=11; $show="EAST TILT [rad]" }
if ( $opt_n ){ $ampl=9; $phald=12; $show="NORTH TILT [rad]" }
if ( $opt_r ){ $sign=-1; }
if ( $opt_R ){ $phaa=180. }
if ( $opt_s ){ $st="s" }
#___________________________________________________________________________
#
# PROCESS FILE 
#
while (<>){
    chomp;
    if ( ! $blq ) # _____________________________________ PROCESS BLQ HEADER
    {
	if ( m|COLUMN ORDER:(.*)| )
	{
	    $tt=$1;
	    $tt=~s/\W(\S+)/$1 /g;
	    $tt=~s/MF/Mf/; 
	    $tt=~s/MM/Mm/;
	    $tt=~s/SSA/Ssa/;
	    @t=( split(/ +/,$tt) )
	}
#
	if ( m|\$\$   TILT E| )
	{ 
	    $tilt=1; 
	    $phal=$phald; $endl=11; 
	    if ( $opt_e || $opt_n ){ $amf=$tiltamf; $amfdm=$tiltamf }
	}
#
	if ( m|\$\$ Ocean tide model: (.*)| )
	{ $model=$1 }
#
	if ( m|\(mGal| )
	{ $amf=$amfdm; print "BLQ Amplitudes were in mGal, using factor $amf\n" }
#
	if ( m|0\.1 urad| )
	{ $tiltamf = 1.e-7 }
#
	if ( m|END HEADER| )
	{ 
	    if ( ( $opt_n || $opt_e ) && ! $tilt )
	    { die "No tilt information in this file\n" }

	    $blq=1; $nr=0; 
	    print "Loading Tide in $show: @t from $model\n"; 
	}
    }
#
    if ( $blq )  #_________________________________________ CONVERT THE DATA 

    {   
	$nr++;
#        print "$nr $_\n";
        if ( $nr > $endl ) # _____________________________ PRINT THE SEGMENT
	{ 
	    print "Station $s long.= $x lat.= $y height= $h  CompCode=$m\n";
            for ($j=1;$j<12;$j++)
	    {
		$i=$indx[$j];
		$pv = -$p[$i]-$phf[$i]*$x+$phaa;
		if ( $pv < -180 ){ $pv+= 360. }
		if ( $pv >  180 ){ $pv-= 360. }
		printf "%-7s: %16.4e%15.2f\n",$st.$t[$i], $a[$i]*$amf*$sign, $pv
	    }
	    $nr=1; 
	    @a = ( );  @p = ( )
	}
#
	if ( m|lon/lat:\s+(\S+)\s+(\S+)\s+(\S+)| )  
	{ $x=$1; $y=$2; $h=$3 }
#
        if ( $nr == 4 )
	{
	    m|^...(\S+)\s+(\S+)\s+|; $m=$1.$2;
	    $m =~ s|Complete(.*)|$1_1S|;
	    $m =~ s|ID:||
	}
#
	if ( $_ =~ /^ / )
	{    
	    if ( $nr == 3 )
	    { 
		m|^..(.*)|; $s=$1
	    }
	    elsif ( $nr == $ampl )
	    {
		@a = ( @a, split(/ +/) )
	    }
	    elsif ( $nr == $phal )
	    {
		@p = ( @p, split(/ +/) )
	    }
	}
    }
}
exit

