#!/usr/bin/perl -w
use strict;
use Readonly;
use Data::Dumper;

# a: calcul L1_ONLY; b: calcul L2_ONLY
Readonly my %EXTENSION_DE => (
                            'L1'   =>  'a',
                            'L2'   =>  'b',
                            );
Readonly my %VU_DE => (
                            'IG01'  =>  'IG01_IG(\d\d)',
                            'SMNE'  =>  'IG(\d\d)_SMNE',
                            );
Readonly my $FIC_SERIE_TEMP =>  'st.ps';

##############################################################
# récupération des noms d'antenne par à partir du station.info
##############################################################
my %antenne_en;
my $commentaire_precedent = '';
open( STATION_INFO, "tables/station.info" );
while ( defined( my $line = <STATION_INFO>)){
    chomp $line;
    if ($line =~ m{\A\ IG(\d\d)\ \  # nom du site
                    (.{16})\ \      # commentaire du site
                    (\d{4})\s+      # annee
                    (\d{1,3})\s+    # jour
                }xm ){
        # pour eviter d'avoir deux fois le nom de la meme antenne
        # ne fonctionne que si la ligne de l'antenne au sud suit
        # la ligne de l'antenne au nord
        if ( $commentaire_precedent ne $2 ){
            $antenne_en{$1}->{$4} = $2;
        }
        $commentaire_precedent = $2;
    }
}
close STATION_INFO;
######################################
# Traitement des résultats des calculs
######################################
my $jours_ref = [];
my %resultats_pour;
foreach my $pivot ( keys( %VU_DE ) ){
    print "*********************************\n";
    print "Lignes de bases vues depuis ${pivot}:\n";
    print "*********************************\n";
    foreach my $freq ( sort keys( %EXTENSION_DE ) ){
# recuperation des valeurs
        print "\t########\n";
        print "\t## $freq ##\n";
        print "\t########\n";
        my @jours;
        my %position_de;
        my @repertoires_jour = glob("2007_???$EXTENSION_DE{$freq}");
        foreach my $repertoire_jour ( @repertoires_jour ) {
            my $annee = substr($repertoire_jour,0,4);
            my $jour = substr($repertoire_jour,5,3);
            push(@jours,$jour);
            my @fichier_o = glob( "$repertoire_jour/ocalba*");
            open( FICHIER_O, "<$fichier_o[0]" );
            while ( defined( my $line = <FICHIER_O>)){
                chomp($line);
                if ($line =~ m{\A$VU_DE{$pivot}\     # ligne de base
                                $annee\.${jour}X\    # ambiguites fixees
                                N\ +([ +-]\d+\.\d+)  # coordonnee Nord
                                \ \+\-\ +(\d+\.\d+)\ # rms nord
                                E\ +([ +-]\d+\.\d+)  # coordonnee Est
                                \ \+\-\ +(\d+\.\d+)\ # rms nord
                                U\ +([ +-]\d+\.\d+)  # coordonnee Up
                                \ \+\-\ +(\d+\.\d+)\ # rms nord
                                L\ +(\d+\.\d+)       # longueur
                                \ \+\-\ +(\d+\.\d+)  # rms nord
                                }xm ){
                    $position_de{$1}{$jour}{'N'}    = $2;
                    $position_de{$1}{$jour}{'Nrms'} = $3;
                    $position_de{$1}{$jour}{'E'}    = $4;
                    $position_de{$1}{$jour}{'Erms'} = $5;
                    $position_de{$1}{$jour}{'U'}    = $6;
                    $position_de{$1}{$jour}{'Urms'} = $7;
                    $position_de{$1}{$jour}{'L'}    = $8;
                    $position_de{$1}{$jour}{'Lrms'} = $9;
                }
            }
            close( FICHIER_O );
        }
# calcul des moyennes
        foreach my $position (keys(%position_de)){
            LOOP1_ON_DAYS:
            foreach my $jour (@jours){
                next LOOP1_ON_DAYS
                    if ( !exists( $position_de{$position}->{$jour} ) );
                $position_de{$position}->{'Nmoy'}
                        += $position_de{$position}->{$jour}{'N'};
                $position_de{$position}->{'Emoy'}
                        += $position_de{$position}->{$jour}{'E'};
                $position_de{$position}->{'Umoy'}
                        += $position_de{$position}->{$jour}{'U'};
                $position_de{$position}->{'Lmoy'}
                        += $position_de{$position}->{$jour}{'L'};
                $position_de{$position}->{'nb_jours'} ++;
            }
            $position_de{$position}->{'Nmoy'} /= $position_de{$position}->{'nb_jours'};
            $position_de{$position}->{'Emoy'} /= $position_de{$position}->{'nb_jours'};
            $position_de{$position}->{'Umoy'} /= $position_de{$position}->{'nb_jours'};
            $position_de{$position}->{'Lmoy'} /= $position_de{$position}->{'nb_jours'};
        }

# remplace les positions par les ecarts a la position moyenne
        foreach my $position (keys(%position_de)){
            LOOP2_ON_DAYS:
            foreach my $jour (@jours){
                next LOOP2_ON_DAYS
                    if ( !exists( $position_de{$position}->{$jour} ) );
                $position_de{$position}->{$jour}{'N'}
                        -= $position_de{$position}->{'Nmoy'};
                $position_de{$position}->{$jour}{'E'}
                        -= $position_de{$position}->{'Emoy'};
                $position_de{$position}->{$jour}{'U'}
                        -= $position_de{$position}->{'Umoy'};
                $position_de{$position}->{$jour}{'L'}
                        -= $position_de{$position}->{'Lmoy'};
            }
        }                
# test
        print '            Nord                 Est        '
                   ."       Up                   Long         jour\n";
        foreach my $position (sort( keys( %position_de ) ) ){
            printf
            "IG%2s      %7.4f              %7.4f              "
                      ."%7.4f              %7.4f\n",
                                $position,
                                $position_de{$position}->{'Nmoy'},
                                $position_de{$position}->{'Emoy'},
                                $position_de{$position}->{'Umoy'},
                                $position_de{$position}->{'Lmoy'};
            LOOP3_ON_DAYS:
            foreach my $jour (@jours){
                next LOOP3_ON_DAYS
                    if ( !exists( $position_de{$position}->{$jour} ) );
                printf "    %7.4f +/- %6.4f   %7.4f +/- %6.4f   "
                          ."%7.4f +/- %6.4f   %7.4f +/- %6.4f   %3s\n",
                        $position_de{$position}->{$jour}{'N'},
                        $position_de{$position}->{$jour}{'Nrms'},
                        $position_de{$position}->{$jour}{'E'},
                        $position_de{$position}->{$jour}{'Erms'},
                        $position_de{$position}->{$jour}{'U'},
                        $position_de{$position}->{$jour}{'Urms'},
                        $position_de{$position}->{$jour}{'L'},
                        $position_de{$position}->{$jour}{'Lrms'},
                        $jour;
            }
        }
        print "\n";
        $jours_ref = \@jours;   # mémorise le tableau des jours
        $resultats_pour{$pivot}{$freq} = \%position_de;
    }   # Frequences
}   # pivot de reference
#########################################################
# Plot des series temporelles par pivot et par composante
#########################################################
foreach my $composante ( 'N','E','U','L' ){
    print STDERR  "prepare serie temporelle sur $composante";
    foreach my $pivot ( keys %resultats_pour ){
        print STDERR  " par rapport a $pivot";
        open( PSTEXT,
            "| pstext -JX29.7/21. -R0./29.7/0./21. -X0 -Y0 -K > ${pivot}_${composante}_${FIC_SERIE_TEMP}"
            );
        print PSTEXT "14.8 20.2 22 0 4 MC Ecarts \344 la position moyenne\n";
        print PSTEXT "14.8 19.5 15 0 4 MC r\345f\345rence ${pivot} composante ${composante}\n";
        close PSTEXT;
        foreach my $freq ( sort keys %{$resultats_pour{$pivot}} ){
            print STDERR  " pour $freq:\n"; 
            my %position_de = %{$resultats_pour{$pivot}->{$freq}};
            my @jours = @$jours_ref;
            my $jour_debut = $jours[0] - 1;
            my $jour_fin = $jours[-1] + 1;
            my $yof = 1.7;
            open( PSTEXT,
                "| pstext -JX -R -X1.5 -G0/0/0 -K -O >> ${pivot}_${composante}_${FIC_SERIE_TEMP}"
                );
            print PSTEXT "6.5175 18.5 12 0 4 MC ${freq}\n";
            print STDERR  "6.5175 18.5 12 0 4 MC ${freq}\n";
            close PSTEXT;
            my $numero_position = 0;
            my $annotation = 'WSne';
            foreach my $position ( sort keys(%position_de)){
                open( PSBASEMAP,
                    "| psbasemap -Ba5f1/a5f1:'${position} (mm)':${annotation} -JX13.035/1.7 "
                    ."-R$jour_debut/$jour_fin/-5/5 -Y$yof -K -O >> ${pivot}_${composante}_${FIC_SERIE_TEMP}");
                close PSBASEMAP;
                $annotation =~ tr/S/s/ if ($numero_position == 0);  # under the first line only
                $annotation =~ tr/WwEe/wWeE/;   # vertical annotation once on left, once on right
                $annotation =~ tr/n/N/ if ($position == '09');  # under the first line only
                # plot the positions with the error bar
                open( PSXY,
                    "| psxy -B -JX -R -Ey.1 -Sd.1 -G0/0/200 -K -O  >> ${pivot}_${composante}_${FIC_SERIE_TEMP}"
                    );
                LOOP4_ON_DAYS:
                foreach my $jour (@jours){
                    next LOOP4_ON_DAYS
                        if ( !defined( $position_de{$position}->{$jour}{$composante} ) );
                    printf PSXY "%d %f %f\n",
                            $jour,
                            $position_de{$position}->{$jour}{$composante} * 1000,
                            $position_de{$position}->{$jour}{$composante.'rms'} * 1000;
                }
                close PSXY;
                # draw a broken line between those points
                open( PSXY,
                    "| psxy -B -JX -R -K -O >> ${pivot}_${composante}_${FIC_SERIE_TEMP}"
                    );
                LOOP5_ON_DAYS:
                foreach my $jour (@jours){
                    next LOOP5_ON_DAYS
                        if ( !exists( $position_de{$position}->{$jour}{$composante} ) );
                    printf PSXY "%d %f\n",
                        $jour,
                        $position_de{$position}->{$jour}{$composante} * 1000;
                }
                close PSXY;
                open( PSTEXT,
                    "| pstext -JX -R -K -O >> ${pivot}_${composante}_${FIC_SERIE_TEMP}"
                    );
                my $antenne_precedente = '';
                # writes the name of the antenna
                LOOP_ON_DAYS1:
                foreach my $jour1 ( sort {$a<=>$b} keys %{$antenne_en{$position}} ){
                    next LOOP_ON_DAYS1
                        if ($antenne_precedente eq $antenne_en{$position}->{$jour1});
                    print PSTEXT 
                        "$jour1 3 08 0 4 TL $antenne_en{$position}->{$jour1}\n";
                    $antenne_precedente = $antenne_en{$position}->{$jour1};
                }
                close PSTEXT;    
                $numero_position++;
            }   # position
            my $recalage_vertical = $yof * $numero_position;
            print STDERR  "recalage vertical = $recalage_vertical\n";
            open( PSTEXT,                              # recalage du second graphique dans la même page
                "| pstext -JX14.535/21 -R0./14.535/0./21 -X13.535 -Y-${recalage_vertical} -K -O >> ${pivot}_${composante}_${FIC_SERIE_TEMP}"
                );
            close PSTEXT;
        }       # frequence
        open( PSTEXT,                                   # cloture du fichier
            "| pstext -JX -R -O >> ${pivot}_${composante}_${FIC_SERIE_TEMP}"
            );
        close PSTEXT;
    }           # pivot
}               # composante
