#!/usr/bin/perl -w
use strict;
use GD;
use GD::Text;
use GD::Graph;
# note this is Andrew's version...
use GD::Graph::histogram;
use GD::Graph::points;
use GD::Graph::lines;


my ($infile, $ASIC)=@ARGV;

local $| = 1;


open INPUT, "<$infile" or die "ERROR - Can't open $infile\n";
my @lines = <INPUT>;
close INPUT;

my $npix = $#lines + 1;
my (@x,@y,@m,@b,@m_norm,@b_norm,@r,@peaks,@fwhms,@res,@fset) = ((),(),(),(),(),(),(),(),(),(),());
my ($ii,$jj,$kk) = (0,0,0);

my @fields = ();

while( $ii< $npix ) {
   @fields = split " ", $lines[$ii];
   chomp $fields[$#fields];
   $x[$ii] = $fields[0];
   $y[$ii] = $fields[1];
   for ($jj=0;$jj<5;$jj++) {
      $peaks[$ii][$jj] = $fields[2 +$jj];
   }
   for ($jj=0;$jj<5;$jj++) {
      $fwhms[$ii][$jj] = $fields[7 +$jj];
      if ($fwhms[$ii][$jj] > 38.0) {$fwhms[$ii][$jj] = 38.0}
   }
   $m[$ii] = $fields[12];
   $b[$ii] = $fields[13];
   $r[$ii] = $fields[14];
   for ($jj=0;$jj<5;$jj++) {
      $res[$ii][$jj] = $fields[15 +$jj];
   }
   for ($jj=0;$jj<16;$jj++) {
      $fset[$ii][$jj] = $fields[20 +$jj];
   }
   $ii++;
}

my ($m_av, $b_av) = (0,0);
my ($m_sdev, $b_sdev) = (0,0);
for ($ii=0;$ii<$npix;$ii++) {
   $m_av += $m[$ii];
   $b_av += $b[$ii];
}
$m_av = $m_av/$npix;
$b_av = $b_av/$npix;
$b_av = $b_av/$m_av;

for ($ii=0;$ii<$npix;$ii++) {
   $m_norm[$ii] = $m[$ii]/$m_av;
#   $m_sdev += (($m_norm[$ii] - 1)*($m_norm[$ii] - 1));
   $m_sdev += (($m[$ii] - $m_av)*($m[$ii] - $m_av));
#   $b_norm[$ii] = $b[$ii]/$b_av;
   $b_norm[$ii] = $b[$ii]/$m_av;
   $b_sdev += (($b_norm[$ii] - $b_av)*($b_norm[$ii] - $b_av));
}
$m_sdev = sqrt($m_sdev/($npix-1));
$b_sdev = sqrt($b_sdev/($npix-1));

my $gainhist = new GD::Graph::histogram(400,400);
my $title = sprintf("%s Gain Dist: mean= %5.3f  sdev= %5.3f",$ASIC,$m_av,$m_sdev);

$gainhist->set( 
                x_label         => 'Gain (Channels/KeV)',
                y_label         => 'Counts/bin',
                title           => $title,
                x_labels_vertical => 0,
		x_tick_number	=> 'auto',
                bar_spacing     => 0,
                shadow_depth    => 1,
                shadowclr       => 'dred',
                transparent     => 0,
		histogram_type  => 'count',
		histogram_bins  => 60,
		uxmin	=> 65,
		uxmax	=> 80
) 
or warn $gainhist->error;

my $gd = $gainhist->plot(\@m) or die $gainhist->error;
#my $filenm = sprintf("%s/gainhist.png",$ASIC);
my $filenm = "png/gainhist.png";
open(IMG, ">$filenm") or die $!;
binmode IMG;
print IMG $gd->png;
close(IMG);

my @pixnum = ();
for ($ii=0;$ii<$npix;$ii++) {
   $pixnum[$ii] = $x[$ii]*32 + $y[$ii];
}
my @gvpdata = (\@pixnum,\@m_norm);

my $gainvpixel = new GD::Graph::points(1000,400);
$title = sprintf("%s Gain vs Pixel#",$ASIC);

$gainvpixel->set( 
                x_label         => 'Pixel #',
                y_label         => 'Normalized Gain',
                title           => $title,
                x_labels_vertical => 0,
		x_long_ticks	=> 1,
		x_tick_number	=> 32,
#		y_tick_number	=> 'auto',
		y_min_value	=> 0.85,
		y_max_value	=> 1.06,
		x_min_value	=> 0,
		x_max_value	=> 1024,
		marker_size	=> 1,
                transparent     => 0
) 
or warn $gainvpixel->error;

$gd = $gainvpixel->plot(\@gvpdata) or die $gainhist->error;
#$filenm = sprintf("%s/gainvpixel.png",$ASIC);
$filenm = "png/gainvpixel.png";
open(IMG, ">$filenm") or die $!;
binmode IMG;
print IMG $gd->png;
close(IMG);


my $bhist = new GD::Graph::histogram(400,400);
$title = sprintf("%s Offset Dist: mean= %5.3f KeV sdev= %5.3f",$ASIC,$b_av,$b_sdev);

$bhist->set( 
                x_label         => 'Offset (b/<m>) (KeV)',
                y_label         => 'Counts/bin',
                title           => $title,
                x_labels_vertical => 0,
		x_tick_number	=> 'auto',
                bar_spacing     => 0,
                shadow_depth    => 1,
                shadowclr       => 'dred',
                transparent     => 0,
		histogram_type  => 'count',
		histogram_bins  => 60,
		uxmin	=> -0.3,
		uxmax	=> 1.0
) 
or warn $bhist->error;

$gd = $bhist->plot(\@b_norm) or die $bhist->error;
#$filenm = sprintf("%s/bhist.png",$ASIC);
$filenm = "png/bhist.png";
open(IMG, ">$filenm") or die $!;
binmode IMG;
print IMG $gd->png;
close(IMG);


my @bvpdata = (\@pixnum,\@b_norm);

my $bvpixel = new GD::Graph::points(1000,400);
$title = sprintf("%s Offset vs Pixel#",$ASIC);

$bvpixel->set( 
                x_label         => 'Pixel #',
                y_label         => 'Offset (b/<m>) (KeV)',
                title           => $title,
                x_labels_vertical => 0,
		x_long_ticks	=> 1,
		x_tick_number	=> 32,
#		y_tick_number	=> 'auto',
		y_min_value	=> -0.2,
		y_max_value	=> 1.0,
		x_min_value	=> 0,
		x_max_value	=> 1024,
		marker_size	=> 1,
                transparent     => 0
) 
or warn $bvpixel->error;

$gd = $bvpixel->plot(\@bvpdata) or die $bhist->error;
#$filenm = sprintf("%s/bvpixel.png",$ASIC);
$filenm = "png/bvpixel.png";
open(IMG, ">$filenm") or die $!;
binmode IMG;
print IMG $gd->png;
close(IMG);

for ($jj=0;$jj<5;$jj++) {
   my @fwhm = ();
   my ($wmean,$wmin,$wmax) = (0,1000,-1000);
   for ($ii=0;$ii<$npix;$ii++) {
      $fwhm[$ii] = $fwhms[$ii][$jj];
      $wmean += $fwhm[$ii];
      if ($fwhm[$ii] <$wmin) { $wmin = $fwhm[$ii]; }
      if ($fwhm[$ii] >$wmax) { $wmax = $fwhm[$ii]; }
   }
   $wmean = $wmean/$npix;

   my @wvpdata = (\@pixnum,\@fwhm);

   my $wvpixel = new GD::Graph::points(1000,400);
   $title = sprintf("%s Peak %d FWHM vs Pixel#. min=%4.1f max=%4.1f mean=%4.1f",$ASIC,$jj, $wmin,$wmax,$wmean);

   $wvpixel->set( 
                x_label         => 'Pixel #',
                y_label         => 'FWHM (Channels)',
                title           => $title,
                x_labels_vertical => 0,
		x_long_ticks	=> 1,
		x_tick_number	=> 32,
#		y_tick_number	=> 'auto',
		y_min_value	=> 0,
		y_max_value	=> 40,
		x_min_value	=> 0,
		x_max_value	=> 1024,
		marker_size	=> 1,
                transparent     => 0
   ) 
   or warn $wvpixel->error;
   $wvpixel->set_y_axis_font(gdSmallFont) or warn $wvpixel->error;

   $gd = $wvpixel->plot(\@wvpdata) or die $bhist->error;
   #$filenm = sprintf("%s/wvpixel%d.png",$ASIC,$jj);
   $filenm = sprintf("png/wvpixel%d.png",$jj);
   open(IMG, ">$filenm") or die $!;
   binmode IMG;
   print IMG $gd->png;
   close(IMG);


   my $whist = new GD::Graph::histogram(400,400);
   $title = sprintf("%s Peak %d FWHM: min=%4.1f max=%4.1f mean=%4.1f",$ASIC,$jj,$wmin,$wmax,$wmean);

   $whist->set( 
                x_label         => 'FWHM (Channels)',
                y_label         => 'Counts/bin',
                title           => $title,
                x_labels_vertical => 0,
		x_tick_number	=> 'auto',
                bar_spacing     => 0,
                shadow_depth    => 1,
                shadowclr       => 'dred',
                transparent     => 0,
		histogram_type  => 'count',
		histogram_bins  => 60,
		uxmin	=> 0,
		uxmax	=> 40
   ) 
   or warn $whist->error;
#   $whist->set_x_axis_font('arial', 18) or warn $whist->error;
   $whist->set_x_axis_font(gdSmallFont) or warn $whist->error;
   $whist->set_y_axis_font(gdSmallFont) or warn $whist->error;

   $gd = $whist->plot(\@fwhm) or die $whist->error;
   #$filenm = sprintf("%s/whist%d.png",$ASIC,$jj);
   $filenm = sprintf("png/whist%d.png",$jj);
   open(IMG, ">$filenm") or die $!;
   binmode IMG;
   print IMG $gd->png;
   close(IMG);
}





my @fsetcap = ();
my $av_fset = 0;
for ($jj=0;$jj<16;$jj++) {
   @fsetcap = ();
   $av_fset = 0;
   for ($ii=0;$ii<$npix;$ii++) {
      $fsetcap[$ii] = $fset[$ii][$jj];
      $av_fset += $fset[$ii][$jj];
   }
   $av_fset = $av_fset/$npix;
   my $capfsethist = new GD::Graph::histogram(200,200);
   $title = sprintf("%s Cap %2d Offsets",$ASIC,$jj);

   $capfsethist->set( 
                x_label         => 'Offset (Channels)',
                y_label         => 'Counts/bin',
                title           => $title,
                x_labels_vertical => 0,
		x_tick_number	=> 'auto',
                bar_spacing     => 0,
                shadow_depth    => 1,
                shadowclr       => 'dred',
                transparent     => 0,
		histogram_type  => 'count',
		histogram_bins  => 60,
		uxmin		=> -30,
		uxmax		=> 30
#		uxmin	=> int($av_fset) -40,
#		uxmax	=> int($av_fset) +40
   ) 
   or warn $capfsethist->error;

   $gd = $capfsethist->plot(\@fsetcap) or die $capfsethist->error;
   #$filenm = sprintf("%s/cap%02dfsethist.png",$ASIC,$jj);
   $filenm = sprintf("png/cap%02dfsethist.png",$jj);
#   print "$filenm\n";
   open(IMG, ">$filenm") or die $!;
   binmode IMG;
   print IMG $gd->png;
   close(IMG);

   my @capfsetvpdata = (\@pixnum,\@fsetcap);
   my $capfsetvpixel = new GD::Graph::points(1000,400);
   $title = sprintf("%s Cap %d Offset vs Pixel#",$ASIC,$jj);

   $capfsetvpixel->set( 
                x_label         => 'Pixel #',
                y_label         => 'Cap Offset',
                title           => $title,
                x_labels_vertical => 0,
		x_long_ticks	=> 1,
		x_tick_number	=> 32,
#		y_tick_number	=> 'auto',
		y_min_value	=> -40,
		y_max_value	=> 40,
		x_min_value	=> 0,
		x_max_value	=> 1024,
		marker_size	=> 1,
                transparent     => 0
   ) 
   or warn $capfsetvpixel->error;
   $capfsetvpixel->set_y_axis_font(gdSmallFont) or warn $capfsetvpixel->error;

   $gd = $capfsetvpixel->plot(\@capfsetvpdata) or die $bhist->error;
   $filenm = sprintf("png/cap%02dfsetvpixel.png",$jj);
   open(IMG, ">$filenm") or die $!;
   binmode IMG;
   print IMG $gd->png;
   close(IMG);
}


my @linx = (1,2,3,4,5);
my @resdata = ();
my @resmax = ();
my $res_xmax = 0.0006;
my $nout = 0;
for($ii=0;$ii<$npix;$ii++) {
   @resdata = ();
   my $rmax = -1000;
   my $rmin = 1000;
   for($jj=0;$jj<5;$jj++) {
      $resdata[$jj] = $res[$ii][$jj];
      if ($resdata[$jj] > $rmax) { $rmax = $resdata[$jj];}
      if ($resdata[$jj] < $rmin) { $rmin = $resdata[$jj];}
#      printf(" %8.6f  ",$resdata[$jj]);
   }
   if (abs($rmin) > $rmax) {
      $resmax[$ii] = abs($rmin);
   } else {
      $resmax[$ii] = $rmax;
   }
   if ($resmax[$ii] > $res_xmax) {
      $nout ++;
   }
#   printf(" %8.6f\n",$resmax[$ii]);

   my @lindata = (\@linx,\@resdata);

   my $resplot = new GD::Graph::lines(200,200);
   $title = sprintf("%s Pixel (%2d,%2d)",$ASIC,$x[$ii],$y[$ii]);

   $resplot->set( 
#                x_label         => undef,
                y_label         => 'Residual',
                title           => $title,
                x_labels_vertical => 0,
#		x_tick_number	=> 'auto',
		y_tick_number	=> 4,
		y_min_value	=> -0.0004,
		y_max_value	=> 0.0004,
		x_min_value	=> 0,
		x_max_value	=> 6,
		marker_size	=> 1,
                transparent     => 0
   ) 
   or warn $resplot->error;
   $resplot->set_y_axis_font(gdSmallFont) or warn $resplot->error;
   
   $gd = $resplot->plot(\@lindata) or die $bhist->error;
   $filenm = sprintf("png/%02d%02dresplot.png",$x[$ii],$y[$ii]);
   open(IMG, ">$filenm") or die $!;
   binmode IMG;
   print IMG $gd->png;
   close(IMG);
}


my $mreshist = new GD::Graph::histogram(400,400);
$title = sprintf("%s Max Linearity Deviation",$ASIC);

$mreshist->set( 
                x_label         => 'Max Deviation from Linearity',
                y_label         => 'Counts/bin',
                title           => $title,
                x_labels_vertical => 1,
		x_tick_number	=> 6,
                bar_spacing     => 0,
                shadow_depth    => 1,
                shadowclr       => 'dred',
                transparent     => 0,
		histogram_type  => 'count',
		histogram_bins  => 60,
		uxmin	=> 0.0,
		uxmax	=> $res_xmax
) 
or warn $mreshist->error;
my $legend = "Outliers: $nout";
my @lkeys = ($legend);
$mreshist->set_legend(@lkeys) or warn $mreshist->error;

my @resmax_r = ();
for ($ii=0;$ii<=$#resmax;$ii++) {
   if ($resmax[$ii] < $res_xmax) {
      $resmax_r[++$#resmax_r] = $resmax[$ii];
   }
}

$gd = $mreshist->plot(\@resmax_r) or die $mreshist->error;
# $filenm = sprintf("%s/mreshist.png",$ASIC);
$filenm = "png/mreshist.png";
open(IMG, ">$filenm") or die $!;
binmode IMG;
print IMG $gd->png;
close(IMG);




my @ctable = ();
my ($r,$g,$b) = (0,0,0);

# Gain map
   my $im = new GD::Image(288,344) || die;
   $ctable[0] = $im->colorAllocate(255,255,255);
   $ctable[1] = $im->colorAllocate(0,0,0);
   for ($ii=2;$ii<256;$ii++) {
      ($r,$g,$b) = colorRamp($ii,0,255);
#      print "$ii $r  $g  $b\n";
     $ctable[$ii] = $im->colorAllocate($r,$g,$b);
   }


$im->filledRectangle(16,16,271,271,$ctable[1]);   
for ($ii=0;$ii<$npix;$ii++) {
      my $clr = int( (($m[$ii] - 75) * 256/10) + 128 );
      if ($clr < 2) {$clr = 2;}
      if ($clr >255) {$clr = 255;}
      my $x1 = $x[$ii] * 8 + 16;
      my $y1 = $y[$ii] * 8 + 16;
      $im->filledRectangle($x1,$y1,$x1+8,$y1+8,$ctable[$clr]);   
}

   # intensity bar
for ($ii=0;$ii<256;$ii++) {
      my $x1 = $ii +16;
      my $y1 = 288;
      $im->filledRectangle($x1,$y1,$x1,$y1+7,$ctable[$ii]);
}
$title = sprintf("%s Gain (channels/keV",$ASIC);
$im->string(gdLargeFont,16,2,$title,$ctable[1]);
my $ibarl0 = "70";
my $ibarl1 = "80";
$im->string(gdMediumBoldFont,16,304,$ibarl0,$ctable[1]);
$im->string(gdMediumBoldFont,264,304,$ibarl1,$ctable[1]);

$filenm = sprintf("png/gainmap.png");
open(IMG, ">$filenm") or die $!;
binmode IMG;
print IMG $im->png;
close(IMG);


# Offset map
   $im = new GD::Image(288,344) || die;
   $ctable[0] = $im->colorAllocate(255,255,255);
   $ctable[1] = $im->colorAllocate(0,0,0);
   for ($ii=2;$ii<256;$ii++) {
      ($r,$g,$b) = colorRamp($ii,0,255);
#      print "$ii $r  $g  $b\n";
     $ctable[$ii] = $im->colorAllocate($r,$g,$b);
   }


$im->filledRectangle(16,16,271,271,$ctable[1]);   
for ($ii=0;$ii<$npix;$ii++) {
      my $clr = int( (($b_norm[$ii] - 0.25) * 256/1.0) + 128 );
      if ($clr < 2) {$clr = 2;}
      if ($clr >255) {$clr = 255;}
      my $x1 = $x[$ii] * 8 + 16;
      my $y1 = $y[$ii] * 8 + 16;
      $im->filledRectangle($x1,$y1,$x1+8,$y1+8,$ctable[$clr]);   
}

   # intensity bar
for ($ii=0;$ii<256;$ii++) {
      my $x1 = $ii +16;
      my $y1 = 288;
      $im->filledRectangle($x1,$y1,$x1,$y1+7,$ctable[$ii]);
}
$title = sprintf("%s Offset (b/m) (keV)",$ASIC);
$im->string(gdLargeFont,16,2,$title,$ctable[1]);
$ibarl0 = "-0.25";
$ibarl1 = "0.75";
$im->string(gdMediumBoldFont,16,304,$ibarl0,$ctable[1]);
$im->string(gdMediumBoldFont,264,304,$ibarl1,$ctable[1]);

$filenm = sprintf("png/fsetmap.png");
open(IMG, ">$filenm") or die $!;
binmode IMG;
print IMG $im->png;
close(IMG);


sub colorRamp {
    my( $v, $vmin, $vmax ) = @_;
    my( $r, $g, $b ) = (1) x 3;

    our @ctable;
    $v = $vmax if $v > $vmax;
    $v = $vmin if $v < $vmin;

## Uncomment below to invert the color range 
## so that small numbers are hot (red) 
## and high numbers are cold (blue)
##    $v = $vmax + $vmin - $v; 

    my $dv = $vmax - $vmin;
    if( $v < ( $vmin + 0.25*$dv ) ) {
        $r = 0;
        $g = 4 * ($v - $vmin) / $dv;
    }
    elsif( $v < ( $vmin + 0.5 * $dv ) ) {
        $r = 0;
        $b = 1 + 4 * ($vmin + 0.25 * $dv - $v) / $dv;
    }
    elsif( $v < ( $vmin + 0.75 * $dv ) ) {
        $r = 4 * ($v - $vmin - 0.5 * $dv) / $dv;
        $b = 0;
    }
    else {
        $g = 1 + 4 * ($vmin + 0.75 * $dv - $v) / $dv;
        $b = 0;
    }
    ## Convert the 0->1 range rgb values to 0->255 
    ## and output as a web/tk style "#rrggbb" hex string
    ## return sprintf "#%02x%02x%02x", $r*255, $g*255, $b*255;

    return int($r*255),int($g*255),int($b*255);
}