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_Plot3DPL.hxx"
29 #include "VISU_CutPlanesPL.hxx"
30 #include "VISU_PipeLineUtils.hxx"
32 #include <vtkAppendPolyData.h>
33 #include <vtkCutter.h>
36 #include <vtkCellDataToPointData.h>
37 #include <vtkGeometryFilter.h>
38 #include <vtkContourFilter.h>
39 #include <vtkWarpScalar.h>
40 #include <vtkOutlineFilter.h>
43 //----------------------------------------------------------------------------
44 vtkStandardNewMacro(VISU_Plot3DPL);
47 //----------------------------------------------------------------------------
50 myCellDataToPointData(vtkCellDataToPointData::New()),
51 myAppendPolyData(vtkAppendPolyData::New()),
52 myGeometryFilter(vtkGeometryFilter::New()),
53 myContourFilter(vtkContourFilter::New()),
54 myWarpScalar(vtkWarpScalar::New()),
55 myOrientation(VISU_CutPlanesPL::YZ),
61 SetIsShrinkable(false);
63 myCellDataToPointData->Delete();
64 myAppendPolyData->Delete();
65 myGeometryFilter->Delete();
66 myContourFilter->Delete();
67 myWarpScalar->Delete();
69 myAngle[0] = myAngle[1] = myAngle[2] = 0.0;
71 SetNumberOfContours(32);
75 //----------------------------------------------------------------------------
81 //----------------------------------------------------------------------------
86 unsigned long int aTime = Superclass::GetMTime();
88 aTime = std::max(aTime, myCellDataToPointData->GetMTime());
89 aTime = std::max(aTime, myAppendPolyData->GetMTime());
90 aTime = std::max(aTime, myGeometryFilter->GetMTime());
91 aTime = std::max(aTime, myContourFilter->GetMTime());
92 aTime = std::max(aTime, myWarpScalar->GetMTime());
98 //----------------------------------------------------------------------------
101 ::DoShallowCopy(VISU_PipeLine *thePipeLine,
104 Superclass::DoShallowCopy(thePipeLine, theIsCopyInput);
106 if(VISU_Plot3DPL *aPipeLine = dynamic_cast<VISU_Plot3DPL*>(thePipeLine)){
107 SetOrientation (aPipeLine->GetPlaneOrientation(),
108 aPipeLine->GetRotateX(), aPipeLine->GetRotateY());
109 SetPlanePosition (aPipeLine->GetPlanePosition(),
110 aPipeLine->IsPositionRelative() );
111 SetScaleFactor( aPipeLine->GetScaleFactor() );
112 SetContourPrs( aPipeLine->GetIsContourPrs() );
113 SetNumberOfContours( aPipeLine->GetNumberOfContours() );
118 //----------------------------------------------------------------------------
119 VISU_CutPlanesPL::PlaneOrientation
121 ::GetOrientation(vtkDataSet* theDataSet)
123 theDataSet->Update();
125 vtkFloatingPointType aBounds[6];
126 theDataSet->GetBounds(aBounds);
127 vtkFloatingPointType aDelta[3] = {aBounds[1] - aBounds[0], aBounds[3] - aBounds[2], aBounds[5] - aBounds[4]};
129 if(aDelta[0] >= aDelta[1] && aDelta[0] >= aDelta[2])
130 if(aDelta[1] >= aDelta[2])
131 return VISU_CutPlanesPL::XY;
133 return VISU_CutPlanesPL::ZX;
135 if(aDelta[1] >= aDelta[0] && aDelta[1] >= aDelta[2])
136 if(aDelta[0] >= aDelta[2])
137 return VISU_CutPlanesPL::XY;
139 return VISU_CutPlanesPL::YZ;
141 if(aDelta[2] >= aDelta[0] && aDelta[2] >= aDelta[1])
142 if(aDelta[0] >= aDelta[1])
143 return VISU_CutPlanesPL::ZX;
145 return VISU_CutPlanesPL::YZ;
147 return VISU_CutPlanesPL::XY;
151 //----------------------------------------------------------------------------
154 ::GetScaleFactor(vtkDataSet* theDataSet)
156 theDataSet->Update();
157 vtkFloatingPointType aLength = theDataSet->GetLength(); // diagonal length
159 vtkFloatingPointType aScalarRange[2];
160 GetSourceRange(aScalarRange);
162 static vtkFloatingPointType EPS = 0.3;
163 vtkFloatingPointType aRange = aScalarRange[1];
165 return aLength / aRange * EPS;
171 //----------------------------------------------------------------------------
178 myOrientation = GetOrientation(GetMergedInput());
179 SetScaleFactor(GetScaleFactor(GetMergedInput()));
183 //----------------------------------------------------------------------------
188 return myAppendPolyData->GetOutput();
192 //----------------------------------------------------------------------------
197 vtkFloatingPointType aPlaneNormal[3];
198 vtkFloatingPointType anOrigin[3];
199 GetBasePlane( anOrigin, aPlaneNormal );
201 vtkPolyData* aPolyData = 0;
202 vtkCutter *aCutPlane = 0;
203 vtkDataSet* aDataSet = GetMergedInput();
205 if ( !IsPlanarInput() )
207 aCutPlane = vtkCutter::New();
208 aCutPlane->SetInput(aDataSet);
210 vtkPlane *aPlane = vtkPlane::New();
211 aPlane->SetOrigin(anOrigin);
212 aPlane->SetNormal(aPlaneNormal);
214 aCutPlane->SetCutFunction(aPlane);
217 aPolyData = aCutPlane->GetOutput();
221 if ( !aPolyData || aPolyData->GetNumberOfCells() == 0 ) {
222 myGeometryFilter->SetInput(aDataSet);
223 aPolyData = myGeometryFilter->GetOutput();
226 if ( !myIsContour ) // surface prs
228 if(VISU::IsDataOnCells(aPolyData)) {
229 myCellDataToPointData->SetInput(aPolyData);
230 myCellDataToPointData->PassCellDataOn();
231 myWarpScalar->SetInput(myCellDataToPointData->GetPolyDataOutput());
233 myWarpScalar->SetInput(aPolyData);
237 if(VISU::IsDataOnCells(aPolyData)) {
238 myCellDataToPointData->SetInput(aPolyData);
239 myCellDataToPointData->PassCellDataOn();
240 myContourFilter->SetInput(myCellDataToPointData->GetOutput());
242 myContourFilter->SetInput(aPolyData);
244 vtkFloatingPointType aScalarRange[2];
245 GetSourceRange(aScalarRange);
247 myContourFilter->GenerateValues(GetNumberOfContours(),aScalarRange);
248 myWarpScalar->SetInput(myContourFilter->GetOutput());
251 VISU_CutPlanesPL::ClearAppendPolyData(myAppendPolyData.GetPointer());
252 myAppendPolyData->AddInput(myWarpScalar->GetPolyDataOutput());
257 myWarpScalar->SetNormal(aPlaneNormal);
259 Superclass::Update();
263 //----------------------------------------------------------------------------
268 unsigned long int aSize = Superclass::GetMemorySize();
270 if(vtkDataObject* aDataObject = myGeometryFilter->GetInput())
271 aSize += aDataObject->GetActualMemorySize() * 1024;
273 if(myCellDataToPointData->GetInput())
274 if(vtkDataSet* aDataSet = myCellDataToPointData->GetOutput())
275 aSize += aDataSet->GetActualMemorySize() * 1024;
277 if(vtkDataObject* aDataObject = myContourFilter->GetInput())
278 aSize += aDataObject->GetActualMemorySize() * 1024;
280 if(vtkDataObject* aDataObject = myWarpScalar->GetInput())
281 aSize += aDataObject->GetActualMemorySize() * 1024;
283 int anEnd = myAppendPolyData->GetNumberOfInputConnections(0);
284 for(int anId = 0; anId < anEnd; anId++){
285 if(vtkDataObject* aDataObject = myAppendPolyData->GetInput(anId))
286 aSize += aDataObject->GetActualMemorySize() * 1024;
293 //----------------------------------------------------------------------------
296 ::SetNumberOfContours(int theNumber)
298 myContourFilter->SetNumberOfContours(theNumber);
302 //----------------------------------------------------------------------------
305 ::GetNumberOfContours()
307 return myContourFilter->GetNumberOfContours();
311 //----------------------------------------------------------------------------
314 ::SetScaleFactor(vtkFloatingPointType theScaleFactor)
316 myScaleFactor = theScaleFactor;
317 myWarpScalar->SetScaleFactor(theScaleFactor);
321 //----------------------------------------------------------------------------
326 return myScaleFactor;
330 //----------------------------------------------------------------------------
333 SetContourPrs(bool theIsContourPrs )
335 if(myIsContour == theIsContourPrs)
338 myIsContour = theIsContourPrs;
343 //----------------------------------------------------------------------------
352 //----------------------------------------------------------------------------
355 ::SetPlanePosition(vtkFloatingPointType thePosition,
358 bool anIsSameValue = VISU::CheckIsSameValue(myIsRelative, theIsRelative);
359 anIsSameValue &= (myPosition == thePosition);
363 myIsRelative = theIsRelative;
364 myPosition = thePosition;
369 //----------------------------------------------------------------------------
372 ::IsPositionRelative()
378 //----------------------------------------------------------------------------
379 VISU_CutPlanesPL::PlaneOrientation
381 ::GetPlaneOrientation()
383 return myOrientation;
387 //----------------------------------------------------------------------------
392 switch(myOrientation){
393 case VISU_CutPlanesPL::XY: return myAngle[0];
394 case VISU_CutPlanesPL::YZ: return myAngle[1];
395 case VISU_CutPlanesPL::ZX: return myAngle[2];
401 //----------------------------------------------------------------------------
405 switch(myOrientation){
406 case VISU_CutPlanesPL::XY: return myAngle[1];
407 case VISU_CutPlanesPL::YZ: return myAngle[2];
408 case VISU_CutPlanesPL::ZX: return myAngle[0];
414 //----------------------------------------------------------------------------
417 SetOrientation(VISU_CutPlanesPL::PlaneOrientation theOrientation,
418 vtkFloatingPointType theXAngle,
419 vtkFloatingPointType theYAngle)
421 bool anIsSameValue = VISU::CheckIsSameValue(GetRotateX(), theXAngle);
422 anIsSameValue &= VISU::CheckIsSameValue(GetRotateY(), theYAngle);
423 anIsSameValue &= (myOrientation == theOrientation);
427 switch(theOrientation){
428 case VISU_CutPlanesPL::XY: myAngle[0] = theXAngle; break;
429 case VISU_CutPlanesPL::YZ: myAngle[1] = theXAngle; break;
430 case VISU_CutPlanesPL::ZX: myAngle[2] = theXAngle; break;
433 switch(theOrientation){
434 case VISU_CutPlanesPL::XY: myAngle[1] = theYAngle; break;
435 case VISU_CutPlanesPL::YZ: myAngle[2] = theYAngle; break;
436 case VISU_CutPlanesPL::ZX: myAngle[0] = theYAngle; break;
439 myOrientation = theOrientation;
444 //----------------------------------------------------------------------------
452 //=======================================================================
453 //function : GetBasePlane
455 //=======================================================================
458 ::GetBasePlane(vtkFloatingPointType theOrigin[3],
459 vtkFloatingPointType theNormal[3],
460 bool theCenterOrigine )
462 VISU_CutPlanesPL::GetDir(theNormal,myAngle,myOrientation);
464 vtkFloatingPointType aPosition = myPosition;
465 vtkFloatingPointType aBounds[6], aBoundPrj[3];
468 GetInput()->GetBounds(aBounds);
469 VISU_CutPlanesPL::GetBoundProject(aBoundPrj,aBounds,theNormal);
470 aPosition = aBoundPrj[0] + aBoundPrj[2]*myPosition;
472 VISU::Mul(theNormal,aPosition,theOrigin);
474 if ( theCenterOrigine ) {
475 // move theOrigin to the center of aBounds projections to the plane
476 GetMergedInput()->GetBounds(aBounds);
477 vtkFloatingPointType boundPoints[8][3] = {
478 {aBounds[0],aBounds[2],aBounds[4]},
479 {aBounds[1],aBounds[2],aBounds[4]},
480 {aBounds[0],aBounds[3],aBounds[4]},
481 {aBounds[1],aBounds[3],aBounds[4]},
482 {aBounds[0],aBounds[2],aBounds[5]},
483 {aBounds[1],aBounds[2],aBounds[5]},
484 {aBounds[0],aBounds[3],aBounds[5]},
485 {aBounds[1],aBounds[3],aBounds[5]}};
486 vtkFloatingPointType newOrigin[3] = { 0,0,0 };
487 for(int i = 0; i < 8; i++) {
488 vtkFloatingPointType proj[3];
489 vtkPlane::ProjectPoint( boundPoints[i], theOrigin, theNormal, proj );
490 newOrigin[0] += proj[0];
491 newOrigin[1] += proj[1];
492 newOrigin[2] += proj[2];
494 theOrigin[0] = newOrigin[0] / 8.;
495 theOrigin[1] = newOrigin[1] / 8.;
496 theOrigin[2] = newOrigin[2] / 8.;
500 //=======================================================================
501 //function : GetMinMaxPosition
502 //purpose : return absolute position range
503 //=======================================================================
506 ::GetMinMaxPosition( vtkFloatingPointType& minPos,
507 vtkFloatingPointType& maxPos )
509 vtkFloatingPointType aBounds[6], aBoundPrj[3], aNormal[3];
510 VISU_CutPlanesPL::GetDir(aNormal,myAngle,myOrientation);
511 GetInput()->GetBounds(aBounds);
512 VISU_CutPlanesPL::GetBoundProject(aBoundPrj,aBounds,aNormal);
513 minPos = aBoundPrj[0];
514 maxPos = aBoundPrj[1];
517 //=======================================================================
518 //function : SetMapScale
520 //=======================================================================
524 ::SetMapScale(vtkFloatingPointType theMapScale)
526 Superclass::SetMapScale(theMapScale);
529 vtkFloatingPointType aRange[2];
530 GetSourceRange(aRange);
531 vtkFloatingPointType aNewRange[] = { aRange[1] - theMapScale*(aRange[1]-aRange[0]), aRange[1] };
532 myContourFilter->GenerateValues(GetNumberOfContours(),aNewRange);
534 myWarpScalar->SetScaleFactor(myScaleFactor*theMapScale);