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.opencascade.org/SALOME/ or email : webmaster.salome@opencascade.org
23 // File: VISU_PipeLine.cxx
24 // Author: Alexey PETROV
28 #include "VISU_StreamLinesPL.hxx"
29 #include "VISU_PipeLineUtils.hxx"
30 #include "VISU_UsedPointsFilter.hxx"
31 #include "VTKViewer_GeometryFilter.h"
36 #include <vtkPointSet.h>
37 #include <vtkStreamLine.h>
40 static int MYDEBUG = 0;
42 static int MYDEBUG = 0;
45 static float EPS = 1.0e-7;
46 static float aMinNbOfSteps = 1.0E+2;
47 //static float aMaxNbOfSteps = 1.0E+3;
48 static float aCoeffOfIntStep = 1.0E+1;
51 vtkStandardNewMacro(VISU_StreamLinesPL);
53 VISU_StreamLinesPL::VISU_StreamLinesPL(){
54 myStream = vtkStreamLine::New();
55 myCenters = vtkCellCenters::New();
56 myGeomFilter = VTKViewer_GeometryFilter::New();
57 myPointsFilter = VISU_UsedPointsFilter::New();
62 VISU_StreamLinesPL::~VISU_StreamLinesPL(){
63 myPointsFilter->UnRegisterAllOutputs();
64 myPointsFilter->Delete();
66 myCenters->UnRegisterAllOutputs();
69 myGeomFilter->UnRegisterAllOutputs();
70 myGeomFilter->Delete();
72 myStream->UnRegisterAllOutputs();
76 void VISU_StreamLinesPL::ShallowCopy(VISU_PipeLine *thePipeLine){
77 if(VISU_StreamLinesPL *aPipeLine = dynamic_cast<VISU_StreamLinesPL*>(thePipeLine)){
78 SetParams(aPipeLine->GetIntegrationStep(),
79 aPipeLine->GetPropagationTime(),
80 aPipeLine->GetStepLength(),
81 aPipeLine->GetSource(),
82 aPipeLine->GetUsedPoints(),
83 aPipeLine->GetDirection());
85 VISU_DeformedShapePL::ShallowCopy(thePipeLine);
89 float VISU_StreamLinesPL::GetNecasseryMemorySize(vtkIdType theNbOfPoints, float theStepLength,
90 float thePropogationTime, float thePercents)
92 static float aStreamPointSize = sizeof(float)*15 + sizeof(vtkIdType)*2;
93 static float aStreamArraySize = aStreamPointSize*1024; // == 69632
95 float aNbCells = thePercents*theNbOfPoints*2.0;
96 float aNbPointsPerCell = thePropogationTime/theStepLength;
97 float aCellsSize = aNbCells*(1+aNbPointsPerCell);
98 float aPointsSize = aCellsSize*3.0*sizeof(float);
100 float aConnectivitySize = aCellsSize*sizeof(vtkIdType);
101 float aTypesSize = aNbCells*sizeof(char);
102 float aLocationsSize = aNbCells*sizeof(int);
103 //float aNbCellsPerPoint = aCellsSize / aNbCells - 1;
104 float aMeshSize = aPointsSize + aConnectivitySize + aTypesSize + aLocationsSize;
106 float anAssignedDataSize = aCellsSize*4.0*sizeof(float);
107 float anOutputDataSetSize = aMeshSize + anAssignedDataSize;
109 float aResult = aStreamArraySize*aNbCells + anOutputDataSetSize;
113 int VISU_StreamLinesPL::FindPossibleParams(vtkPointSet* theDataSet, float& theStepLength,
114 float& thePropogationTime, float& thePercents)
116 static float aPercentsDecrease = 3.0, aStepLengthIncrease = 9.0;
117 vtkIdType aNbOfPoints = theDataSet->GetNumberOfPoints();
118 float aSize = GetNecasseryMemorySize(aNbOfPoints,theStepLength,thePropogationTime,thePercents);
119 int isPoss = CheckAvailableMemory(aSize);
121 float aMaxStepLength = max(GetMaxStepLength(theDataSet),thePropogationTime);
122 float aMinStepLength = GetMinStepLength(theDataSet);
123 float aDeltaStepLength = (aMaxStepLength - aMinStepLength)/aStepLengthIncrease;
124 for(int i = 2, aStepChanged = 1, aPerecentsChanged = 1; aStepChanged || aPerecentsChanged; i++){
125 float aStepLength = theStepLength + aDeltaStepLength;
126 if(aStepLength < aMaxStepLength) theStepLength = aStepLength;
127 else if(aStepChanged){
128 aStepLength = aMaxStepLength;
131 float aPercents = thePercents /= aPercentsDecrease;
132 if(aPercents*aNbOfPoints > 1) thePercents = aPercents;
133 else if(aPerecentsChanged) {
134 thePercents = 1.1 / aNbOfPoints;
135 aPerecentsChanged = 0;
137 aSize = GetNecasseryMemorySize(aNbOfPoints,theStepLength,thePropogationTime,thePercents);
138 if(CheckAvailableMemory(aSize)){
144 if(MYDEBUG) MESSAGE("FindPossibleParams - aSize = "<<aSize<<"; isPoss = "<<isPoss);
149 int VISU_StreamLinesPL::SetParams(float theIntStep,
150 float thePropogationTime,
152 vtkPointSet* theSource,
157 vtkPointSet* aDataSet = theSource? theSource: myFieldTransform->GetUnstructuredGridOutput();
159 vtkIdType aNbOfPoints = aDataSet->GetNumberOfPoints();
160 vtkPointSet* aPointSet = myExtractor->GetOutput();
161 if(thePercents*aNbOfPoints < 1) thePercents = 2.0/aNbOfPoints;
162 theIntStep = CorrectIntegrationStep(theIntStep,aPointSet,thePercents);
163 thePropogationTime = CorrectPropagationTime(thePropogationTime,aPointSet);
164 theStepLength = CorrectStepLength(theStepLength,aPointSet);
165 int isAccepted = FindPossibleParams(aPointSet,theStepLength,thePropogationTime,thePercents);
166 if((!isOnlyTry && isAccepted) || (isOnlyTry && isAccepted == 1)){
167 mySource = theSource;
168 myPercents = thePercents;
169 if(GetInput2()->GetCellData()->GetNumberOfArrays()){
170 myCenters->SetInput(aDataSet);
171 myCenters->VertexCellsOn();
172 aDataSet = myCenters->GetOutput();
174 myPointsFilter->SetInput(aDataSet);
175 myPointsFilter->SetPercentsOfUsedPoints(thePercents);
176 aDataSet = myPointsFilter->GetOutput();
177 myStream->SetSource(aDataSet);
178 myStream->SetIntegrationStepLength(theIntStep);
179 myStream->SetMaximumPropagationTime(thePropogationTime);
180 myStream->SetStepLength(theStepLength);
181 myStream->SetSavePointInterval(theIntStep*aMinNbOfSteps);
182 myStream->SetIntegrationDirection(theDirection);
183 myStream->Modified();
190 vtkPointSet* VISU_StreamLinesPL::GetSource() {
194 float VISU_StreamLinesPL::GetUsedPoints() {
198 vtkDataSet* VISU_StreamLinesPL::GetStreamerSource(){
199 return myStream->GetSource();
202 float VISU_StreamLinesPL::GetVelocityCoeff(){
203 return GetVelocityCoeff(myExtractor->GetOutput());
206 float VISU_StreamLinesPL::GetVelocityCoeff(vtkPointSet* theDataSet){
207 float* aScalarRange = theDataSet->GetScalarRange();
208 return (fabs(aScalarRange[1]) + fabs(aScalarRange[0]))/2.0;
212 int VISU_StreamLinesPL::IsPossible(vtkPointSet* theDataSet, float thePercents){
213 float aStepLength = GetBaseStepLength(theDataSet);
214 float aBasePropTime = GetBasePropagationTime(theDataSet);
215 VISU_UsedPointsFilter *aPointsFilter = VISU_UsedPointsFilter::New();
216 aPointsFilter->SetInput(theDataSet);
217 vtkPointSet* aDataSet = aPointsFilter->GetOutput();
219 int aRes = FindPossibleParams(aDataSet,aStepLength,aBasePropTime,thePercents);
220 aPointsFilter->UnRegisterAllOutputs();
221 aPointsFilter->Delete();
226 float VISU_StreamLinesPL::GetIntegrationStep(){
227 return myStream->GetIntegrationStepLength();
229 float VISU_StreamLinesPL::GetStepLength() {
230 return myStream->GetStepLength();
232 float VISU_StreamLinesPL::GetPropagationTime() {
233 return myStream->GetMaximumPropagationTime();
235 int VISU_StreamLinesPL::GetDirection(){
236 return myStream->GetIntegrationDirection();
240 float VISU_StreamLinesPL::GetMinIntegrationStep(vtkPointSet* theDataSet, float thePercents) {
241 if(!theDataSet) return -1.0;
244 theDataSet->Update();
245 float* aBounds = theDataSet->GetBounds();
246 for(int i = 0, j = 0; i < 3; ++i, j = 2*i){
247 float tmp = aBounds[j+1] - aBounds[j];
253 if (degree < 1) return 0.0; // absolutely empty object
254 float anStepLength = GetMaxIntegrationStep(theDataSet)/aCoeffOfIntStep;
255 float aBasePropTime = GetBasePropagationTime(theDataSet)/GetVelocityCoeff(theDataSet);
257 vtkIdType aNbOfPoints = theDataSet->GetNumberOfPoints();
258 float aSize = GetNecasseryMemorySize(aNbOfPoints,anStepLength,aBasePropTime,thePercents);
259 float aRealSize = GetAvailableMemory(aSize);
260 float anAverageVolume = aVolume / aRealSize;
261 float aStep = pow(double(anAverageVolume), double(1.0/double(degree)));
264 float VISU_StreamLinesPL::GetMinIntegrationStep(){
265 return GetMinIntegrationStep(myExtractor->GetOutput(),GetUsedPoints());
269 float VISU_StreamLinesPL::GetMaxIntegrationStep(vtkPointSet* theDataSet) {
270 if(!theDataSet) return -1.0;
271 theDataSet->Update();
272 float aLength = theDataSet->GetLength();
273 float* aBounds = theDataSet->GetBounds();
274 float aMaxSizeY = (aBounds[3]-aBounds[2])/aLength;
275 float aMaxSizeZ = (aBounds[5]-aBounds[4])/aLength;
276 float aMinMax = (aBounds[1] - aBounds[0])/aLength;
277 if (aMinMax < EPS || (aMaxSizeY < aMinMax && aMaxSizeY > EPS)) aMinMax = aMaxSizeY;
278 if (aMinMax < EPS || (aMaxSizeZ < aMinMax && aMaxSizeZ > EPS)) aMinMax = aMaxSizeZ;
279 return aMinMax*aLength/2.0;
281 float VISU_StreamLinesPL::GetMaxIntegrationStep(){
282 return GetMaxIntegrationStep(myExtractor->GetOutput());
285 float VISU_StreamLinesPL::GetBaseIntegrationStep(vtkPointSet* theDataSet, float thePercents){
286 theDataSet->Update();
287 float aMinIntegrationStep = GetMinIntegrationStep(theDataSet,thePercents);
288 float aMaxIntegrationStep = GetMaxIntegrationStep(theDataSet);
289 float anIntegrationStep = aMaxIntegrationStep / aCoeffOfIntStep;
290 float aMinMax = theDataSet->GetLength()/theDataSet->GetNumberOfPoints();
291 if(aMinMax > anIntegrationStep)
292 anIntegrationStep = (anIntegrationStep*aCoeffOfIntStep*0.9+aMinMax)/aCoeffOfIntStep;
293 if(aMinIntegrationStep > anIntegrationStep)
294 anIntegrationStep = aMinIntegrationStep;
295 return anIntegrationStep;
298 float VISU_StreamLinesPL::CorrectIntegrationStep(float theStep, vtkPointSet* theDataSet, float thePercents){
299 theDataSet->Update();
300 float aMinIntegrationStep = GetMinIntegrationStep(theDataSet,thePercents);
301 float aMaxIntegrationStep = GetMaxIntegrationStep(theDataSet);
302 if(aMinIntegrationStep > theStep)
303 theStep = aMinIntegrationStep;
304 if(aMaxIntegrationStep < theStep)
305 theStep = aMaxIntegrationStep;
310 float VISU_StreamLinesPL::GetMinPropagationTime(vtkPointSet* theDataSet){
311 if(!theDataSet) return -1.0;
312 return GetMinStepLength(theDataSet);
314 float VISU_StreamLinesPL::GetMinPropagationTime(){
315 return GetMinPropagationTime(myExtractor->GetOutput());
318 float VISU_StreamLinesPL::GetMaxPropagationTime(vtkPointSet* theDataSet){
319 if(!theDataSet) return -1.0;
320 return GetBasePropagationTime(theDataSet)*aMinNbOfSteps;
322 float VISU_StreamLinesPL::GetMaxPropagationTime(){
323 return GetMaxPropagationTime(myExtractor->GetOutput());
326 float VISU_StreamLinesPL::CorrectPropagationTime(float thePropagationTime, vtkPointSet* theDataSet){
327 float aMinPropagationTime = GetMinPropagationTime(theDataSet);
328 float aMaxPropagationTime = GetMaxPropagationTime(theDataSet);
329 if(aMinPropagationTime > thePropagationTime)
330 thePropagationTime = aMinPropagationTime;
331 if(aMaxPropagationTime < thePropagationTime)
332 thePropagationTime = aMaxPropagationTime;
333 return thePropagationTime;
335 float VISU_StreamLinesPL::GetBasePropagationTime(vtkPointSet* theDataSet){
336 if(!theDataSet) return -1.0;
337 theDataSet->Update();
338 float aPropagationTime = theDataSet->GetLength()/GetVelocityCoeff(theDataSet);
339 return aPropagationTime;
341 float VISU_StreamLinesPL::GetBasePropagationTime(){
342 return GetBasePropagationTime(myExtractor->GetOutput());
346 float VISU_StreamLinesPL::GetMinStepLength(vtkPointSet* theDataSet){
347 static float aNbOfStepsOfIntStep = 1.0E+1;
348 float anIntStep = GetMinIntegrationStep(theDataSet);
349 float aStepLength = anIntStep*aNbOfStepsOfIntStep/GetVelocityCoeff(theDataSet);
352 float VISU_StreamLinesPL::GetMinStepLength(){
353 return GetMinStepLength(myExtractor->GetOutput());
356 float VISU_StreamLinesPL::GetMaxStepLength(vtkPointSet* theDataSet){
357 float aStepLength = GetBasePropagationTime(theDataSet);
360 float VISU_StreamLinesPL::GetMaxStepLength(){
361 return GetMaxStepLength(myExtractor->GetOutput());
364 float VISU_StreamLinesPL::CorrectStepLength(float theStep, vtkPointSet* theDataSet){
365 float aMinStep = GetMinStepLength(theDataSet);
366 if(theStep < aMinStep) theStep = aMinStep;
367 float aMaxStep = GetMaxStepLength(theDataSet);
368 if(theStep > aMaxStep) theStep = aMaxStep;
371 float VISU_StreamLinesPL::GetBaseStepLength(vtkPointSet* theDataSet){
372 static float anAvgNbOfSteps = 1.0E+2;
373 float aPropagationTime = GetBasePropagationTime(theDataSet);
374 float aStepLength = aPropagationTime/anAvgNbOfSteps;
375 aStepLength = CorrectStepLength(aStepLength,theDataSet);
380 void VISU_StreamLinesPL::Init(){
381 VISU_ScalarMapPL::Init();
382 vtkPointSet* aDataSet = myExtractor->GetOutput();
383 float anIntStep = GetBaseIntegrationStep(aDataSet);
384 float aPropagationTime = GetBasePropagationTime(aDataSet);
385 float aStepLength = GetBaseStepLength(aDataSet);
386 SetParams(anIntStep,aPropagationTime,aStepLength);
389 VISU_ScalarMapPL::THook* VISU_StreamLinesPL::DoHook(){
390 GetInput2()->Update();
391 VISU::CellDataToPoint(myStream,myCellDataToPointData,GetInput2(),myFieldTransform);
392 float *aBounds = GetInput2()->GetBounds();
393 myGeomFilter->SetExtent(aBounds);
394 myGeomFilter->ExtentClippingOn();
395 myGeomFilter->SetInput(myStream->GetOutput());
396 return myGeomFilter->GetOutput();
399 void VISU_StreamLinesPL::Update(){
400 VISU_ScalarMapPL::Update();
403 void VISU_StreamLinesPL::SetMapScale(float theMapScale){
404 VISU_ScalarMapPL::SetMapScale(theMapScale);