#include <CylBase.h>

Inheritance diagram for CylBase:

Collaboration diagram for CylBase:

Public Member Functions
	CylBase (double x, double y, double z, double dia, double length, unsigned int numDivs)

	CylBase ()

vector< double >	getCoordinates (const CylBase &parent, unsigned int entry) const

double	getDia () const

double	getDiffusionArea (const CylBase &parent, unsigned int index) const

bool	getIsCylinder () const

double	getLength () const

double	getMiddleArea (const CylBase &parent, unsigned int index) const
	Return cross-section area of middle of specified voxel. More...

unsigned int	getNumDivs () const

double	getVoxelLength () const
	Return length of voxel. All are equal. More...

double	getX () const

double	getY () const

double	getZ () const

void	matchCubeMeshEntries (const ChemCompt *other, const CylBase &parent, unsigned int startIndex, double granularity, vector< VoxelJunction > &ret, bool useCylinderCurve, bool useCylinderCap) const

double	nearest (double x, double y, double z, const CylBase &parent, double &linePos, double &r) const

double	nearest (double x, double y, double z, const CylBase &parent, unsigned int &index) const

double	selectGridSize (double h, double dia1, double granularity) const

void	setDia (double v)

void	setIsCylinder (bool v)

void	setLength (double v)

void	setNumDivs (unsigned int v)

void	setX (double v)

void	setY (double v)

void	setZ (double v)

double	volume (const CylBase &parent) const
	Returns vol of current node. Usually needs to refer to parent. More...

double	voxelVolume (const CylBase &parent, unsigned int fid) const

Private Attributes
double	dia_

bool	isCylinder_
	Number of subdivisions of cylinder. More...

double	length_
	Diameter of node end. More...

unsigned int	numDivs_
	Length of compartment. More...

double	x_
	end of the node. The start is given by parent coords. More...

double	y_

double	z_

Detailed Description

Base class for cylinder calculations.

Definition at line 17 of file CylBase.h.

Constructor & Destructor Documentation

CylBase::CylBase	(	double	x,
		double	y,
		double	z,
		double	dia,
		double	length,
		unsigned int	numDivs
	)

Definition at line 26 of file CylBase.cpp.

     :
         x_( x ),
         y_( y ),
         z_( z ),
         dia_( dia ),
         length_( length ),
         numDivs_( numDivs ),
         isCylinder_( false )
 {
     ;
 }

CylBase::CylBase ( )

Definition at line 40 of file CylBase.cpp.

     :
         x_( 0.0 ),
         y_( 0.0 ),
         z_( 0.0 ),
         dia_( 1.0 ),
         length_( 1.0 ),
         numDivs_( 1.0 ),
         isCylinder_( false )
 {
     ;
 }

Member Function Documentation

vector< double > CylBase::getCoordinates	(	const CylBase &	parent,
		unsigned int	fid
	)		const

Returns a vector with the coordinates of the specified MeshEntry. The MeshEntry (entry argument) is a voxel within this CylBase, 'entry' is indexed from 0 to numDivs. The parent is the parent CylBase with its own set of coords. The function does interpolation from the ends of the Node. These are detailed coordinates representative of the shape of the MeshEntry. For cylinders it is x1y1z1 x2y2z2 r0 r1 phi0 phi1 For spheres it is xyzdia.

Virtual function to return coords of mesh Entry. For Cylindrical mesh, coords are x1y1z1 x2y2z2 r0 r1 phi0 phi1

Definition at line 172 of file CylBase.cpp.

References dia_, numDivs_, x_, y_, and z_.

Referenced by NeuroMesh::vGetVoxelMidpoint().

 {
     assert( numDivs_ > fid );
     double frac0 = ( static_cast< double >( fid ) ) /
                 static_cast< double >( numDivs_ );
     double frac1 = ( static_cast< double >( fid + 1 ) ) /
                 static_cast< double >( numDivs_ );
 
     double r0 = 0.5 * ( parent.dia_ * ( 1.0 - frac0 ) + dia_ * frac0 );
     // equivalent: double r0 = parent.dia_ + frac0 * ( dia_ - parent.dia_ );
     double r1 = 0.5 * ( parent.dia_ * ( 1.0 - frac1 ) + dia_ * frac1 );
 
     vector< double > ret( 10, 0.0 );
     ret[0] = parent.x_ + frac0 * ( x_ - parent.x_ );
     ret[1] = parent.y_ + frac0 * ( y_ - parent.y_ );
     ret[2] = parent.z_ + frac0 * ( z_ - parent.z_ );
     ret[3] = parent.x_ + frac1 * ( x_ - parent.x_ );
     ret[4] = parent.y_ + frac1 * ( y_ - parent.y_ );
     ret[5] = parent.z_ + frac1 * ( z_ - parent.z_ );
     ret[6] = r0;
     ret[7] = r1;
     ret[8] = 0;
     ret[9] = 0;
 
     return ret;
 }

Here is the caller graph for this function:

double CylBase::getDia ( ) const

Definition at line 88 of file CylBase.cpp.

References dia_.

Referenced by NeuroNode::findStartNode(), SpineEntry::psdCoords(), SpineEntry::rootArea(), SpineEntry::setVolume(), and NeuroMesh::vSetVolumeNotRates().

 {
     return dia_;
 }

Here is the caller graph for this function:

double CylBase::getDiffusionArea	(	const CylBase &	parent,
		unsigned int	fid
	)		const

Returns diffusion cross-section from specified index to next. For index 0, this is cross-section to parent. For index numDivs-1, it is the cross-section from the second-last to the last voxel in this CylBase. Thus there is no valid value for (index == numDivs - 1), it has to be computed external to the CylBase, typically by calling the getDiffusionArea for the child CylBase.

Returns diffusion cross-section from specified index to previous. For index 0, this is cross-section area of parent. For index numDivs-1, it is the cross-section from the second-last to the last voxel in this CylBase. For index numDivs it is the area of this CylBase. Thus there is no valid value for (index > numDivs), it has to be computed external to the CylBase, typically by calling the getDiffusionArea for the child CylBase.

Definition at line 210 of file CylBase.cpp.

References dia_, isCylinder_, numDivs_, and PI.

Referenced by NeuroMesh::getAdx().

 {
     assert( fid < numDivs_ + 1 );
     if ( isCylinder_ )
             return PI * dia_ * dia_ / 4.0;
     double frac0 = ( static_cast< double >( fid ) ) /
                 static_cast< double >( numDivs_ );
     double r0 = 0.5 * ( parent.dia_ * ( 1.0 - frac0 ) + dia_ * frac0 );
     return PI * r0 * r0;
 }

Here is the caller graph for this function:

bool CylBase::getIsCylinder ( ) const

Definition at line 118 of file CylBase.cpp.

References isCylinder_.

 {
     return isCylinder_;
 }

double CylBase::getLength ( ) const

Definition at line 98 of file CylBase.cpp.

References length_.

Referenced by NeuroNode::calculateLength(), SpineEntry::diffusionLength(), NeuroMesh::getAdx(), SpineEntry::setVolume(), and NeuroMesh::vSetVolumeNotRates().

 {
     return length_;
 }

Here is the caller graph for this function:

double CylBase::getMiddleArea	(	const CylBase &	parent,
		unsigned int	index
	)		const

Return cross-section area of middle of specified voxel.

Return the cross section area of the middle of the specified voxel.

Definition at line 223 of file CylBase.cpp.

References dia_, isCylinder_, numDivs_, and PI.

Referenced by NeuroMesh::updateCoords().

 {
     assert( fid < numDivs_ );
     if ( isCylinder_ )
             return PI * dia_ * dia_ / 4.0;
     double frac0 = ( 0.5 + static_cast< double >( fid ) ) /
                 static_cast< double >( numDivs_ );
     double r0 = 0.5 * ( parent.dia_ * ( 1.0 - frac0 ) + dia_ * frac0 );
     return PI * r0 * r0;
 }

Here is the caller graph for this function:

unsigned int CylBase::getNumDivs ( ) const

Definition at line 108 of file CylBase.cpp.

References numDivs_.

Referenced by NeuroMesh::getAdx(), NeuroMesh::indexToSpace(), NeuroNode::isDummyNode(), NeuroMesh::nearest(), NeuroMesh::updateCoords(), and NeuroMesh::vGetVoxelMidpoint().

 {
     return numDivs_;
 }

Here is the caller graph for this function:

double CylBase::getVoxelLength ( ) const

Return length of voxel. All are equal.

Definition at line 235 of file CylBase.cpp.

References length_, and numDivs_.

Referenced by NeuroMesh::updateCoords().

 {
     return length_ / numDivs_;
 }

Here is the caller graph for this function:

double CylBase::getX ( ) const

Definition at line 58 of file CylBase.cpp.

References x_.

Referenced by NeuroNode::calculateLength(), NeuroMesh::indexToSpace(), SpineEntry::mid(), SpineEntry::positionShaftBase(), SpineEntry::psdCoords(), SpineEntry::setVolume(), and NeuroMesh::updateShaftParents().

 {
     return x_;
 }

Here is the caller graph for this function:

double CylBase::getY ( ) const

Definition at line 68 of file CylBase.cpp.

References y_.

Referenced by NeuroNode::calculateLength(), NeuroMesh::indexToSpace(), SpineEntry::mid(), SpineEntry::positionShaftBase(), SpineEntry::psdCoords(), SpineEntry::setVolume(), and NeuroMesh::updateShaftParents().

 {
     return y_;
 }

Here is the caller graph for this function:

double CylBase::getZ ( ) const

Definition at line 78 of file CylBase.cpp.

References z_.

Referenced by NeuroNode::calculateLength(), NeuroMesh::indexToSpace(), SpineEntry::mid(), SpineEntry::positionShaftBase(), SpineEntry::psdCoords(), SpineEntry::setVolume(), and NeuroMesh::updateShaftParents().

 {
     return z_;
 }

Here is the caller graph for this function:

void CylBase::matchCubeMeshEntries	(	const ChemCompt *	other,
		const CylBase &	parent,
		unsigned int	startIndex,
		double	granularity,
		vector< VoxelJunction > &	ret,
		bool	useCylinderCurve,
		bool	useCylinderCap
	)		const

Definition at line 320 of file CylBase.cpp.

References Vec::a0(), Vec::a1(), Vec::a2(), dia_, EPSILON, fillPointsOnCircle(), fillPointsOnDisc(), CubeMesh::getDx(), ChemCompt::getNumEntries(), isCylinder_, length_, numDivs_, Vec::orthogonalAxes(), selectGridSize(), x_, y_, and z_.

Referenced by SpineEntry::matchCubeMeshEntries(), PsdMesh::matchCubeMeshEntries(), NeuroMesh::matchCubeMeshEntries(), SpineEntry::matchCubeMeshEntriesToHead(), and SpineEntry::matchCubeMeshEntriesToPSD().

 {
     const CubeMesh* other = dynamic_cast< const CubeMesh* >( compt );
     assert( other );
     const double EPSILON = 1e-18;
     Vec a( parent.x_ - x_, parent.y_ - y_, parent.z_ - z_ );
     Vec u;
     Vec v;
     a.orthogonalAxes( u, v );
 
     double h = selectGridSize( other->getDx(), parent.dia_/2, granularity );
     double lambda = length_ / numDivs_;
 
     unsigned int num = floor( 0.1 + lambda / h );
     // March along axis of cylinder.
     // q is the location of the point along axis.
     double rSlope = ( dia_ - parent.dia_ ) * 0.5 / length_;
     for ( unsigned int i = 0; i < numDivs_; ++i ) {
         vector< double >area( other->getNumEntries(), 0.0 );
         if ( useCylinderCurve ) {
             for ( unsigned int j = 0; j < num; ++j ) {
                 unsigned int m = i * num + j;
                 double frac = ( m * h + h/2.0 ) / length_;
                 double q0 = x_ + a.a0() * frac;
                 double q1 = y_ + a.a1() * frac;
                 double q2 = z_ + a.a2() * frac;
                 // get radius of cylinder at this point.
                 double r = dia_/2.0;
                 if ( !isCylinder_ ) // Use the more complicated conic value
                 r = parent.dia_/2.0 + frac * rSlope;
                 fillPointsOnCircle( u, v, Vec( q0, q1, q2 ),
                             h, r, area, other );
             }
         }
         if ( useCylinderCap && i == numDivs_ - 1 ) {
             fillPointsOnDisc( u, v, Vec( x_, y_, z_ ),
                             h, dia_/2.0, area, other );
         }
         // Go through all cubeMesh entries and compute diffusion
         // cross-section. Assume this is through a membrane, so the
         // only factor relevant is area. Not the distance.
         for ( unsigned int k = 0; k < area.size(); ++k ) {
             if ( area[k] > EPSILON ) {
                 ret.push_back( VoxelJunction( i + startIndex, k, area[k] ));
             }
         }
     }
 }

Here is the call graph for this function:

Here is the caller graph for this function:

double CylBase::nearest	(	double	x,
		double	y,
		double	z,
		const CylBase &	parent,
		double &	linePos,
		double &	r
	)		const

Definition at line 378 of file CylBase.cpp.

References dia_, Vec::distance(), EPSILON, numDivs_, Vec::pointOnLine(), x_, y_, and z_.

Referenced by nearest(), NeuroMesh::nearest(), and NeuroMesh::updateShaftParents().

 {
     const double EPSILON = 1e-8; // 10 nm is too small for any compt.
     // Consider a = parent and b = self, and c = x,y,z.
     // Fraction along cylinder axis ab = k
     // k = (c - a).(b - a)/((b-a).(b-a))
     //
     // Then, point p along line from parent to self is
     // p = k( self-parent) + parent.
     //
     // So distance from c to p is what we want.
     //
     // If k is
     //
     Vec a( parent.x_, parent.y_, parent.z_ );
     Vec b( x_, y_, z_ );
     Vec c( x, y, z );
 
     double dist = b.distance( a );
     assert( dist > EPSILON );
     double k = ( b - a ).dotProduct( c - a );
     k /= dist * dist;
     Vec pt = a.pointOnLine( b, k );
 
     double ret = c.distance(pt);
     linePos = k;
     // For now I disable the positioning with a tapering cylinder.
     // double rSlope = 0.5 * ( dia_ - parent.dia_ ) / numDivs_;
     double rSlope = 0;
     r = dia_/2.0 + k * numDivs_ * rSlope;
     return ret;
 }

Here is the call graph for this function:

Here is the caller graph for this function:

double CylBase::nearest	(	double	x,
		double	y,
		double	z,
		const CylBase &	parent,
		unsigned int &	index
	)		const

Definition at line 414 of file CylBase.cpp.

References nearest(), and numDivs_.

 {
     double k = 0.0;
     double r;
     double ret = nearest( x, y, z, parent, k, r );
     if ( k < 0.0 ) {
         ret = -ret;
         index = 0;
     } else if ( k > 1.0 ) {
         ret = -ret;
         index = numDivs_ - 1;
     } else { // Inside length of cylinder, now is it inside radius?
         index = k * numDivs_;
         if ( index >= numDivs_ )
             index = numDivs_ - 1;
         // double ri = r0_ + (index + 0.5) * rSlope_;
         if ( ret > r * 1.01 ) // Handle roundoff
             ret = -ret;
     }
     return ret;
 }

Here is the call graph for this function:

double CylBase::selectGridSize	(	double	h,
		double	dia1,
		double	granularity
	)		const

Utility commands for building interface to cube mesh

Definition at line 242 of file CylBase.cpp.

References dia_, length_, and numDivs_.

Referenced by matchCubeMeshEntries().

 {
     if ( length_ < 1e-7 && numDivs_ == 1 ) { // It is a disc, not a cylinder
         return granularity * dia_ / 2.0;
     }
 
     double lambda = length_ / numDivs_;
     if ( h > lambda )
         h = lambda;
     if ( h > dia_ / 2.0 )
         h = dia_ / 2.0;
     if ( h > dia1/2.0 )
         h = dia1/2.0;
     h *= granularity;
     unsigned int num = ceil( lambda / h );
     h = lambda / num;
 
     return h;
 }

Here is the caller graph for this function:

void CylBase::setDia ( double v )

Definition at line 83 of file CylBase.cpp.

References dia_.

Referenced by NeuroNode::NeuroNode(), SpineEntry::setVolume(), SpineEntry::SpineEntry(), and NeuroMesh::vSetVolumeNotRates().

 {
     dia_ = v;
 }

Here is the caller graph for this function:

void CylBase::setIsCylinder ( bool v )

Definition at line 113 of file CylBase.cpp.

References isCylinder_.

Referenced by NeuroMesh::insertSingleDummy(), and SpineEntry::SpineEntry().

 {
     isCylinder_ = v;
 }

Here is the caller graph for this function:

void CylBase::setLength ( double v )

Definition at line 93 of file CylBase.cpp.

References length_.

Referenced by NeuroNode::calculateLength(), NeuroNode::NeuroNode(), SpineEntry::setVolume(), SpineEntry::SpineEntry(), and NeuroMesh::vSetVolumeNotRates().

 {
     length_ = v;
 }

Here is the caller graph for this function:

void CylBase::setNumDivs ( unsigned int v )

Definition at line 103 of file CylBase.cpp.

References numDivs_.

Referenced by NeuroMesh::insertSingleDummy(), and SpineEntry::SpineEntry().

 {
     numDivs_ = v;
 }

Here is the caller graph for this function:

void CylBase::setX ( double v )

Definition at line 53 of file CylBase.cpp.

References x_.

Referenced by NeuroMesh::insertSingleDummy(), NeuroNode::NeuroNode(), SpineEntry::positionShaftBase(), SpineEntry::setVolume(), and SpineEntry::SpineEntry().

 {
     x_ = v;
 }

Here is the caller graph for this function:

void CylBase::setY ( double v )

Definition at line 63 of file CylBase.cpp.

References y_.

Referenced by NeuroMesh::insertSingleDummy(), NeuroNode::NeuroNode(), SpineEntry::positionShaftBase(), SpineEntry::setVolume(), and SpineEntry::SpineEntry().

 {
     y_ = v;
 }

Here is the caller graph for this function:

void CylBase::setZ ( double v )

Definition at line 73 of file CylBase.cpp.

References z_.

Referenced by NeuroMesh::insertSingleDummy(), NeuroNode::NeuroNode(), SpineEntry::positionShaftBase(), SpineEntry::setVolume(), and SpineEntry::SpineEntry().

 {
     z_ = v;
 }

Here is the caller graph for this function:

double CylBase::volume ( const CylBase & parent ) const

Returns vol of current node. Usually needs to refer to parent.

The entire volume for a truncated cone is given by: V = 1/3 pi.length.(r0^2 + r0.r1 + r1^2) where the length is the length of the cone r0 is radius at base r1 is radius at top. Note that this converges to volume of a cone if r0 or r1 is zero, and to the volume of a cylinder if r0 == r1.

Definition at line 135 of file CylBase.cpp.

References dia_, isCylinder_, length_, and PI.

Referenced by SpineEntry::setVolume(), SpineEntry::volume(), and NeuroMesh::vSetVolumeNotRates().

 {
     if ( isCylinder_ )
             return length_ * dia_ * dia_ * PI / 4.0;
     double r0 = parent.dia_/2.0;
     double r1 = dia_/2.0;
     return length_ * ( r0*r0 + r0 *r1 + r1 * r1 ) * PI / 3.0;
 }

Here is the caller graph for this function:

double CylBase::voxelVolume	(	const CylBase &	parent,
		unsigned int	fid
	)		const

Returns volume of specified voxel (MeshEntry)

Returns volume of MeshEntry. This isn't the best subdivision of the cylinder from the viewpoint of keeping the length constants all the same for different volumes. Ideally the thinner segments should have a smaller length. But this is simple and so is the diffusion calculation, so leave it. Easy to fine-tune later by modifying how one computes frac0 and frac1.

Definition at line 152 of file CylBase.cpp.

References dia_, isCylinder_, length_, numDivs_, PI, and s0.

Referenced by NeuroMesh::updateCoords().

 {
     assert( numDivs_ > fid );
     if ( isCylinder_ )
             return length_ * dia_ * dia_ * PI / ( 4.0 * numDivs_ );
 
     double frac0 = ( static_cast< double >( fid ) ) /
                 static_cast< double >( numDivs_ );
     double frac1 = ( static_cast< double >( fid + 1 ) ) /
                 static_cast< double >( numDivs_ );
     double r0 = 0.5 * ( parent.dia_ * ( 1.0 - frac0 ) + dia_ * frac0 );
     double r1 = 0.5 * ( parent.dia_ * ( 1.0 - frac1 ) + dia_ * frac1 );
     double s0 = length_ * frac0;
     double s1 = length_ * frac1;
 
     return (s1 - s0) * ( r0*r0 + r0 *r1 + r1 * r1 ) * PI / 3.0;
 }