view release on metacpan or  search on metacpan

bin/MorphometryI  view on Meta::CPAN

{
   my $w = shift;

   # input: xy array
   # output in array context: (area,perimeter,formfactor,xcen,ycen,fxcen,fycen)
   # initialization:
   my ($xCalib, $yCalib) = ( $w-> {ini}-> {XCalibration}, $w-> {ini}-> {YCalibration});

   # algorithm
   my $xflag = 1;
   my( @x, @y);
   for (@_) {
      push @x, $_ if $xflag;
      push @y, $_ unless $xflag;
      $xflag = !$xflag;
   }
   return () unless @x;
   unless ($x[$#x] == $x[0] && $y[$#y] == $y[0]) {
      push @x, $x[0];
      push @y, $y[0];
   }
   my ($xyCalib,$xxCalib,$yyCalib) = ($xCalib*$yCalib,$xCalib*$xCalib,$yCalib*$yCalib);
   my ($area,$perimeter,$xcen,$ycen,$fxcen,$fycen,$ff) = (0,0,0,0,0,0,0,0,0);
   for my $i ( 1..$#x) {
      $area += $xyCalib * ($x[$i-1] * $y[$i] - $x[$i] * $y[$i-1]);
      my $dx = $x[$i-1] - $x[$i];
      my $dy = $y[$i-1] - $y[$i];
      $perimeter += sqrt( $xxCalib * $dx * $dx + $yyCalib * $dy *$dy);
      $xcen += $x[$i];
      $ycen += $y[$i];
      $fxcen += $xCalib * $x[$i];
      $fycen += $yCalib * $y[$i];
   }
   $area = abs( $area / 2);
   $ff = 4 * PI * $area / $perimeter / $perimeter
      if $perimeter > 0;
   $xcen /= @x;
   $ycen /= @y;
   $fxcen /= @x;
   $fycen /= @y;
   return ($area,$perimeter,$ff,$xcen,$ycen,$fxcen,$fycen);
}


sub win_saveframe
{
   my $w = $_[0];
   my $xmlname = $w-> win_extname( $w-> {file});

   return 1 unless $w-> {modified};

   if ( open F, "> $xmlname") {
      my $waitPtr = $::application-> pointer;
      $::application-> pointer( cr::Wait);
      $w-> sb_text("saving $xmlname");
      $w-> win_validate(1);
      $w-> win_syncrecdata;
      my $image = $w-> IV-> image;
      if ( $w->{ini}->{CalcBrightness} && $w->{ini}->{EqualBrightness}) {
         # subtracting low frequencies
         $w-> sb_text("Equalizing background ...");
         my $i1 = Prima::IPA::Global::butterworth( $image, 
            low        => 1,
            homomorph  => 0,
            power      => 2,
            cutoff     => 20,
            boost      => 0.7,
            spatial    => 1,
            lowquality => 1,
         );
         $i1-> type( $image-> type);
         $image = Prima::IPA::Point::subtract( $image, $i1);
         $w-> sb_text("saving $xmlname");
      }   
         

      my ( $iname, $ix, $iy, $path, $datestr, $xc, $yc, $objects, $i) = (
         $w->{file}, $image-> size, $w->{ini}->{path}, scalar(gmtime(time)),
         $w->{ini}->{XCalibration}, $w->{ini}->{YCalibration}, 0
      );
      $iname =~ m{[/\\]([^/\\]*)$};
      $iname = $1;

      for ( $i = 0; $i < 2; $i++) {
         my $ptr = $w-> {lineStorage}->[$i];
         next unless defined $ptr;
         $objects += scalar @$ptr;
      }
      $objects += scalar @{$w-> {points}} / 2 if defined $w-> {points};
      my $objCount = 0;

print F <<HEADER;
<?xml version="1.0"?>
<!DOCTYPE morphology_data SYSTEM "morphology_data.dtd">
<!-- This is a generated file.  Do not edit! -->
<morphology_data
  imagename    = "$iname"
  imagewidth   = "$ix"
  imageheight  = "$iy"
  directory    = "$path"
  creator      = "MorphologyI"
  creationdate = "$datestr"
  xcalib       = "$xc"
  ycalib       = "$yc"
  objects      = "$objects"
>

HEADER

      for ( $i = 0; $i < 2; $i++) {
         my $ptr = $w-> {lineStorage}->[$i];
         next unless defined $ptr;
         my $type = $w-> menu-> text( $i);
         $type =~ s[\~][]g;
         $type = lc $type;
         for ( @$ptr) {
             $objCount++;
             my ( $xs, $ys) = ( '', '', '', '');
             my $j;
             my $pp = $_;
             for ( $j = 0; $j < scalar @$pp; $j += 2) {

bin/MorphometryI  view on Meta::CPAN

      $d-> insert( CheckList => 
         origin => [ 10, 76],
         size   => [ 300, 110],
         name   => 'CheckList',
         vector => $w-> {ini}->{CalcOptions},
         vScroll=> 1,
         items  => [
            'List of files',
            'Morphometric parameters',
            'Brightness histogram',
            'Subtract backround from histogram',
         ],
      );
      $w-> dlg_okcancel( $d);
      $w->{CalcOptDialog} = $d;
   }
   return if $d-> execute != mb::OK;
   return $w-> {ini}-> {CalcOptions} = $d-> CheckList-> vector;
}

sub opt_statistics
{
   my $w = $_[0];
   if ( defined $w-> {file}) {
      return unless $w-> win_saveframe;
   }
   my $d = $w-> dlg_file(
      cwd       => 1,
      directory => $w-> {ini}-> {path},
      filter    => [
         ['Data files' => '*.xml'],
         ['All files' => '*']
      ],
      multiSelect => 1,
      filterIndex => 0,
   );
   my @res = $d-> execute;
   return unless @res;

   my $dres = $w-> opt_calcdlg;
   return unless defined $dres;

   my %dataCommon = ();
   my %dataFore   = ();
   my %dataBack   = ();
   my @dataHist;
   my %disabled = map {$_=>1} qw(x y type brightnessN brightnessSum brightnessSum2
      xcentroid ycentroid fxcentroid fycentroid harea histogram
      convex_xcentroid convex_ycentroid convex_fxcentroid convex_fycentroid
      );
   my $index = 0;
   my $statCount = 4;
   my @resName;


   for my $xmlname ( @res) {
      my ( $bN, $bSum, @hist);
      push( @resName, "$xmlname : ");
      my $state = $w-> win_xmlload( $xmlname, onObject => sub {
         my ( $attrs, $n0, $n1) = @_;
         # collect sum of backgrounds
         if ( $attrs->{type} eq $n0 || $attrs->{type} eq $n1) {
            my $stref = $attrs->{type} eq $n0 ? \%dataFore : \%dataBack;
            if ( 3 == grep {exists $attrs->{$_}} qw{brightnessN brightnessSum brightnessSum2}) {
               $stref->{brightness} = [(0) x $statCount] unless defined $stref->{brightness};
               $stref->{brightness}-> [0] += $attrs->{brightnessN};
               $stref->{brightness}-> [1] += $attrs->{brightnessSum};
               $stref->{brightness}-> [2] += $attrs->{brightnessSum2};
               # XXX warning - bogus data for backgrounds !!
               $stref->{brightness}-> [$statCount + $index] = $attrs->{brightnessSum} / $attrs->{brightnessN};
               # collect local sum of background data
               if ( $attrs->{type} eq $n1) {
                  $bN   += $attrs->{brightnessN};
                  $bSum += $attrs->{brightnessSum};
               }
            }
            return if $attrs->{type} eq $n1; # no background relevant data further

            push @dataHist, exists( $attrs->{histogram}) ? [split(' ', $attrs->{histogram})] : undef;
            push( @hist, $dataHist[-1]) if $dataHist[-1]; # push the histogram ref for the post-processing
            
            if ( exists $attrs->{area} && exists $attrs->{harea}) {
               my $nh = 0;
               my $tharea = 0;
               my $minArea = $attrs->{area} * $w->{ini}->{HolesPercent} / 100;
               my @holz = split ' ', $attrs->{harea};
               for my $harea ( @holz) {
                  next if $harea < $minArea;
                  $nh++;
                  $tharea += $harea;
               }
               $dataCommon{process_index} = [(0) x $statCount]
                  unless defined $dataCommon{process_index};
               $dataCommon{process_index}-> [0]++;
               $dataCommon{process_index}-> [1] += $nh;
               $dataCommon{process_index}-> [2] += $nh * $nh;
               $dataCommon{process_index}-> [$statCount + $index] = $nh;
               $dataCommon{process_domain} = [(0) x $statCount]
                  unless defined $dataCommon{process_domain};
               $dataCommon{process_domain}-> [0]++;
               $dataCommon{process_domain}-> [1] += $tharea;
               $dataCommon{process_domain}-> [2] += $tharea * $tharea;
               $dataCommon{process_domain}-> [$statCount + $index] = $tharea;
            }
            for ( keys %{$attrs}) {
               next if $disabled{$_};
               $dataCommon{$_} = [(0) x $statCount] unless defined $dataCommon{$_};
               $dataCommon{$_}-> [0]++;
               $dataCommon{$_}-> [1] += $attrs->{$_};
               $dataCommon{$_}-> [2] += $attrs->{$_} * $attrs->{$_};
               $dataCommon{$_}-> [$statCount + $index] = $attrs->{$_};
            }
            $index++;
            $resName[-1] .= "$index ";
         }
      } 
      );
      return unless $state;
      # subtract backgrounds
      if ( vec( $dres, ocq::Equalize, 1)) {
          if (!$bN && @hist) {
             # noting to subtract, nulling histograms
             splice( @dataHist, -(@hist), (@hist), (undef)x scalar(@hist));
          } elsif ( @hist) {
             my $b = int($bSum / $bN); # average brightness for background(s) within file
             if ( $b > 255 || $b < 0) {
                $@ = "Somehow the average brightness is not in range of 0..255.\n";
                $w-> win_xmlerror( $xmlname);
             }   
             for ( @hist) { # shift array left $b times
                splice( @$_, 0, 0, (0)x(255-$b)); # add space for negative brightnesses
                push( @$_, (0)x($b)); 
             }
          }   
      }
   }

   $dataCommon{brightness} = $dataFore{brightness} if defined $dataFore{brightness};
   $dataCommon{backingbrightness} = $dataBack{brightness} if defined $dataBack{brightness};
   for ( keys %dataCommon) {
      $dataCommon{$_}->[1] /= $dataCommon{$_}->[0];
      $dataCommon{$_}->[2] = undef, next if $dataCommon{$_}->[0] < 2;
      $dataCommon{$_}->[2] = ($dataCommon{$_}->[2] - $dataCommon{$_}->[0] *
         $dataCommon{$_}->[1] * $dataCommon{$_}->[1]) /
         ($dataCommon{$_}->[0] - 1);
      $dataCommon{$_}->[3] = $dataCommon{$_}->[2] / $dataCommon{$_}->[0];
      $dataCommon{$_}->[2] = ( $dataCommon{$_}->[2] > 0) ? sqrt( $dataCommon{$_}->[2]) : 0;
      $dataCommon{$_}->[3] = ( $dataCommon{$_}->[3] > 0) ? sqrt( $dataCommon{$_}->[3]) : 0;
   }
   delete $dataCommon{backingbrightness};

   my @fmt = ();
   my @order = qw(area perimeter formfactor length width breadth
                  brightness convex_area
                  convex_perimeter convex_formfactor convex_width
                  process_index process_domain length_width
                  convex_length_width spreading_index);
   my @names = ();
   my %order = ();
   for (@order) {
      if ( defined $dataCommon{$_}) {
         push @names, $_;
         $order{$_} = 1;
      }
   }
   for ( sort keys %dataCommon) {
      push @names, $_ unless exists $order{$_};
   }
   for ( @names) {
      push @fmt, 13;
      $fmt[-1] = length( $_) if $fmt[-1] < length( $_);
   }
   my $colsformat2 = " %3s ";
   my $colsformat  = " $colsformat2| ";
   my $exformat   = '';
   $exformat .= '%'.$_.'s | ' for @fmt;
   $exformat .= "\n";
   my $exformat2 = " $exformat";
   $exformat2 =~ s/\s\|//g;

   my $divline    = 7 + 3 * scalar @names;
   $divline += $_ for @fmt;
   $divline = '-' x $divline;
   $divline .= "\n";

   my $hdr        = sprintf( "    # | $exformat", @names);
   my $texts = '';
   my $tm = localtime;
   $texts = <<HR;
Statistics morphology data
Calculated on $tm
HR
   if ( vec( $dres, ocq::Files, 1)) {
      my $howmany = ( scalar @resName > 1) ? "s were" : " was";
      $texts .= "The following file$howmany used for data gathering:\n";
      $texts .= "   $_\n" for @resName;
   }
   $texts .= sprintf("Average background brightness: %.7g +/- %.7g (averaged by %d pixels)\n",
      $dataBack{brightness}->[1],
      $dataBack{brightness}->[3],
      $dataBack{brightness}->[0],
   ) if defined $dataBack{brightness};
   my ( $title, $meta1, $meta2) = ( "$texts\n", '', '');

   if ( vec( $dres, ocq::Basics, 1)) 
   {
      {
         (my $xhdr = $hdr) =~ s/\s\|//g;
         my $xdivline = '-' x length $xhdr;
         $meta1 = $meta2 = "$xdivline\n$xhdr$xdivline\n";
      }

      $texts .="\n$divline$hdr$divline";
      my $i;
      my @avail = ( 1 ) x ( scalar @names);

      for ( $i = 0; $i < $index; $i++) {
         my @data = map { $dataCommon{$_}->[$i + $statCount]} @names;
         my $j = 0;
         for ( @data) {
            $avail[ $j++] &= defined ($_) ? 1 : 0;
            $_ = defined $_ ? sprintf( '%13.7g', $_) : "N/A";
         }
         $texts .= sprintf( "$colsformat$exformat",   $i + 1, @data);
         $meta1 .= sprintf( "$colsformat2$exformat2", $i + 1, @data);
      }
      $texts .= $divline . $hdr . $divline;
      my @elabels = qw( N Ave SD SEM);
      for $i ( 1, 3, 2, 0) {
         my @data = map { $dataCommon{$_}->[$i]} @names;
         my $j = 0;
         for ( @data) {
            $_ = defined $_ ? ( sprintf( '%1'.($avail[$j] ? '3' : '2').'.7g', $_).($avail[$j] ? '' : '*')) : "N/A";
            $j++;
         }
         $texts .= sprintf( "$colsformat$exformat",  $elabels[$i], @data);
         $meta2 .= sprintf( "$colsformat2$exformat2", $elabels[$i], @data);
      }
      $texts .= $divline;
   }
   
   if ( vec( $dres, ocq::Histogram, 1)) {
       my $maxSpace = 10;
       my ( $cnt, $min, $max) = ( 0, 255, 0);
       for ( @dataHist) {
          next unless defined $_;
          $cnt++;
          my $j = 0;
          $j++ while $j < @$_ && $$_[$j] == 0;
          $min = $j if $j < $min;
          $j = @$_ - 1;
          $j-- while $j && $$_[$j] == 0;
          $max = $j if $j > $max;
       }   
       for ( @dataHist) {
          splice( @$_, $max);
          splice( @$_, 0, $min);
       }