1 // VISU OBJECT : interactive object for VISU entities implementation
3 // Copyright (C) 2003 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
4 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
6 // This library is free software; you can redistribute it and/or
7 // modify it under the terms of the GNU Lesser General Public
8 // License as published by the Free Software Foundation; either
9 // version 2.1 of the License.
11 // This library is distributed in the hope that it will be useful,
12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
14 // Lesser General Public License for more details.
16 // You should have received a copy of the GNU Lesser General Public
17 // License along with this library; if not, write to the Free Software
18 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
20 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
23 // File: VISU_PipeLine.cxx
24 // Author: Alexey PETROV
28 #include "VISU_CutPlanesPL.hxx"
29 #include "VISU_PipeLineUtils.hxx"
30 #include "VTKViewer_GeometryFilter.h"
32 #include <vtkAppendPolyData.h>
33 #include <vtkCutter.h>
38 static vtkFloatingPointType EPS = 1.0E-3;
40 vtkStandardNewMacro(VISU_CutPlanesPL);
45 myAppendPolyData = vtkAppendPolyData::New();
46 myIsShrinkable = false;
53 myDisplacement[0] = myDisplacement[1] = 0.5;
55 myAng[0][0] = myAng[0][1] = myAng[0][2] = 0.0;
56 myAng[1][0] = myAng[1][1] = myAng[1][2] = 0.0;
62 myAppendPolyData->Delete();
67 ::ShallowCopy(VISU_PipeLine *thePipeLine)
69 if(VISU_CutPlanesPL *aPipeLine = dynamic_cast<VISU_CutPlanesPL*>(thePipeLine)){
70 SetOrientation(aPipeLine->GetPlaneOrientation(),
71 aPipeLine->GetRotateX(),aPipeLine->GetRotateY());
72 SetDisplacement(aPipeLine->GetDisplacement());
73 SetNbParts(aPipeLine->GetNbParts());
74 for (int i = 0, iend = GetNbParts(); i < iend; i++)
75 if(!aPipeLine->IsPartDefault(i)) SetPartPosition(i, aPipeLine->GetPartPosition(i));
77 VISU_ScalarMapPL::ShallowCopy(thePipeLine);
84 VISU_ScalarMapPL::Init();
88 myDisplacement[0] = 0.5;
89 myAng[0][0] = myAng[0][1] = myAng[0][2] = 0.0;
92 VISU_ScalarMapPL::THook*
96 return myAppendPolyData->GetOutput();
103 ClearAppendPolyData(myAppendPolyData);
105 vtkFloatingPointType aDir[3];
106 GetDir(aDir,myAng[0],myBasePlane[0]);
107 vtkFloatingPointType aBounds[6];
108 GetInput2()->GetBounds(aBounds);
109 vtkDataSet* aDataSet = myFieldTransform->GetUnstructuredGridOutput();
110 CutWithPlanes(myAppendPolyData,aDataSet,myNbParts,aDir,aBounds,
111 myPartPosition,myPartCondition,myDisplacement[0]);
113 VISU_ScalarMapPL::Update();
118 ::SetPartPosition(int theNum)
120 for(int i = 0; i < myNbParts; i++)
121 myPartPosition[i] = GetPartPosition(i,theNum);
126 ::ClearAppendPolyData(vtkAppendPolyData *theAppendPolyData)
128 int iEnd = theAppendPolyData->GetNumberOfInputs();
129 for(int i = iEnd-1; i >= 0; i--)
130 theAppendPolyData->RemoveInput(theAppendPolyData->GetInput(i));
133 vtkFloatingPointType*
135 GetRx(vtkFloatingPointType theRx[3][3],
136 vtkFloatingPointType thaAng)
138 theRx[0][0] = 1.0; theRx[0][1] = 0.0; theRx[0][2] = 0.0;
139 theRx[1][0] = 0.0; theRx[1][1] = cos(thaAng); theRx[1][2] = -sin(thaAng);
140 theRx[2][0] = 0.0; theRx[2][1] = sin(thaAng); theRx[2][2] = cos(thaAng);
145 vtkFloatingPointType*
147 ::GetRy(vtkFloatingPointType theRy[3][3],
148 vtkFloatingPointType thaAng)
150 theRy[0][0] = cos(thaAng); theRy[0][1] = 0.0; theRy[0][2] = sin(thaAng);
151 theRy[1][0] = 0.0; theRy[1][1] = 1.0; theRy[1][2] = 0.0;
152 theRy[2][0] = -sin(thaAng); theRy[2][1] = 0.0; theRy[2][2] = cos(thaAng);
157 vtkFloatingPointType*
159 ::GetRz(vtkFloatingPointType theRz[3][3],
160 vtkFloatingPointType thaAng)
162 theRz[0][0] = cos(thaAng); theRz[0][1] = -sin(thaAng); theRz[0][2] = 0.0;
163 theRz[1][0] = sin(thaAng); theRz[1][1] = cos(thaAng); theRz[1][2] = 0.0;
164 theRz[2][0] = 0.0; theRz[2][1] = 0.0; theRz[2][2] = 1.0;
171 ::CorrectPnt(vtkFloatingPointType thePnt[3],
172 const vtkFloatingPointType BoundPrj[6])
174 for(int i = 0, j = 0; i < 3; ++i, j=2*i){
175 if(thePnt[i] < BoundPrj[j]) thePnt[i] = BoundPrj[j];
176 if(thePnt[i] > BoundPrj[j+1]) thePnt[i] = BoundPrj[j+1];
182 ::GetBoundProject(vtkFloatingPointType BoundPrj[3],
183 const vtkFloatingPointType BoundBox[6],
184 const vtkFloatingPointType Dir[3])
186 vtkFloatingPointType BoundPoints[8][3] = { {BoundBox[0],BoundBox[2],BoundBox[4]},
187 {BoundBox[1],BoundBox[2],BoundBox[4]},
188 {BoundBox[0],BoundBox[3],BoundBox[4]},
189 {BoundBox[1],BoundBox[3],BoundBox[4]},
190 {BoundBox[0],BoundBox[2],BoundBox[5]},
191 {BoundBox[1],BoundBox[2],BoundBox[5]},
192 {BoundBox[0],BoundBox[3],BoundBox[5]},
193 {BoundBox[1],BoundBox[3],BoundBox[5]}};
194 BoundPrj[0] = vtkMath::Dot(Dir,BoundPoints[0]), BoundPrj[1] = BoundPrj[0];
195 for(int i = 1; i < 8; i++){
196 vtkFloatingPointType tmp = vtkMath::Dot(Dir,BoundPoints[i]);
197 if(BoundPrj[1] < tmp) BoundPrj[1] = tmp;
198 if(BoundPrj[0] > tmp) BoundPrj[0] = tmp;
200 BoundPrj[2] = BoundPrj[1] - BoundPrj[0];
201 BoundPrj[1] = BoundPrj[0] + (1.0 - EPS)*BoundPrj[2];
202 BoundPrj[0] = BoundPrj[0] + EPS*BoundPrj[2];
203 BoundPrj[2] = BoundPrj[1] - BoundPrj[0];
209 ::SetOrientation(const VISU_CutPlanesPL::PlaneOrientation& theOrient,
210 vtkFloatingPointType theXAng,
211 vtkFloatingPointType theYAng,
214 myBasePlane[theNum] = theOrient;
215 switch(myBasePlane[theNum]){
216 case XY: myAng[theNum][0] = theXAng; break;
217 case YZ: myAng[theNum][1] = theXAng; break;
218 case ZX: myAng[theNum][2] = theXAng; break;
220 switch(myBasePlane[theNum]){
221 case XY: myAng[theNum][1] = theYAng; break;
222 case YZ: myAng[theNum][2] = theYAng; break;
223 case ZX: myAng[theNum][0] = theYAng; break;
228 const VISU_CutPlanesPL::PlaneOrientation&
230 ::GetPlaneOrientation(int theNum)
232 return myBasePlane[theNum];
237 ::GetRotateX(int theNum)
239 switch(myBasePlane[theNum]){
240 case XY: return myAng[theNum][0];
241 case YZ: return myAng[theNum][1];
242 case ZX: return myAng[theNum][2];
249 ::GetRotateY(int theNum)
251 switch(myBasePlane[theNum]){
252 case XY: return myAng[theNum][1];
253 case YZ: return myAng[theNum][2];
254 case ZX: return myAng[theNum][0];
262 ::SetNbParts(int theNb)
265 myPartPosition.resize(myNbParts);
266 myPartCondition.resize(myNbParts,1);
273 ::SetPartPosition(int thePartNumber,
274 vtkFloatingPointType thePartPosition)
276 if(thePartNumber >= myNbParts) return;
277 myPartPosition[thePartNumber] = thePartPosition;
278 myPartCondition[thePartNumber] = 0;
284 ::GetPartPosition(int thePartNumber,
287 if(thePartNumber >= myNbParts) return 0;
288 vtkFloatingPointType aPosition = myPartPosition[thePartNumber];
289 if(myPartCondition[thePartNumber]){
290 vtkFloatingPointType aDir[3], aBounds[6], aBoundPrj[3];
291 GetInput2()->GetBounds(aBounds);
292 GetDir(aDir,myAng[theNum],myBasePlane[theNum]);
293 GetBoundProject(aBoundPrj,aBounds,aDir);
295 vtkFloatingPointType aDBoundPrj = aBoundPrj[2]/(myNbParts - 1);
296 vtkFloatingPointType aDisplacement = aDBoundPrj * myDisplacement[theNum];
297 vtkFloatingPointType aStartPosition = aBoundPrj[0] - 0.5*aDBoundPrj + aDisplacement;
298 aPosition = aStartPosition + thePartNumber*aDBoundPrj;
300 aPosition = aBoundPrj[0] + aBoundPrj[2]*myDisplacement[theNum];
308 ::SetPartDefault(int thePartNumber)
310 if(thePartNumber >= myNbParts) return;
311 myPartPosition[thePartNumber] = GetPartPosition(thePartNumber);
312 myPartCondition[thePartNumber] = 1;
317 ::IsPartDefault(int thePartNumber)
319 if(thePartNumber >= myNbParts) return 1;
320 return myPartCondition[thePartNumber];
326 ::GetDir(vtkFloatingPointType theDir[3],
327 const vtkFloatingPointType theAng[3],
328 const PlaneOrientation& theBasePlane)
331 vtkFloatingPointType aRx[3][3], aRy[3][3], aRz[3][3], aRotation[3][3];
332 switch(theBasePlane){
334 if(fabs(theAng[0]) > EPS) GetRx(aRx,theAng[0]); else vtkMath::Identity3x3(aRx);
335 if(fabs(theAng[1]) > EPS) GetRy(aRy,theAng[1]); else vtkMath::Identity3x3(aRy);
336 vtkMath::Multiply3x3(aRx,aRy,aRotation);
340 if(fabs(theAng[1]) > EPS) GetRy(aRy,theAng[1]); else vtkMath::Identity3x3(aRy);
341 if(fabs(theAng[2]) > EPS) GetRz(aRz,theAng[2]); else vtkMath::Identity3x3(aRz);
342 vtkMath::Multiply3x3(aRy,aRz,aRotation);
346 if(fabs(theAng[2]) > EPS) GetRz(aRz,theAng[2]); else vtkMath::Identity3x3(aRz);
347 if(fabs(theAng[0]) > EPS) GetRx(aRx,theAng[0]); else vtkMath::Identity3x3(aRx);
348 vtkMath::Multiply3x3(aRz,aRx,aRotation);
352 for(int i = 0; i < 3; i++)
353 theDir[i] = aRotation[i][iPlane];
359 ::CutWithPlane(vtkAppendPolyData* theAppendPolyData,
360 vtkDataSet* theDataSet,
361 vtkFloatingPointType theDir[3],
362 vtkFloatingPointType theOrig[3])
364 vtkCutter *aCutPlane = vtkCutter::New();
365 aCutPlane->SetInput(theDataSet);
366 vtkPlane *aPlane = vtkPlane::New();
367 aPlane->SetOrigin(theOrig);
369 aPlane->SetNormal(theDir);
370 aCutPlane->SetCutFunction(aPlane);
372 theAppendPolyData->AddInput(aCutPlane->GetOutput());
373 aCutPlane->Register(theAppendPolyData);
380 ::CutWithPlanes(vtkAppendPolyData* theAppendPolyData,
381 vtkDataSet* theDataSet,
383 vtkFloatingPointType theDir[3],
384 vtkFloatingPointType theBounds[6],
385 const vector<vtkFloatingPointType>& thePlanePosition,
386 const vector<int>& thePlaneCondition,
387 vtkFloatingPointType theDisplacement)
389 vtkFloatingPointType aBoundPrj[3], aOrig[3], aPosition;
390 GetBoundProject(aBoundPrj, theBounds, theDir);
392 vtkFloatingPointType aDBoundPrj = aBoundPrj[2]/(theNbPlanes - 1);
393 vtkFloatingPointType aDisplacement = aDBoundPrj*theDisplacement;
394 vtkFloatingPointType aStartPosition = aBoundPrj[0] - 0.5*aDBoundPrj + aDisplacement;
395 for (int i = 0; i < theNbPlanes; i++){
396 aPosition = aStartPosition + i*aDBoundPrj;
397 if(thePlaneCondition[i]){
398 aPosition = aStartPosition + i*aDBoundPrj;
400 aPosition = thePlanePosition[i];
401 VISU::Mul(theDir,aPosition,aOrig);
402 CutWithPlane(theAppendPolyData,theDataSet,theDir,aOrig);
405 if(thePlaneCondition[0])
406 aPosition = aBoundPrj[0] + aBoundPrj[2]*theDisplacement;
408 aPosition = thePlanePosition[0];
409 VISU::Mul(theDir,aPosition,aOrig);
410 CutWithPlane(theAppendPolyData,theDataSet,theDir,aOrig);
412 vtkPolyData *aPolyData = theAppendPolyData->GetOutput();
414 theAppendPolyData->Update();