plane.H 6.49 KB
Newer Older
1
2
3
4
/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     |
5
    \\  /    A nd           | Copyright (C) 2011-2016 OpenFOAM Foundation
6
     \\/     M anipulation  | Copyright (C) 2017 OpenCFD Ltd.
7
8
9
10
-------------------------------------------------------------------------------
License
    This file is part of OpenFOAM.

11
12
13
14
    OpenFOAM is free software: you can redistribute it and/or modify it
    under the terms of the GNU General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.
15
16
17
18
19
20
21

    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
    for more details.

    You should have received a copy of the GNU General Public License
22
    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
23
24
25
26
27

Class
    Foam::plane

Description
28
    Geometric class that creates a 3D plane and can return the intersection
29
    point between a line and the plane.
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51

SourceFiles
    plane.C

\*---------------------------------------------------------------------------*/

#ifndef plane_H
#define plane_H

#include "point.H"
#include "scalarList.H"
#include "dictionary.H"
#include "line.H"

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

namespace Foam
{

// Forward declaration of friend functions and operators

class plane;
52
Ostream& operator<<(Ostream& os, const plane& pln);
53
54
55


/*---------------------------------------------------------------------------*\
56
                            Class plane Declaration
57
58
59
60
61
62
\*---------------------------------------------------------------------------*/

class plane
{
public:

63
64
65
66
67
68
69
    //- Side of the plane
    enum side
    {
        NORMAL,
        FLIP
    };

70
71
72
73
74
    //- A direction and a reference point
    class ray
    {
        point refPoint_;
        vector dir_;
75

76
    public:
77

78
79
80
81
82
        ray(const point& refPoint, const vector& dir)
        :
            refPoint_(refPoint),
            dir_(dir)
        {}
83

84
85
86
87
        const point& refPoint() const
        {
            return refPoint_;
        }
88

89
90
91
92
93
        const vector& dir() const
        {
            return dir_;
        }
    };
94
95
96
97
98
99


private:

    // Private data

100
101
        //- Normal
        vector normal_;
102

103
104
        //- Reference point
        point point_;
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125


    // Private Member Functions

        //- Calculates basePoint and normal vector given plane coefficients
        void calcPntAndVec(const scalarList& C);

        //- Calculates basePoint and normal vector given three points
        //- Normal vector determined using right hand rule
        void calcPntAndVec
        (
            const point& point1,
            const point& point2,
            const point& point3
        );


public:

    // Constructors

126
127
128
        //- Construct null. Does not set normal and point.
        plane();

129
        //- Construct from normal vector through the origin
130
        explicit plane(const vector& normalVector);
131
132

        //- Construct from normal vector and point in plane
133
        explicit plane
134
135
136
137
138
        (
            const point& basePoint,
            const vector& normalVector,
            const bool normalise = true
        );
139
140

        //- Construct from three points in plane
141
        plane
142
143
144
145
146
        (
            const point& point1,
            const point& point2,
            const point& point3
        );
147
148
149

        //- Construct from coefficients for the
        //  plane equation: ax + by + cz + d = 0
150
        explicit plane(const scalarList& C);
151
152

        //- Construct from dictionary
153
        explicit plane(const dictionary& planeDict);
154
155

        //- Construct from Istream. Assumes the base + normal notation.
156
        explicit plane(Istream& is);
157
158
159
160
161
162
163


    // Member Functions

        //- Return plane normal
        const vector& normal() const;

164
        //- Return plane base point
165
166
167
168
        const point& refPoint() const;

        //- Return coefficients for the
        //  plane equation: ax + by + cz + d = 0
169
        FixedList<scalar, 4> planeCoeffs() const;
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196

        //- Return nearest point in the plane for the given point
        point nearestPoint(const point& p) const;

        //- Return distance from the given point to the plane
        scalar distance(const point& p) const;

        //- Return cut coefficient for plane and line defined by
        //  origin and direction
        scalar normalIntersect(const point& pnt0, const vector& dir) const;

        //- Return cut coefficient for plane and ray
        scalar normalIntersect(const ray& r) const
        {
            return normalIntersect(r.refPoint(), r.dir());
        }

        //- Return the cutting point between the plane and
        // a line passing through the supplied points
        template<class Point, class PointRef>
        scalar lineIntersect(const line<Point, PointRef>& l) const
        {
            return normalIntersect(l.start(), l.vec());
        }

        //- Return the cutting line between this plane and another.
        //  Returned as direction vector and point line goes through.
197
        ray planeIntersect(const plane& plane2) const;
198
199

        //- Return the cutting point between this plane and two other planes
200
201
202
203
204
205
206
207
208
        point planePlaneIntersect
        (
            const plane& plane2,
            const plane& plane3
        )
        const;

        //- Return a point somewhere on the plane, a distance from the base
        point somePointInPlane(const scalar dist = 1e-3) const;
209

210
211
212
213
        //- Return the side of the plane that the point is on.
        //  If the point is on the plane, then returns NORMAL.
        side sideOfPlane(const point& p) const;

214
215
216
        //- Mirror the supplied point in the plane. Return the mirrored point.
        point mirror(const point& p) const;

217
        //- Write to dictionary
218
        void writeDict(Ostream& os) const;
219

220
221
222
223

    // IOstream Operators

        //- Write plane properties
224
225
        friend Ostream& operator<<(Ostream& os, const plane& pln);

226
227
228
};


229
230
231
232
233
234
// Global Operators

bool operator==(const plane& a, const plane& b);
bool operator!=(const plane& a, const plane& b);


235
236
237
238
239
240
241
242
243
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

} // End namespace Foam

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

#endif

// ************************************************************************* //