| 1 | #!/usr/bin/env perl
|
|---|
| 2 | # ========================================================================
|
|---|
| 3 | # compareSP3Plot.pl
|
|---|
| 4 | # ========================================================================
|
|---|
| 5 | #
|
|---|
| 6 | # Purpose: Script using gnuplot to plot BNC's SP3 Comparison Results
|
|---|
| 7 | #
|
|---|
| 8 | # Comment: stored values:
|
|---|
| 9 | # To compare satellite clocks provided by the two files, BNC first converts
|
|---|
| 10 | # coordinate differences dX,dY,dZ into along track, out-of-plane, and radial
|
|---|
| 11 | # components.
|
|---|
| 12 | # It then corrects the clock differences for the radial components of
|
|---|
| 13 | # coordinate differences. RMS values of clock differences are
|
|---|
| 14 | # finally calculated after introducing
|
|---|
| 15 | # - at first one offset 'per epoch for all satellites' and
|
|---|
| 16 | # - secondly one offset 'per satellite for all epochs'.
|
|---|
| 17 | # ========================================================================
|
|---|
| 18 | # Uses
|
|---|
| 19 | use strict;
|
|---|
| 20 | use warnings;
|
|---|
| 21 | use Getopt::Long;
|
|---|
| 22 | use Chart::Gnuplot;
|
|---|
| 23 | use Data::Dumper qw(Dumper);
|
|---|
| 24 | use diagnostics;
|
|---|
| 25 | use File::Basename;
|
|---|
| 26 | use Date::Manip;
|
|---|
| 27 | use Time::Piece 1.30;
|
|---|
| 28 | use List::MoreUtils qw( minmax );
|
|---|
| 29 | use experimental 'smartmatch';
|
|---|
| 30 |
|
|---|
| 31 | # -----------------------------------------------------------------------------
|
|---|
| 32 | # required to create a pdf
|
|---|
| 33 | # -----------------------------------------------------------------------------
|
|---|
| 34 | use PDF::API2;
|
|---|
| 35 | use constant {mm => 25.4 / 72, inch => 1 / 72,pt => 1,
|
|---|
| 36 | };# There are 72 postscript points in an inch and there are 25.4 millimeters in an inch.
|
|---|
| 37 | # -----------------------------------------------------------------------------
|
|---|
| 38 | # Options
|
|---|
| 39 | # -----------------------------------------------------------------------------
|
|---|
| 40 | # init options
|
|---|
| 41 | my $help = 0;
|
|---|
| 42 | my @plotTypes = ();
|
|---|
| 43 | my @logFiles = ();
|
|---|
| 44 | # read options from command line
|
|---|
| 45 | GetOptions(
|
|---|
| 46 | 'help' => \$help,
|
|---|
| 47 | 'plotTypes=s' => \@plotTypes,
|
|---|
| 48 | 'logFiles=s' => \@logFiles,
|
|---|
| 49 | );
|
|---|
| 50 | help() if $help;
|
|---|
| 51 | @plotTypes = qw(NEU) unless (@plotTypes);# Default values
|
|---|
| 52 | @plotTypes = map {uc} split ( /[, ]/, join ( ',', @plotTypes ) );# Reformat options
|
|---|
| 53 | @logFiles = split ( /[, ]/, join ( ',', @logFiles ) );
|
|---|
| 54 | unless (@logFiles) {
|
|---|
| 55 | print "ERROR: logfiles missing\n";
|
|---|
| 56 | exit;
|
|---|
| 57 | }
|
|---|
| 58 | print " plotTpes: @plotTypes\n logfiles: @logFiles\n\n";
|
|---|
| 59 | # -----------------------------------------------------------------------------
|
|---|
| 60 | # generate data sets for gnuplot
|
|---|
| 61 | # -----------------------------------------------------------------------------
|
|---|
| 62 | my $epoSec = 0;
|
|---|
| 63 | # for pdf generation
|
|---|
| 64 | my ($png, $page, $headline, $headline_text);
|
|---|
| 65 | my $y0 = 190/mm;
|
|---|
| 66 | my ($x, $y, $width, $height) = (40/mm, $y0, 130/mm, 70/mm);
|
|---|
| 67 | my $dy = $height+7/mm;
|
|---|
| 68 | # file by file ..
|
|---|
| 69 | foreach my $file (@logFiles) {
|
|---|
| 70 | my (%COMPARISON);
|
|---|
| 71 | my (@array);
|
|---|
| 72 | my ($sys, $sat, $epo, $radial, $along, $out, $drange, $clk, $clkRed);
|
|---|
| 73 | # -----------------------------------------------------------------------------
|
|---|
| 74 | # create pdf for plot results
|
|---|
| 75 | # -----------------------------------------------------------------------------
|
|---|
| 76 | my($inputFilename, $inputDir, $inputSuffix) = fileparse($file, '\..*');
|
|---|
| 77 | my $pdf_name = sprintf("%s.pdf", $inputFilename);
|
|---|
| 78 | my $pdf = PDF::API2->new(-file => "$inputDir$pdf_name");
|
|---|
| 79 | my $font1 = $pdf->corefont('Helvetica-Bold');
|
|---|
| 80 | # -----------------------------------------------------------------------------
|
|---|
| 81 | # read file and save data
|
|---|
| 82 | # -----------------------------------------------------------------------------
|
|---|
| 83 | open ( INPUT, "<$file" ) || die $!;
|
|---|
| 84 | print "Parse logfile $file ...\n";
|
|---|
| 85 | while (<INPUT>) {#56985.000000 G01 0.0160 -0.0242 0.0961 0.0082 -0.0077 1
|
|---|
| 86 | if ( $_ !~ /^!/) { #print "$_\n";
|
|---|
| 87 | @array = split(/\s+/, $_);
|
|---|
| 88 | $epo = $array[0]; #print "$epo\n";
|
|---|
| 89 | $epo =~s/(\.\d*)//; #print "$epo\n";
|
|---|
| 90 | $sys = substr $array[1], 0, 1;
|
|---|
| 91 | $sat = $array[1];#print "$sat\n";
|
|---|
| 92 | $radial = $array[2];#print "$radial\n";
|
|---|
| 93 | $along = $array[3];#print "$along\n";
|
|---|
| 94 | $out = $array[4];#print "$out\n";
|
|---|
| 95 | $clk = $array[5];#print "$clk\n";
|
|---|
| 96 | $clkRed = $array[6];#print "$clkRed\n";
|
|---|
| 97 | if ($out !~ /^\.$/) {
|
|---|
| 98 | push @{$COMPARISON{$sys}{$sat}{EPOCH }}, $epo;
|
|---|
| 99 | push @{$COMPARISON{$sys}{$sat}{RADIAL}}, $radial;
|
|---|
| 100 | push @{$COMPARISON{$sys}{$sat}{ALONG }}, $along;
|
|---|
| 101 | push @{$COMPARISON{$sys}{$sat}{OUT }}, $out;
|
|---|
| 102 | push @{$COMPARISON{$sys}{$sat}{CLK }}, $clk;
|
|---|
| 103 | push @{$COMPARISON{$sys}{$sat}{CLKRED}}, $clkRed;
|
|---|
| 104 | }
|
|---|
| 105 | }
|
|---|
| 106 | }
|
|---|
| 107 | close INPUT;
|
|---|
| 108 |
|
|---|
| 109 | # -----------------------------------------------------------------------------
|
|---|
| 110 | # plot several data sets
|
|---|
| 111 | # -----------------------------------------------------------------------------
|
|---|
| 112 | ######### RAO #####################
|
|---|
| 113 | if (("ALL" ~~ @plotTypes ) || ("RAO" ~~ @plotTypes)) {
|
|---|
| 114 | print "Plot RADIAL, ALONG TRACK and OUT OF PLANE COMPONENTS ...\n";
|
|---|
| 115 | my $yrange = .6;
|
|---|
| 116 | $page = $pdf->page(); $page->mediabox('A4');
|
|---|
| 117 | $headline = sprintf("Orbit Differences in radial, along-track and out-of-plane component (%s)", $inputFilename);
|
|---|
| 118 | $headline_text = $page->text;
|
|---|
| 119 | $headline_text->font($font1, 11/pt);
|
|---|
| 120 | $headline_text->translate(15/mm, 280/mm);
|
|---|
| 121 | $headline_text->text($headline);
|
|---|
| 122 | $y = $y0+$dy;
|
|---|
| 123 | #SYSTEM
|
|---|
| 124 | foreach my $key_sys (sort keys %COMPARISON) {
|
|---|
| 125 | print " for $key_sys\n";
|
|---|
| 126 | my (@datasets_r, @datasets_a, @datasets_o); # init datasets
|
|---|
| 127 | my $pngName_r = sprintf ("%s_RADIAL_%s.png", $inputFilename, $key_sys);
|
|---|
| 128 | my $pngName_a = sprintf ("%s_ALONG_%s.png", $inputFilename, $key_sys);
|
|---|
| 129 | my $pngName_o = sprintf ("%s_OUT_%s.png", $inputFilename, $key_sys);
|
|---|
| 130 | my $chart_r = Chart::Gnuplot->new(
|
|---|
| 131 | output => $pngName_r, terminal => 'png', title => "Radial component",
|
|---|
| 132 | ylabel => "Orbit Difference [m]", yrange => ["-$yrange", "$yrange"],
|
|---|
| 133 | xlabel => "Time [h]", timeaxis => 'x', xtics => {labelfmt => '%H:%M', rotate => '-270',},
|
|---|
| 134 | legend => { position => "outside right",},
|
|---|
| 135 | grid => 'on',);
|
|---|
| 136 | my $chart_a = Chart::Gnuplot->new(
|
|---|
| 137 | output => $pngName_a, terminal => 'png', title => "Along track component",
|
|---|
| 138 | ylabel => "Orbit Difference [m]", yrange => ["-$yrange", "$yrange"],
|
|---|
| 139 | xlabel => "Time [h]", timeaxis => 'x', xtics => {labelfmt => '%H:%M', rotate => '-270',},
|
|---|
| 140 | legend => { position => "outside right",},
|
|---|
| 141 | grid => 'on',);
|
|---|
| 142 | my $chart_o = Chart::Gnuplot->new(
|
|---|
| 143 | output => $pngName_o, terminal => 'png', title => "Out of plane component",
|
|---|
| 144 | ylabel => "Orbit Difference [m]", yrange => ["-$yrange", "$yrange"],
|
|---|
| 145 | xlabel => "Time [h]", timeaxis => 'x', xtics => {labelfmt => '%H:%M', rotate => '-270',},
|
|---|
| 146 | legend => { position => "outside right",},
|
|---|
| 147 | grid => 'on',);
|
|---|
| 148 | #SATELLITE
|
|---|
| 149 | foreach my $key_sat (sort keys %{$COMPARISON{$key_sys}}) {
|
|---|
| 150 | my $colour = sprintf("%d", (substr $key_sat, -2)+1);
|
|---|
| 151 | #RADIAL
|
|---|
| 152 | my $dataset_r = Chart::Gnuplot::DataSet->new(
|
|---|
| 153 | xdata => \@{$COMPARISON{$key_sys}{$key_sat}{EPOCH}},
|
|---|
| 154 | ydata => \@{$COMPARISON{$key_sys}{$key_sat}{RADIAL}},
|
|---|
| 155 | title => "$key_sat",
|
|---|
| 156 | timefmt => '%Y-%m-%d_%H:%M:%S',
|
|---|
| 157 | linetype => $colour,
|
|---|
| 158 | #style => " linespoints ",);
|
|---|
| 159 | style => " lines ",);
|
|---|
| 160 | push (@datasets_r, $dataset_r);
|
|---|
| 161 | # ALONG TRACK
|
|---|
| 162 | my $dataset_a = Chart::Gnuplot::DataSet->new(
|
|---|
| 163 | xdata => \@{$COMPARISON{$key_sys}{$key_sat}{EPOCH}},
|
|---|
| 164 | ydata => \@{$COMPARISON{$key_sys}{$key_sat}{ALONG}},
|
|---|
| 165 | title => "$key_sat",
|
|---|
| 166 | timefmt => '%Y-%m-%d_%H:%M:%S',
|
|---|
| 167 | linetype => $colour,
|
|---|
| 168 | style => " lines ",);
|
|---|
| 169 | push (@datasets_a, $dataset_a);
|
|---|
| 170 | # OUT OF PLANE
|
|---|
| 171 | my $dataset_o = Chart::Gnuplot::DataSet->new(
|
|---|
| 172 | xdata => \@{$COMPARISON{$key_sys}{$key_sat}{EPOCH}},
|
|---|
| 173 | ydata => \@{$COMPARISON{$key_sys}{$key_sat}{OUT}},
|
|---|
| 174 | title => "$key_sat",
|
|---|
| 175 | timefmt => '%Y-%m-%d_%H:%M:%S',
|
|---|
| 176 | linetype => $colour,
|
|---|
| 177 | style => " lines",);
|
|---|
| 178 | push (@datasets_o, $dataset_o);
|
|---|
| 179 | }
|
|---|
| 180 | $chart_r->plot2d(@datasets_r);
|
|---|
| 181 | $chart_a->plot2d(@datasets_a);
|
|---|
| 182 | $chart_o->plot2d(@datasets_o);
|
|---|
| 183 | # RADIAL
|
|---|
| 184 | $y = $y - $dy;
|
|---|
| 185 | if($y < 30/mm) {
|
|---|
| 186 | $page = $pdf->page(); $page->mediabox('A4');
|
|---|
| 187 | $y = $y0;
|
|---|
| 188 | }
|
|---|
| 189 | $png = $page->gfx(); die("Unable to find image file: $!") unless -e $pngName_r;
|
|---|
| 190 | $png->image($pdf->image_png($pngName_r), $x, $y, $width, $height);
|
|---|
| 191 | # ALONG
|
|---|
| 192 | $y = $y - $dy;
|
|---|
| 193 | if($y < 30/mm) {
|
|---|
| 194 | $page = $pdf->page(); $page->mediabox('A4');
|
|---|
| 195 | $y = $y0;
|
|---|
| 196 | }
|
|---|
| 197 | $png = $page->gfx(); die("Unable to find image file: $!") unless -e $pngName_a;
|
|---|
| 198 | $png->image($pdf->image_png($pngName_a), $x, $y, $width, $height);
|
|---|
| 199 | # OUT
|
|---|
| 200 | $y = $y - $dy;
|
|---|
| 201 | if($y < 30/mm) {
|
|---|
| 202 | $page = $pdf->page(); $page->mediabox('A4');
|
|---|
| 203 | $y = $y0;
|
|---|
| 204 | }
|
|---|
| 205 | $png = $page->gfx(); die("Unable to find image file: $!") unless -e $pngName_o;
|
|---|
| 206 | $png->image($pdf->image_png($pngName_o), $x, $y, $width, $height);
|
|---|
| 207 | }
|
|---|
| 208 | }
|
|---|
| 209 | ######## CLK / CLKRED #####################
|
|---|
| 210 | if (("ALL" ~~ @plotTypes ) || ("CLK" ~~ @plotTypes)) {
|
|---|
| 211 | print "Plot CLOCK and CLOCK REDUCED BY RADIAL COMPONENT ...\n";
|
|---|
| 212 | my $yrange = 0.6;
|
|---|
| 213 | $page = $pdf->page(); $page->mediabox('A4');
|
|---|
| 214 | $headline = sprintf("Clock residuals and reduced clock residuals(%s)", $inputFilename);
|
|---|
| 215 | $headline_text = $page->text;
|
|---|
| 216 | $headline_text->font($font1, 11/pt);
|
|---|
| 217 | $headline_text->translate(15/mm, 280/mm);
|
|---|
| 218 | $headline_text->text($headline);
|
|---|
| 219 | $y = $y0+$dy;
|
|---|
| 220 | #SYSTEM
|
|---|
| 221 | foreach my $key_sys (sort keys %COMPARISON) {
|
|---|
| 222 | print " for $key_sys\n";
|
|---|
| 223 | my (@datasets_c, @datasets_r); # init datasets
|
|---|
| 224 | my $pngName_c = sprintf ("%s_CLK_%s.png", $inputFilename, $key_sys);
|
|---|
| 225 | my $pngName_r = sprintf ("%s_CLKRED_%s.png", $inputFilename, $key_sys);
|
|---|
| 226 | my $chart_c = Chart::Gnuplot->new(
|
|---|
| 227 | output => $pngName_c, terminal => 'png',title => "Clock only",
|
|---|
| 228 | ylabel => "Clock Residuals [m]", yrange => ["-$yrange", "$yrange"],
|
|---|
| 229 | xlabel => "Time [h]", timeaxis => 'x', xtics => {labelfmt => '%H:%M', rotate => '-270',},
|
|---|
| 230 | legend => { position => "outside right",},
|
|---|
| 231 | grid => 'on',);
|
|---|
| 232 | my $chart_r = Chart::Gnuplot->new(
|
|---|
| 233 | output => $pngName_r, terminal => 'png',title => "Clock minus radial component",
|
|---|
| 234 | ylabel => "Clock Residuals [m]", yrange => ["-$yrange", "$yrange"],
|
|---|
| 235 | xlabel => "Time [h]", timeaxis => 'x', xtics => {labelfmt => '%H:%M', rotate => '-270',},
|
|---|
| 236 | legend => { position => "outside right",},
|
|---|
| 237 | grid => 'on',);
|
|---|
| 238 | #SATELLITE
|
|---|
| 239 | foreach my $key_sat (sort keys %{$COMPARISON{$key_sys}}) {
|
|---|
| 240 | my $colour = sprintf("%d", (substr $key_sat, -2)+1);
|
|---|
| 241 | # CLK
|
|---|
| 242 | my $dataset_c = Chart::Gnuplot::DataSet->new(
|
|---|
| 243 | xdata => \@{$COMPARISON{$key_sys}{$key_sat}{EPOCH}},
|
|---|
| 244 | ydata => \@{$COMPARISON{$key_sys}{$key_sat}{CLK}},
|
|---|
| 245 | title => "$key_sat",
|
|---|
| 246 | timefmt => '%Y-%m-%d_%H:%M:%S',
|
|---|
| 247 | linetype => $colour,
|
|---|
| 248 | style => " lines ",);
|
|---|
| 249 | push (@datasets_c, $dataset_c);
|
|---|
| 250 | my ( $min, $max ) = minmax @{$COMPARISON{$key_sys}{$key_sat}{CLK}};
|
|---|
| 251 | if ($max > 1.0) {print "$key_sat: max clk = $max\n"};
|
|---|
| 252 | # CLKRED
|
|---|
| 253 | my $dataset_r = Chart::Gnuplot::DataSet->new(
|
|---|
| 254 | xdata => \@{$COMPARISON{$key_sys}{$key_sat}{EPOCH}},
|
|---|
| 255 | ydata => \@{$COMPARISON{$key_sys}{$key_sat}{CLKRED}},
|
|---|
| 256 | title => "$key_sat",
|
|---|
| 257 | timefmt => '%Y-%m-%d_%H:%M:%S',
|
|---|
| 258 | linetype => $colour,
|
|---|
| 259 | style => " lines ",);
|
|---|
| 260 | push (@datasets_r, $dataset_r);
|
|---|
| 261 | }
|
|---|
| 262 | $chart_c->plot2d(@datasets_c);
|
|---|
| 263 | $chart_r->plot2d(@datasets_r);
|
|---|
| 264 | # CLK
|
|---|
| 265 | $y = $y - $dy;
|
|---|
| 266 | if($y < 30/mm) {
|
|---|
| 267 | $page = $pdf->page(); $page->mediabox('A4');
|
|---|
| 268 | $y = $y0;
|
|---|
| 269 | }
|
|---|
| 270 | $png = $page->gfx(); die("Unable to find image file: $!") unless -e $pngName_c;
|
|---|
| 271 | $png->image($pdf->image_png($pngName_c), $x, $y, $width, $height);
|
|---|
| 272 | # CLKRED
|
|---|
| 273 | $y = $y - $dy;
|
|---|
| 274 | if($y < 30/mm) {
|
|---|
| 275 | $page = $pdf->page(); $page->mediabox('A4');
|
|---|
| 276 | $y = $y0;
|
|---|
| 277 | }
|
|---|
| 278 | $png = $page->gfx(); die("Unable to find image file: $!") unless -e $pngName_r;
|
|---|
| 279 | $png->image($pdf->image_png($pngName_r), $x, $y, $width, $height);
|
|---|
| 280 | }
|
|---|
| 281 | }
|
|---|
| 282 | $pdf->save();
|
|---|
| 283 | $pdf->end();
|
|---|
| 284 | system("rm *png");
|
|---|
| 285 | #system("evince $inputDir/$pdf_name&");
|
|---|
| 286 | }
|
|---|
| 287 |
|
|---|
| 288 |
|
|---|
| 289 | #########################################
|
|---|
| 290 | sub help {
|
|---|
| 291 | print <<EOI_HILFE;
|
|---|
| 292 | CALL:
|
|---|
| 293 | compareSp3Plot.pl [paramter]
|
|---|
| 294 |
|
|---|
| 295 | PARAMETER:
|
|---|
| 296 |
|
|---|
| 297 | --plotTypes possible: ALL => RAO and CLOCK; CLK => CLOCK and CLOCKRED; RAO => RAO only
|
|---|
| 298 | --logFiles comma separated list of BNC's combination logfiles
|
|---|
| 299 | --help show help contents
|
|---|
| 300 |
|
|---|
| 301 | DESCRIPTION:
|
|---|
| 302 |
|
|---|
| 303 | Script to plot BNC's SP3 Comparison Results
|
|---|
| 304 |
|
|---|
| 305 | EOI_HILFE
|
|---|
| 306 | exit;
|
|---|
| 307 | }
|
|---|
