Документ взят из кэша поисковой машины. Адрес оригинального документа : http://www.naic.edu/~phil/hardware/pdev/fpga/gx/sp/test/stokesplot Дата изменения: Thu Jun 26 04:27:55 2008 Дата индексирования: Sat Sep 6 20:05:52 2008 Кодировка: Поисковые слова: http www.astronomy.ru forum index.php topic 4644.0.html

#!/usr/bin/perl

# Jeff Mock
# 2030 Gough St.
# San Francisco, CA 94109
# jeff@mock.com
#
# Copyright 2005,2006
#
# $URL: https://www.mock.com/svn/pdev/trunk/gx/sp/test/stokesplot $
# $Id: stokesplot 177 2006-06-16 00:14:54Z jeff $

use Getopt::Long;
use Math::BigInt;

$opt_width = 18; # datapath width
$opt_n = 32; # length of fft
$opt_start = undef;
$opt_stop = undef;
$opt_fn = "plot";
$opt_debug = 1;
$opt_pack = 0;

$0 =~ /(.*)\/.*/;
$opt_dir = $1 eq "" ? "." : $1;

%opts = (
'width=o' => \$opt_width,
'n=o' => \$opt_n,
'start=f' => \$opt_start,
'stop=f' => \$opt_stop,
'fn=s' => \$opt_fn,
'debug' => \$opt_debug,
'pack' => \$opt_pack,
);

if (!GetOptions(%opts)) {
print STDERR "
Plot accumulation output of gx chip

stokesplot [options]
[--n=n] Length of FFT ($opt_n)
[--width=n] Datapath width of FFT ($opt_width)
[--start=f] Start frequency ($opt_start) [-n/2 .. n/2]
[--stop=f] Stop frequency ($opt_stop) [-n/2 .. n/2]
[--fn=s] Output filename for plot
[--debug] Leave intermediate files around
[--pack] Use dma.out file for plotting

\n";
exit 1;
}

sub bitrev {
my $n = shift;
my $bits = shift;

my $u = $n & 1;
while(--$bits) {
$n >>= 1;
$u <<= 1;
$u += $n & 1;
}
return $u;
}

sub log2 {
my $v = shift;
return int(log(2*$v-1)/log(2));
}

sub readpfb {
my $fn = shift;

my $wsize = Math::BigInt->bone();
$wsize->blsft($opt_width);
my $wsizem1 = $wsize->copy;
$wsizem1->bdec();
# printf "wsizem1 is %s\n", $wsizem1->as_hex;

my $highbit = Math::BigInt->bone; (1 << ($opt_width-1));
$highbit->blsft($opt_width-1);

my $highbitp1 = Math::BigInt->bone; (1 << ($opt_width));
$highbitp1->blsft($opt_width);

my $bits = log2($opt_n)-1;
my @retlist = ();

open $fd, "< $fn" or die "Cannot open $fn: $!";
while (<$fd>) {
my @s0 = ();
my @s1 = ();
my @s2 = ();
my @s3 = ();
my %acc = ();
for my $i (0 .. $opt_n/2-1) {
die "Premature end of pfb file" unless ($line=<$fd>);
my @f = split(' ', $line);
die "Bad Line marker $f[0]" unless $f[0] == $i;

my $col=0;
for $v (@f) {
die "x's found in data..."
if $v =~ /x/;
# $v = hex($v) & ($wsize-1);
# $v = $v - $wsize if $v >= $highbit && $opt_width<32;
# $v = ($v + 0.0) / ($highbit + 0.0);
$v = Math::BigInt->new("0x" . $v);
$v->band($wsizem1);

# columns 1,2,5,6 are unsigned
if ($col==1 || $col==2 || $col==5 || $col==6) {
$v = $v->numify / ($highbitp1->numify + 0.0);
} else {
$v->bsub($wsize) if $v->bcmp($highbit) >= 0;
$v = $v->numify / ($highbit->numify + 0.0);
}
$col++;
}
my $idx = bitrev($i, $bits);

# The X output has all of the positive freq bins
$s0[$idx + $opt_n/2] = $f[1];
$s1[$idx + $opt_n/2] = $f[2];
$s2[$idx + $opt_n/2] = $f[3];
$s3[$idx + $opt_n/2] = $f[4];

# The Y output has all of the negative freq bins
$s0[$idx] = $f[5];
$s1[$idx] = $f[6];
$s2[$idx] = $f[7];
$s3[$idx] = $f[8];
}
$acc{si} = [];
$acc{s0} = \@s0;
$acc{s1} = \@s1;
$acc{s2} = \@s2;
$acc{s3} = \@s3;
$acc{n} = $opt_n;
$acc{bits} = 32;
$acc{type} = 2;
push @retlist, \%acc;
}
close $fd;
return \@retlist;
}

sub myval_us {
my $str = shift;
my $sz = shift;

my $highbitp1 = Math::BigInt->bone;
$highbitp1->blsft($sz);

my $v = Math::BigInt->new("0x" . $str);
$v = $v->numify / ($highbitp1->numify + 0.0);
return $v;
}

sub myval_s {
my $str = shift;
my $sz = shift;

my $highbit = Math::BigInt->bone;
$highbit->blsft($sz-1);
my $wsize = Math::BigInt->bone();
$wsize->blsft($sz);

my $v = Math::BigInt->new("0x" . $str);
$v->bsub($wsize) if $v->bcmp($highbit) >= 0;
$v = $v->numify / ($highbit->numify + 0.0);
return $v;
}

sub readpack {
my $fn = shift;
my @retlist = ();

open $fd, "< $fn" or die "Cannot open $fn: $!";
($sig, $pack_n, $pack_wid, $pack_type, $pack_dumpstart, $pack_dumpstop)
= split ' ', <$fd>;
die "bad signature $sig for dma.out" unless $sig eq "PACK";

my $bits = log2($pack_n)-1;

$pack_dumpstart = 0 if $pack_dumpstart<0;
$pack_dumpstop = $pack_n-1 if $pack_dumpstop>=$pack_n;

my $pack_act = ($pack_dumpstop-$pack_dumpstart+1);

# Number of lines depends on format and size
#
if ($pack_wid==32) {
if ($pack_type == 0) { # I
$pack_ni = int(($pack_act+1)/2);
} elsif ($pack_type == 1) { # s0, s1
$pack_ni = $pack_act;
} elsif ($pack_type == 2) {
$pack_ni = $pack_act*2;
}
} elsif ($pack_wid==16) {
if ($pack_type == 0) {
$pack_ni = int(($pack_act+3)/4);
} elsif ($pack_type == 1) {
$pack_ni = int(($pack_act+1)/2);
} elsif ($pack_type == 2) {
$pack_ni = $pack_act;
}
} elsif ($pack_wid==8) {
if ($pack_type==0) {
$pack_ni = int(($pack_act+7)/8);
} elsif ($pack_type==1) {
$pack_ni = int(($pack_act+3)/4);
} elsif ($pack_type==2) {
$pack_ni = int(($pack_act+1)/2);
}
}
$pack_ni++; # For status word at end of each frame

# Map index into packed output to FFT bin number
# First make a bit reverse map as if all bins are dumped.
#
my @rmapx = ();
for (0..$pack_n-1) {
$rmapx[$_] = ($_ < $pack_n/2) ? bitrev($_, $bits) :
bitrev($_ - $pack_n/2, $bits) + $pack_n/2;
}
# Go through list and remove entries for bins not sent
my @rmap = ();
for my $v (@rmapx) {
push @rmap,$v if $v>=$pack_dumpstart && $v<=$pack_dumpstop;
}

while (<$fd>) {
my @si = ();
my @s0 = ();
my @s1 = ();
my @s2 = ();
my @s3 = ();
my %acc = ();
my @upack = ();

for (0..$pack_n-1) {
push @si, 0;
push @s0, 0;
push @s1, 0;
push @s2, 0;
push @s3, 0;
}

# Read packed accumulation into unpack array using
# correct size for items.
#
for my $i (0 .. $pack_ni-1) {
my $v = <$fd>;
# die "x's found in data..." if $v =~ /x/;
($idxc, $pack) = split ' ', $v;
die "Index doesn't compare" if $i != $idxc;
if ($pack_wid==8) {
push @upack, substr($pack, 2*$_, 2) for (0..7);
} elsif ($pack_wid==16) {
push @upack, substr($pack, 4*$_, 4) for (0..3);
} elsif ($pack_wid==32) {
push @upack, substr($pack, 8*$_, 8) for (0..1);
}
}
my $pack_status = $pack;

# Fill si,s0,s1,s2,s3 arrays with values according
# to type and size of dump.
#
# The data is in weirdo bit reverse order. The first n/2
# bins are the negative freq bins in bit reverse order. The
# next n/2 items are the positive bins in bit reverse order.
# So, it's like the high order bit indexing into the unpacked
# area is not bit reverses but the lower n-1 bits are reversed.
#
if ($pack_type==0) {
for my $i (0..$pack_act-1) {
$si[$rmap[$i]] = myval_us($upack[$i], $pack_wid);
}
} elsif ($pack_type==1) {
for my $i (0..$pack_act-1) {
$s0[$rmap[$i]] = myval_us($upack[$i*2+0], $pack_wid);
$s1[$rmap[$i]] = myval_us($upack[$i*2+1], $pack_wid);
}
} elsif ($pack_type==2) {
if ($pack_wid==32) {
# special case for the four pass dump
# In this case, s0,s1 are packed in the first n/2 values
# s2,s3 are packed in the second n/2 values.
#
for my $i (0..$pack_act-1) {
$s0[$rmap[$i]] = myval_us($upack[$i*2+0], $pack_wid);
$s1[$rmap[$i]] = myval_us($upack[$i*2+1], $pack_wid);
$s2[$rmap[$i]] = myval_s ($upack[($i+$pack_act)*2+0],
$pack_wid);
$s3[$rmap[$i]] = myval_s ($upack[($i+$pack_act)*2+1],
$pack_wid);
}
} else {
for my $i (0..$pack_act-1) {
$s0[$rmap[$i]] = myval_us($upack[$i*4+0], $pack_wid);
$s1[$rmap[$i]] = myval_us($upack[$i*4+1], $pack_wid);
$s2[$rmap[$i]] = myval_s ($upack[$i*4+2], $pack_wid);
$s3[$rmap[$i]] = myval_s ($upack[$i*4+3], $pack_wid);
}
}
}

$acc{si} = \@si;
$acc{s0} = \@s0;
$acc{s1} = \@s1;
$acc{s2} = \@s2;
$acc{s3} = \@s3;
$acc{n} = $pack_n;
$acc{bits} = $pack_wid;
$acc{type} = $pack_type;
$acc{info} = $pack_status;
print "Info is $acc{info}\n";
push @retlist, \%acc;
}
close $fd;
return \@retlist;
}

sub doplot {
my $fn = shift;
my $tag = shift;
my $ymin = shift;
my $ymax = shift;
my $info = shift;

open $gp, "> plot.cmd" or die "Cannot open plot.cmd";
# print $gp "set terminal png color\n";
print $gp "set terminal png small\n";
print $gp "set output \"${fn}\"\n";
print $gp "set grid lt 0\n";
printf $gp "set yrange [%f:%f]\n", $ymin, $ymax;
# print $gp "set ylabel \"dB\"\n";
print $gp "set xlabel \" ${d_n} FFT Bins ${tag} \"\n";
# print $gp "set xtics 1\n" if ($stop-$start)<=33;
printf $gp "set xrange [%f:%f]\n", $opt_start, $opt_stop;
print $gp "plot ";
for my $acc (0 .. $accn-1) {
my $col = $acc+2;
my $color = $acc+1;
my $iter = hex(substr($$info[$acc], 12, 4));
my $acc_tot = hex(substr($$info[$acc], 8, 4));
my $err = substr($$info[$acc], 0, 7);
my $cal = hex(substr($$info[$acc], 7, 1)) & 0x1;
my $title = sprintf("Acc %d, int %d", $iter, $acc_tot);
$title .= " cal-on" if $cal;
$title .= " error $err" if $err ne "0000000";
print $gp "\"plot.dat\" using 1:${col} title \"${title}\" with lines ${color} ";
print $gp "," if $acc != $accn-1;
}
print $gp "\n";
close $gp;
die "gnu plot failed" if
system("gnuplot plot.cmd");
}


if ($opt_pack) {
$x = readpack("dma.out");
$accn = scalar @$x;
print "$accn integrations.\n";

$d_type = ${$$x[0]}{type};
$d_n = ${$$x[0]}{n};
$d_bits = ${$$x[0]}{bits};

print "$d_n points, $d_type type.\n";
} else {
$x = readpfb("pfb.out");
$accn = scalar @$x;
# print "$accn integrations.\n";

$ar = $$x[$i];
$s0 = $$ar{s0};
$d_type=2;
$d_n = $opt_n;
$d_bits = 32;
# print "$d_n points.\n";
}

$opt_start = -$d_n/2 if !defined($opt_start);
$opt_stop = $d_n/2-1 if !defined($opt_stop);


# Find ymin/ymax
#
$ymax = 0;
$ymin = 0;
for my $ac (0 .. $accn-1) {
$ar = $$x[$ac];
$si = $$ar{si};
$s0 = $$ar{s0};
$s1 = $$ar{s1};
$s2 = $$ar{s2};
$s3 = $$ar{s3};
if ($d_type==0) {
for my $i (0 .. $d_n-1) {
$stoke_i = $$si[$i];
$ymin = $stoke_i if $stoke_i < $ymin;
$ymax = $stoke_i if $stoke_i > $ymax;
}
} elsif ($d_type==1) {
for my $i (0 .. $d_n-1) {
$stoke_i = ($$s0[$i] + $$s1[$i])/2.0;
$stoke_q = ($$s0[$i] - $$s1[$i])/2.0;
$ymin = $stoke_i if $stoke_i < $ymin;
$ymin = $stoke_q if $stoke_q < $ymin;
$ymax = $stoke_i if $stoke_i > $ymax;
$ymax = $stoke_q if $stoke_q > $ymax;
}
} elsif ($d_type==2) {
for my $i (0 .. $d_n-1) {
$stoke_i = ($$s0[$i] + $$s1[$i])/2.0;
$stoke_q = ($$s0[$i] - $$s1[$i])/2.0;
$stoke_u = $$s2[$i];
$stoke_v = $$s3[$i];
$ymin = $stoke_i if $stoke_i < $ymin;
$ymin = $stoke_q if $stoke_q < $ymin;
$ymin = $stoke_u if $stoke_u < $ymin;
$ymin = $stoke_v if $stoke_v < $ymin;
$ymax = $stoke_i if $stoke_i > $ymax;
$ymax = $stoke_q if $stoke_q > $ymax;
$ymax = $stoke_u if $stoke_u > $ymax;
$ymax = $stoke_v if $stoke_v > $ymax;
}
}
}

# Put 0 in the middle of the y-axis
$ymin = -$ymax if -$ymax < $ymin;
$ymax = -$ymin if -$ymin > $ymax;

# collect the 64-bit info word from the accumulations. This contains
# the sequence number and the error terms.
#
my @info = ();
for my $ac (0 .. $accn-1) {
$ar = $$x[$ac];
push @info, $$ar{info};
}

if ($d_type==0) {
open $gd, "> plot.dat" or die "Cannot open plot.dat";
for my $i (0 .. $d_n-1) {
printf $gd "%d ", $i-$d_n/2;
for my $ac (0 .. $accn-1) {
$ar = $$x[$ac];
$si = $$ar{si};
printf $gd "%g ", $$si[$i];
}
printf $gd "\n";
}
close $gd;
doplot($opt_fn . "-i.png", "I stokes parameter, ${d_bits}-bits",
$ymin, $ymax, \@info);
}

if ($d_type==1 || $d_type==2) {
open $gd, "> plot.dat" or die "Cannot open plot.dat";
for my $i (0 .. $d_n-1) {
printf $gd "%d ", $i-$d_n/2;
for my $ac (0 .. $accn-1) {
$ar = $$x[$ac];
$s0 = $$ar{s0};
$s1 = $$ar{s1};
$s2 = $$ar{s2};
$s3 = $$ar{s3};
printf $gd "%g ", ($$s0[$i] + $$s1[$i])/2.0;
# printf $gd "%g ", $$s0[$i];
}
printf $gd "\n";
}
close $gd;
doplot($opt_fn . "-i.png", "I stokes parameter, ${d_bits}-bits",
$ymin, $ymax, \@info);

open $gd, "> plot.dat" or die "Cannot open plot.dat";
for my $i (0 .. $d_n-1) {
printf $gd "%d ", $i-$d_n/2;
for my $ac (0 .. $accn-1) {
$ar = $$x[$ac];
$s0 = $$ar{s0};
$s1 = $$ar{s1};
$s2 = $$ar{s2};
$s3 = $$ar{s3};
printf $gd "%g ", ($$s0[$i] - $$s1[$i])/2.0;
# printf $gd "%g ", $$s1[$i];
}
printf $gd "\n";
}
close $gd;
doplot($opt_fn . "-q.png", "Q stokes parameter, ${d_bits}-bits",
$ymin, $ymax, \@info);
}

if ($d_type==2) {
open $gd, "> plot.dat" or die "Cannot open plot.dat";
for my $i (0 .. $d_n-1) {
printf $gd "%d ", $i-$d_n/2;
for my $ac (0 .. $accn-1) {
$ar = $$x[$ac];
$s0 = $$ar{s0};
$s1 = $$ar{s1};
$s2 = $$ar{s2};
$s3 = $$ar{s3};
printf $gd "%g ", $$s2[$i];
}
printf $gd "\n";
}
close $gd;
doplot($opt_fn . "-u.png", "U stokes parameter, ${d_bits}-bits",
$ymin, $ymax, \@info);

open $gd, "> plot.dat" or die "Cannot open plot.dat";
for my $i (0 .. $d_n-1) {
printf $gd "%d ", $i-$d_n/2;
for my $ac (0 .. $accn-1) {
$ar = $$x[$ac];
$s0 = $$ar{s0};
$s1 = $$ar{s1};
$s2 = $$ar{s2};
$s3 = $$ar{s3};
printf $gd "%g ", $$s3[$i];
}
printf $gd "\n";
}
close $gd;
doplot($opt_fn . "-v.png", "V stokes parameter, ${d_bits}-bits",
$ymin, $ymax, \@info);
}

exit(0);