#!/opt/perl/bin/perl
use GPS::IGS::Sum;
use GPS::Defaults;
use File::Basename;
use File::Path;
use File::Copy;
use Cwd;
use GPS::DATES;
my $progname = basename($0);
my $usage=qq(
  $progname provides information about the daily fit of the satellites
  resulting from the weighted average combination by the IGS.
  
  $progname date orb_type [sat_flag] [sat_id]
  
  Options:
   date: date of the combination (required)
   orb_type: rapid|final (required)
   sat_flag: (optional) [a|x|s|e|m] show selected satellites:
             a: all
             x: excluded from the RMS,WRMS combination (default)
             e: eclipsing by Earth
             m: eclipsing by Moon
             s: eclipsing by Earth or Moon (sats in Shadow file)
   sat_id: prn or gps (default)
   
);
die $usage if ! @ARGV;
my $date = get_date(shift,'cal') ;
my $orb = shift || '';
my $ac = shift || 'igs';
if ( ! $date || $orb !~ /(rapid|final)/i ) {
 die $usage;
}
my $w = get_date($date,'gpsw');
my $d = get_date($date,'gdow');
my %par = get_defaults();
my $cal = get_date($date,'cal');
my $flags = lc shift || 'x';
my $pg = lc shift || 'gps';

if ( $pg ne 'gps' && $pg ne 'prn' ) {
  die "Invalid sat_id:$pg\n$usage";
}
my $AC = uc $ac;
if ( $ac == 'ig1') {
  $w = sprintf("%04d",$w);
}
my $dir = $par{$AC."_ORBITS_DIR"}."/$w";

my $file = $orb eq 'rapid' ? "igr$w$d.sum" : "$ac${w}7.sum";
my $here = cwd;
if ( -e "$par{CACHE_DIR}/orbits/$ac/$cal/$file" ){
  copy ("$par{CACHE_DIR}/orbits/$ac/$cal/$file",$here);
}
elsif ( -e "$dir/$file.Z" ) {
  `gzip -dc $dir/$file.Z > $file`;
  if ( ! -d "$par{CACHE_DIR}/orbits/igs/$cal" ) {
    mkpath "$par{CACHE_DIR}/orbits/igs/$cal";
  }
  copy $file,"$par{CACHE_DIR}/orbits/igs/$cal";
}

@flags = split(/,/,uc $flags);
$flags = join('|',@flags);
#$flags = $flags =~ /(s|e|m)/i ? 'E|M' : 
#         $flags =~ /x/i ? 'X' : '.*';

my $s = GPS::IGS::Sum->new($file);
my %seen=();
my @sats = grep { ! $seen{$_}++ }
           map { $_->{$pg} } 
           grep { $_->{flags} =~ /$flags/i } $s->sats($date);
%seen=();
my @h_sats =  grep { ! $seen{$_}++ }
              map { [$_->{$pg},$_->{rms}] } 
              grep { $_->{flags} !~ /X/i } $s->sats($date);
my $m = median( map { $_->[1] } @h_sats );           
@h_sats = map { $_->[0] } grep { ($_->[1]/$m - 1)*100 >= 60 } @h_sats;
#push @sats,@h_sats if @h_sats;
print join (',',@sats),"\n";
#print join (',',@h_sats),"\n";


sub median {
  my @v = sort { $a <=> $b } @_;
  my $n = @v;
  my $half= int(@v/2);
  return  @v % 2 ? $v[$half] : ($v[$half-1]+$v[$half])/2;
}
