#!perl -w

=pod

=head1 NAME

Neuron::Morphometrics::Compartment - Information on compartments

=head1 SYNOPSIS

Neuron::Morphometrics::Compartment provides different functions to return
morphometric values concerning a neuron's compartments, like average diameter,
area, volume etc.

=head1 FUNCTIONS

In all the following functions, the optional parameters I<min> and I<max> define
limits in the parameter specified by I<which>: branch order ("order"), path
distance ("distance") or euclidian distance ("edistance"); for example,
numSegments(\@myneuron, 2, 4, order) would only consider segments of branch
orders 2 through 4, inclusively. If parameters are given at all, both I<min>
and I<max> must be given.

=over 0

=item *
I<$n> = avgArea(I<neuron> [, I<min>, I<max>, I<which>]);

returns average compartment area.

=item *
I<$n> = avgVolume(I<neuron> [, I<min>, I<max>, I<which>]);

returns average compartment volume.

=item *
I<$n> = (min|max|avg)BurkeTaperRate(I<neuron> [, I<min>, I<max>, I<which>]);

returns minimal, maximal or average Burke's Taper Rate of compartments.

=item *
I<$n> = (min|max|avg)CompartmentEuclidianDistance(I<neuron> [, I<min>, I<max>, 
I<which>]);

returns minimal, maximal or average euclidian distance of compartments from
soma.

=item *
I<$n> = (min|max|avg)CompartmentLength(I<neuron> [, I<min>, I<max>, I<which>]);

returns minimal, maximal or average length (path length equals euclidian length)
of compartments.

=item *
I<$n> = (min|max|avg)CompartmentPathDistance(I<neuron> [, I<min>, I<max>,
I<which>]);

returns minimal, maximal or average path distance of compartments.

=item *
I<$n> = (min|max|avg)Diameter(I<neuron> [, I<min>, I<max>, I<which>]);

returns minimal, maximal or average diameter of compartments.

=back

=head1 AUTHOR

Gunther Schmidl <gschmidl@gmx.at>

=cut

package Neuron::Morphometrics::Compartment;
use strict;
use warnings;
no warnings 'recursion';

use Neuron::Morphometrics::EntireNeuron;
use Neuron::Morphometrics::Tools;

# prepare variable export

use vars qw(@ISA @EXPORT @EXPORT_OK $VERSION);
require Exporter;
@ISA = qw(Exporter);
@EXPORT = qw(minDiameter maxDiameter avgDiameter avgVolume avgArea
    minCompartmentPathDistance avgCompartmentPathDistance
    maxCompartmentPathDistance minCompartmentLength avgCompartmentLength
    maxCompartmentLength minCompartmentEuclidianDistance
    avgCompartmentEuclidianDistance maxCompartmentEuclidianDistance
    minBurkeTaperRate avgBurkeTaperRate maxBurkeTaperRate);
@EXPORT_OK = qw();
$VERSION = 0.01;






#==============================================================================
#
# generic min/max subroutine for branch order limitation
#
# usage: lminmax(array, variable, value, lolimit, uplimit, which)
#
# variable    the hash key to look up
# value       0 for minimum, 1 for maximum, 2 for average
# lolimit     lower limit
# hilimit     higher limit
# which       order, distance, euclidian distance
#
#==============================================================================

sub lminmax
{
    my @arr = @{$_[0]};
    my $min = (defined($_[3]) ? $_[3] : 0);
    my $max = (defined($_[4]) ? $_[4] : maxBranchOrder(\@arr));
    my $wch = (defined($_[5]) ? $_[5] : "order");

    my $num = numCompartments(\@arr);
    my $count = 0;                          # for average
    my $comp;                               # workhorse value

    SWITCH:
    {
        $comp = 65000, last SWITCH if ($_[2] == 0);
        $comp = -1,    last SWITCH if ($_[2] == 1);
        $comp = 0,     last SWITCH if ($_[2] == 2);
    }

    for (my $i = 1; $i <= $num; $i++)
    {
        # for compartment statistics, we include the stem if order is 0
        if ((($arr[$i]->{$wch} >= $min) && ($arr[$i]->{$wch} <= $max))
            || ($min == 0 && $i == 1))
        {
            if ($_[2] == 0)                 # minimum
            {
                if ($arr[$i]->{"$_[1]"} < $comp)
                {
                    $comp = $arr[$i]->{"$_[1]"};
                }
            }
            elsif ($_[2] == 1)              # maximum
            {
                if ($arr[$i]->{"$_[1]"} > $comp)
                {
                    $comp = $arr[$i]->{"$_[1]"};
                }
            }
            else                            # average
            {
                $count++;
                $comp = $comp + $arr[$i]->{"$_[1]"};
            }
        }
    }

    if ($_[2] == 2)
    {
        $comp = (($count == 0) ? 0 : ($comp / $count));
    }

    return $comp;
}

#==============================================================================
#
# (min|max|avg)Diameter(neuron [, min, max, which]) - return min/max/average 
# diameter.
#
#==============================================================================

sub minDiameter
{
    return lminmax(\@{$_[0]}, "d", 0, $_[1], $_[2], $_[3]);
}

sub maxDiameter
{
    return lminmax(\@{$_[0]}, "d", 1, $_[1], $_[2], $_[3]);
}

sub avgDiameter
{
    return lminmax(\@{$_[0]}, "d", 2, $_[1], $_[2], $_[3]);
}






#==============================================================================
#
# avgArea(neuron [, min, max, which]) - return average area.
#
#==============================================================================

sub avgArea
{
    return lminmax(\@{$_[0]}, "area", 2, $_[1], $_[2], $_[3]);
}






#==============================================================================
#
# avgVolume(neuron [, min, max, which]) - return average volume.
#
#==============================================================================

sub avgVolume
{
    return lminmax(\@{$_[0]}, "volume", 2, $_[1], $_[2], $_[3]);
}






#==============================================================================
#
# (min|max|avg)CompartmentPathDistance(neuron [, min, max, which]) - return 
# min/max/average path distance of compartment from soma.
#
#==============================================================================

sub minCompartmentPathDistance
{
    return lminmax(\@{$_[0]}, "distance", 0, $_[1], $_[2], $_[3]);
}

sub maxCompartmentPathDistance
{
    return lminmax(\@{$_[0]}, "distance", 1, $_[1], $_[2], $_[3]);
}

sub avgCompartmentPathDistance
{
    return lminmax(\@{$_[0]}, "distance", 2, $_[1], $_[2], $_[3]);
}






#==============================================================================
#
# (min|max|avg)CompartmentLength(neuron [, min, max, which]) - return 
# min/max/average compartment length. (euclidian length == path length for a 
# compartment)
#
#==============================================================================

sub minCompartmentLength
{
    return lminmax(\@{$_[0]}, "length", 0, $_[1], $_[2], $_[3]);
}

sub maxCompartmentLength
{
    return lminmax(\@{$_[0]}, "length", 1, $_[1], $_[2], $_[3]);
}

sub avgCompartmentLength
{
    return lminmax(\@{$_[0]}, "length", 2, $_[1], $_[2], $_[3]);
}






#==============================================================================
#
# (min|max|avg)CompartmentEuclidianDistance(neuron [, min, max, which]) - 
# return min/max/average euclidian distance of compartment from soma.
#
#==============================================================================

sub minCompartmentEuclidianDistance
{
    return lminmax(\@{$_[0]}, "edistance", 0, $_[1], $_[2], $_[3]);
}

sub maxCompartmentEuclidianDistance
{
    return lminmax(\@{$_[0]}, "edistance", 1, $_[1], $_[2], $_[3]);
}

sub avgCompartmentEuclidianDistance
{
    return lminmax(\@{$_[0]}, "edistance", 2, $_[1], $_[2], $_[3]);
}






#==============================================================================
#
# (min|max|avg)BurkeTaperRate(neuron [, min, max, which]) - return 
# min/max/average of Burke's Taper Rate, which is defined as 
# ((parent diameter - daughter diameter) / (parent length + daughter length)).
#
#==============================================================================

sub minBurkeTaperRate
{
    return lminmax(\@{$_[0]}, "taper", 0, $_[1], $_[2], $_[3]);
}

sub maxBurkeTaperRate
{
    return lminmax(\@{$_[0]}, "taper", 1, $_[1], $_[2], $_[3]);
}

sub avgBurkeTaperRate
{
    return lminmax(\@{$_[0]}, "taper", 2, $_[1], $_[2], $_[3]);
}
