FreeFOAM The Cross-Platform CFD Toolkit
Hosted by SourceForge:
Get FreeFOAM at SourceForge.net.
            Fast, secure and Free Open Source software downloads

polyLine.C

Go to the documentation of this file.
00001 /*---------------------------------------------------------------------------*\
00002   =========                 |
00003   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
00004    \\    /   O peration     |
00005     \\  /    A nd           | Copyright (C) 1991-2010 OpenCFD Ltd.
00006      \\/     M anipulation  |
00007 -------------------------------------------------------------------------------
00008 License
00009     This file is part of OpenFOAM.
00010 
00011     OpenFOAM is free software: you can redistribute it and/or modify it
00012     under the terms of the GNU General Public License as published by
00013     the Free Software Foundation, either version 3 of the License, or
00014     (at your option) any later version.
00015 
00016     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
00017     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
00018     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
00019     for more details.
00020 
00021     You should have received a copy of the GNU General Public License
00022     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
00023 
00024 \*---------------------------------------------------------------------------*/
00025 
00026 #include <OpenFOAM/error.H>
00027 #include "polyLine.H"
00028 
00029 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
00030 
00031 void Foam::polyLine::calcParam()
00032 {
00033     param_.setSize(points_.size());
00034 
00035     if (param_.size())
00036     {
00037         param_[0] = 0.0;
00038 
00039         for (label i=1; i < param_.size(); i++)
00040         {
00041             param_[i] = param_[i-1] + mag(points_[i] - points_[i-1]);
00042         }
00043 
00044         // normalize on the interval 0-1
00045         lineLength_ = param_[param_.size()-1];
00046         for (label i=1; i < param_.size() - 1; i++)
00047         {
00048             param_[i] /= lineLength_;
00049         }
00050         param_[param_.size()-1] = 1.0;
00051     }
00052     else
00053     {
00054         lineLength_ = 0.0;
00055     }
00056 }
00057 
00058 
00059 
00060 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
00061 
00062 Foam::polyLine::polyLine(const pointField& ps, const bool)
00063 :
00064     points_(ps),
00065     lineLength_(0.0),
00066     param_(0)
00067 {
00068     calcParam();
00069 }
00070 
00071 
00072 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
00073 
00074 const Foam::pointField& Foam::polyLine::points() const
00075 {
00076     return points_;
00077 }
00078 
00079 
00080 Foam::label Foam::polyLine::nSegments() const
00081 {
00082     return points_.size()-1;
00083 }
00084 
00085 
00086 Foam::label Foam::polyLine::localParameter(scalar& lambda) const
00087 {
00088     // check endpoints
00089     if (lambda < SMALL)
00090     {
00091         lambda = 0;
00092         return 0;
00093     }
00094     else if (lambda > 1 - SMALL)
00095     {
00096         lambda = 1;
00097         return nSegments();
00098     }
00099 
00100     // search table of cumulative distances to find which line-segment
00101     // we are on. Check the upper bound.
00102 
00103     label segmentI = 1;
00104     while (param_[segmentI] < lambda)
00105     {
00106         segmentI++;
00107     }
00108     segmentI--;   // we want the corresponding lower bound
00109 
00110     // the local parameter [0-1] on this line segment
00111     lambda =
00112     (
00113         ( lambda - param_[segmentI] )
00114       / ( param_[segmentI+1] - param_[segmentI] )
00115     );
00116 
00117     return segmentI;
00118 }
00119 
00120 
00121 Foam::point Foam::polyLine::position(const scalar mu) const
00122 {
00123     // check endpoints
00124     if (mu < SMALL)
00125     {
00126         return points_[0];
00127     }
00128     else if (mu > 1 - SMALL)
00129     {
00130         return points_[points_.size()-1];
00131     }
00132 
00133 
00134     scalar lambda = mu;
00135     label segment = localParameter(lambda);
00136     return position(segment, lambda);
00137 }
00138 
00139 
00140 Foam::point Foam::polyLine::position
00141 (
00142     const label segment,
00143     const scalar mu
00144 ) const
00145 {
00146     // out-of-bounds
00147     if (segment < 0)
00148     {
00149         return points_[0];
00150     }
00151     else if (segment > nSegments())
00152     {
00153         return points_[points_.size()-1];
00154     }
00155 
00156     const point& p0 = points()[segment];
00157     const point& p1 = points()[segment+1];
00158 
00159     // special cases - no calculation needed
00160     if (mu <= 0.0)
00161     {
00162         return p0;
00163     }
00164     else if (mu >= 1.0)
00165     {
00166         return p1;
00167     }
00168     else
00169     {
00170         // linear interpolation
00171         return points_[segment] + mu * (p1 - p0);
00172     }
00173 }
00174 
00175 
00176 Foam::scalar Foam::polyLine::length() const
00177 {
00178     return lineLength_;
00179 }
00180 
00181 
00182 // ************************ vim: set sw=4 sts=4 et: ************************ //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines