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_Convertor_impl.cxx
24 // Author : Alexey PETROV
27 #include "VISU_Convertor_impl.hxx"
28 #include "VISU_ConvertorUtils.hxx"
30 #include <vtkIdList.h>
31 #include <vtkCellType.h>
32 #include <vtkIntArray.h>
33 #include <vtkCellArray.h>
34 #include <vtkFloatArray.h>
35 #include <vtkUnsignedCharArray.h>
36 #include <vtkPointData.h>
37 #include <vtkCellData.h>
38 #include <vtkCellLinks.h>
40 #include <vtkMergeDataObjectFilter.h>
43 #include <qfileinfo.h>
51 static float ERR_SIZE_CALC = 1.00;
53 static int MYVTKDEBUG = 0;
56 static int MYDEBUG = 0;
57 static int MYDEBUGWITHFILES = 0;
59 static int MYDEBUG = 0;
60 static int MYDEBUGWITHFILES = 0;
67 std::string dtos(const std::string& fmt, T val){
68 static QString aString;
69 aString.sprintf(fmt.c_str(),val);
70 return aString.latin1();
73 enum ECoordName{eX, eY, eZ, eNone};
74 typedef VISU::TCoord (*TGetCoord)(const VISU::TMeshImpl::TPointsCoord&, int);
76 template<ECoordName TheCoordId>
78 GetCoord(const VISU::TMeshImpl::TPointsCoord& thePointsCoord,
81 return thePointsCoord[theStartPos+TheCoordId];
86 GetCoord<eNone>(const VISU::TMeshImpl::TPointsCoord& thePointsCoord,
93 TGetCoord aXYZGetCoord[3] = {
100 TGetCoord aXYGetCoord[3] = {
106 TGetCoord aYZGetCoord[3] = {
112 TGetCoord aXZGetCoord[3] = {
119 TGetCoord aXGetCoord[3] = {
125 TGetCoord aYGetCoord[3] = {
131 TGetCoord aZGetCoord[3] = {
139 const VISU::TMeshImpl::TPointsCoord& myPointsCoord;
140 TGetCoord* myGetCoord;
142 TCoordHelper(const VISU::TMeshImpl::TPointsCoord& thePointsCoord,
143 TGetCoord* theGetCoord):
144 myPointsCoord(thePointsCoord),
145 myGetCoord(theGetCoord)
147 virtual ~TCoordHelper(){}
149 GetCoord(int theStartPos, int theCoodId)
151 return (*myGetCoord[theCoodId])(myPointsCoord,theStartPos);
154 typedef std::auto_ptr<TCoordHelper> TCoordHelperPtr;
156 void GetPoints(VISU::TVTKSource& theStorage, VISU::PMeshImpl theMesh)
158 vtkPoints* aPoints = theMesh->myPoints.GetPointer();
160 aPoints = vtkPoints::New();
161 TCoordHelperPtr aCoordHelperPtr;
162 const VISU::TMeshImpl::TPointsCoord& anArray = theMesh->myPointsCoord;
164 int aMeshDimension = theMesh->myDim;
165 bool anIsDimPresent[3] = {false, false, false};
166 for(int iDim = 0; iDim < aMeshDimension; iDim++){
167 string aDimName = theMesh->myPointsDim[iDim];
168 if(aDimName == "x" || aDimName == "X")
169 anIsDimPresent[eX] = true;
170 else if(aDimName == "y" || aDimName == "Y")
171 anIsDimPresent[eY] = true;
172 else if(aDimName == "z" || aDimName == "Z")
173 anIsDimPresent[eZ] = true;
176 switch(aMeshDimension){
178 aCoordHelperPtr.reset(new TCoordHelper(anArray,aXYZGetCoord));
181 if(anIsDimPresent[eY] && anIsDimPresent[eZ])
182 aCoordHelperPtr.reset(new TCoordHelper(anArray,aYZGetCoord));
183 else if(anIsDimPresent[eX] && anIsDimPresent[eZ])
184 aCoordHelperPtr.reset(new TCoordHelper(anArray,aXZGetCoord));
186 aCoordHelperPtr.reset(new TCoordHelper(anArray,aXYGetCoord));
189 if(anIsDimPresent[eY])
190 aCoordHelperPtr.reset(new TCoordHelper(anArray,aYGetCoord));
191 else if(anIsDimPresent[eZ])
192 aCoordHelperPtr.reset(new TCoordHelper(anArray,aZGetCoord));
194 aCoordHelperPtr.reset(new TCoordHelper(anArray,aXGetCoord));
199 if(MYVTKDEBUG) aPoints->DebugOn();
200 vtkIdType iEnd = theMesh->myPointsCoord.size();
201 vtkIdType aNbPoints = iEnd / theMesh->myDim;
202 aPoints->SetNumberOfPoints(aNbPoints);
203 MSG(MYDEBUG,"GetPoints - aNbPoints = "<<aNbPoints<<"; myDim = "<<theMesh->myDim);
204 for (vtkIdType i = 0, j = 0; i < iEnd; i += theMesh->myDim, j++)
206 aCoordHelperPtr->GetCoord(i,eX),
207 aCoordHelperPtr->GetCoord(i,eY),
208 aCoordHelperPtr->GetCoord(i,eZ));
209 theMesh->myPoints = aPoints;
212 theStorage->SetPoints(aPoints);
216 inline void PrintCells(int& theStartId,
217 vtkCellArray* theConnectivity,
218 const VISU::TMeshOnEntityImpl::TConnect& theVector)
220 vtkIdList *anIdList = vtkIdList::New();
221 int kEnd = theVector.size();
222 anIdList->SetNumberOfIds(kEnd);
223 for(int k = 0; k < kEnd; k++)
224 anIdList->SetId(k,theVector[k]);
225 theConnectivity->InsertNextCell(anIdList);
229 void GetCellsOnEntity(VISU::TVTKSource& theStorage,
230 const VISU::PMeshOnEntityImpl theMeshOnEntity,
231 const string& theFamilyName)
233 //Check on existing family
234 PFamilyImpl aFamily = GetFamily(theMeshOnEntity,theFamilyName);
236 pair<int,int> aCellsDim = theMeshOnEntity->GetCellsDims(theFamilyName);
237 int aNbCells = aCellsDim.first, aCellsSize = aCellsDim.second;
238 vtkCellArray* aConnectivity = vtkCellArray::New();
239 aConnectivity->Allocate(aCellsSize,0);
240 vtkUnsignedCharArray* aCellTypesArray = vtkUnsignedCharArray::New();
241 aCellTypesArray->SetNumberOfComponents(1);
242 aCellTypesArray->SetNumberOfTuples(aNbCells);
243 MSG(MYDEBUG,"GetCellsOnEntity - isFamilyPresent = "<<bool(aFamily));
244 const VISU::TMeshOnEntityImpl::TCellsConn &aCellsConn = theMeshOnEntity->myCellsConn;
245 VISU::TMeshOnEntityImpl::TCellsConn::const_iterator aCellsConnIter = aCellsConn.begin();
246 for(int i = 0, j = 0; aCellsConnIter != aCellsConn.end(); aCellsConnIter++){
247 const VISU::TMeshOnEntityImpl::TConnForCellType& anArray = aCellsConnIter->second;
248 int aVtkType = aCellsConnIter->first;
249 MSG(MYDEBUG,"GetCellsOnEntity - aVtkType = "<<aVtkType<<"; anArray.size() = "<<anArray.size());
251 for(int k = 0, kEnd = anArray.size(); k < kEnd; k++, i++){
252 PrintCells(i,aConnectivity,anArray[k]);
253 aCellTypesArray->SetValue(j++,(unsigned char)aVtkType);
256 const VISU::TFamilyImpl::TSubMesh& aSubMesh = aFamily->mySubMesh;
258 EXCEPTION(runtime_error,"GetCells >> There is no elements on the family !!!");
259 VISU::TFamilyImpl::TSubMesh::const_iterator aSubMeshIter = aSubMesh.find(aVtkType);
260 if(aSubMeshIter == aSubMesh.end()) continue;
261 const VISU::TFamilyImpl::TSubMeshOnCellType& aSubMeshOnCellType = aSubMeshIter->second;
262 MSG(MYDEBUG,"GetCellsOnEntity - aSubMeshOnCellType.size() = "<<aSubMeshOnCellType.size());
263 VISU::TFamilyImpl::TSubMeshOnCellType::const_iterator aSubMeshOnCellTypeIter = aSubMeshOnCellType.begin();
264 for(; aSubMeshOnCellTypeIter != aSubMeshOnCellType.end(); aSubMeshOnCellTypeIter++, i++){
265 PrintCells(i,aConnectivity,anArray[*aSubMeshOnCellTypeIter]);
266 aCellTypesArray->SetValue(j++,(unsigned char)aVtkType);
270 vtkIdType *pts = 0, npts = 0;
271 vtkIntArray* aCellLocationsArray = vtkIntArray::New();
272 aCellLocationsArray->SetNumberOfComponents(1);
273 aCellLocationsArray->SetNumberOfTuples(aNbCells);
274 aConnectivity->InitTraversal();
275 for(int i=0; aConnectivity->GetNextCell(npts,pts); i++)
276 aCellLocationsArray->SetValue(i,aConnectivity->GetTraversalLocation(npts));
277 theStorage->SetCells(aCellTypesArray,aCellLocationsArray,aConnectivity);
278 if(MYVTKDEBUG) aConnectivity->DebugOn();
279 aCellLocationsArray->Delete();
280 aCellTypesArray->Delete();
281 aConnectivity->Delete();
285 void GetCellsOnGroup(VISU::TVTKSource& theStorage,
286 VISU::PMeshImpl theMesh,
287 const VISU::TFamilyAndEntitySet& theFamilyAndEntitySet)
289 //Calculate dimentions of the group
290 int aNbCells = 0, aCellsSize = 0;
291 VISU::TFamilyAndEntitySet::const_iterator aFamilyAndEntitySetIter = theFamilyAndEntitySet.begin();
292 for(; aFamilyAndEntitySetIter != theFamilyAndEntitySet.end(); aFamilyAndEntitySetIter++){
293 const VISU::TFamilyAndEntity& aFamilyAndEntity = *aFamilyAndEntitySetIter;
294 const string& aFamilyName = aFamilyAndEntity.first;
295 const VISU::TEntity& anEntity = aFamilyAndEntity.second;
296 const VISU::PMeshOnEntity aMeshOnEntity = theMesh->myMeshOnEntityMap[anEntity];
297 pair<int,int> aCellsDim = aMeshOnEntity->GetCellsDims(aFamilyName);
298 aNbCells += aCellsDim.first;
299 aCellsSize += aCellsDim.second;
301 vtkCellArray* aConnectivity = vtkCellArray::New();
302 aConnectivity->Allocate(aCellsSize,0);
303 vtkUnsignedCharArray* aCellTypesArray = vtkUnsignedCharArray::New();
304 aCellTypesArray->SetNumberOfComponents(1);
305 aCellTypesArray->SetNumberOfTuples(aNbCells);
306 aFamilyAndEntitySetIter = theFamilyAndEntitySet.begin();
307 for(int i = 0, j = 0; aFamilyAndEntitySetIter != theFamilyAndEntitySet.end(); aFamilyAndEntitySetIter++){
308 const VISU::TFamilyAndEntity& aFamilyAndEntity = *aFamilyAndEntitySetIter;
309 const string& aFamilyName = aFamilyAndEntity.first;
310 const VISU::TEntity& anEntity = aFamilyAndEntity.second;
311 PMeshOnEntityImpl aMeshOnEntity = theMesh->myMeshOnEntityMap[anEntity];
312 PFamilyImpl aFamily = GetFamily(aMeshOnEntity,aFamilyName);
313 const VISU::TMeshOnEntityImpl::TCellsConn &aCellsConn = aMeshOnEntity->myCellsConn;
314 VISU::TMeshOnEntityImpl::TCellsConn::const_iterator aCellsConnIter = aCellsConn.begin();
315 for(; aCellsConnIter != aCellsConn.end(); aCellsConnIter++){
316 const VISU::TMeshOnEntityImpl::TConnForCellType& anArray = aCellsConnIter->second;
317 int aVtkType = aCellsConnIter->first;
318 MSG(MYDEBUG,"GetCellsOnGroup - aVtkType = "<<aVtkType<<"; anArray.size() = "<<anArray.size());
319 const VISU::TFamilyImpl::TSubMesh& aSubMesh = aFamily->mySubMesh;
321 EXCEPTION(runtime_error,"GetCells >> There is no elements on the family !!!");
322 VISU::TFamilyImpl::TSubMesh::const_iterator aSubMeshIter = aSubMesh.find(aVtkType);
323 if(aSubMeshIter == aSubMesh.end()) continue;
324 const VISU::TFamilyImpl::TSubMeshOnCellType& aSubMeshOnCellType = aSubMeshIter->second;
325 MSG(MYDEBUG,"GetCellsOnGroup - aSubMeshOnCellType.size() = "<<aSubMeshOnCellType.size());
326 VISU::TFamilyImpl::TSubMeshOnCellType::const_iterator aSubMeshOnCellTypeIter = aSubMeshOnCellType.begin();
327 for(; aSubMeshOnCellTypeIter != aSubMeshOnCellType.end(); aSubMeshOnCellTypeIter++, i++){
328 PrintCells(i,aConnectivity,anArray[*aSubMeshOnCellTypeIter]);
329 aCellTypesArray->SetValue(j++,(unsigned char)aVtkType);
333 vtkIdType *pts = 0, npts = 0;
334 vtkIntArray* aCellLocationsArray = vtkIntArray::New();
335 aCellLocationsArray->SetNumberOfComponents(1);
336 aCellLocationsArray->SetNumberOfTuples(aNbCells);
337 aConnectivity->InitTraversal();
338 for(int i = 0; aConnectivity->GetNextCell(npts,pts); i++)
339 aCellLocationsArray->SetValue(i,aConnectivity->GetTraversalLocation(npts));
340 theStorage->SetCells(aCellTypesArray,aCellLocationsArray,aConnectivity);
341 aCellLocationsArray->Delete();
342 aCellTypesArray->Delete();
343 aConnectivity->Delete();
347 void InitProfile(VISU::TVTKExtractFilter& theFilter,
348 PMeshOnEntityImpl theMeshOnEntity,
349 PValForTimeImpl theValForTime)
351 const VISU::TValForTimeImpl::TValForCells& aValForCells = theValForTime->myValForCells;
352 const VISU::TMeshOnEntityImpl::TCellsConn &aCellsConn = theMeshOnEntity->myCellsConn;
353 VISU::TMeshOnEntityImpl::TCellsConn::const_iterator aCellsConnIter = aCellsConn.begin();
354 for(; aCellsConnIter != aCellsConn.end(); aCellsConnIter++){
355 const vtkIdType& aCellType = aCellsConnIter->first;
356 if(aValForCells.find(aCellType) == aValForCells.end())
357 theFilter->RemoveCellsWithType(aCellType);
362 void GetValsOnTimeStamp(vtkFloatArray *theFloatArray,
363 const vtkIdType& theNumberOfTuples,
364 const std::string& theFieldName,
365 VISU::PFieldImpl theField,
366 VISU::PValForTimeImpl theValForTime)
368 //theFloatArray->DebugOn();
369 theFloatArray->SetNumberOfTuples(theNumberOfTuples);
370 theFloatArray->SetName(theFieldName.c_str());
371 MSG(MYDEBUG,"GetValsOnTimeStamp - theNumberOfTuples = "<<theNumberOfTuples);
372 const VISU::TValForTimeImpl::TValForCells& aValForCells = theValForTime->myValForCells;
373 VISU::TValForTimeImpl::TValForCells::const_iterator aValForCellsIter = aValForCells.begin();
374 for(int k = 0; aValForCellsIter != aValForCells.end(); aValForCellsIter++) {
375 const VISU::TValForTimeImpl::TValForCellsWithType& anArray = aValForCellsIter->second;
376 int iEnd = anArray.size()/theField->myNbComp;
377 int aVtkType = aValForCellsIter->first;
378 MSG(MYDEBUG,"GetValsOnTimeStamp - iEnd = "<<iEnd<<"; aVtkType = "<<aVtkType);
379 switch(theField->myNbComp) {
381 for (int i = 0; i < iEnd; i++)
382 theFloatArray->SetTuple1(k++,anArray[i]);
385 for (int i = 0, ji = 0; i < iEnd; ++i, ji = i*2)
386 theFloatArray->SetTuple3(k++,anArray[ji],anArray[ji+1],0.0);
389 for (int i = 0, ji = 0; i < iEnd; ++i, ji = i*3)
390 theFloatArray->SetTuple3(k++,anArray[ji],anArray[ji+1],anArray[ji+2]);
393 for (int i = 0, ji = 0; i < iEnd; ++i, ji = i*4)
394 theFloatArray->SetTuple3(k++,anArray[ji],anArray[ji+1],0.0);
397 for (int i = 0, ji = 0; i < iEnd; ++i, ji = i*6)
398 theFloatArray->SetTuple3(k++,anArray[ji],anArray[ji+1],anArray[ji+2]);
401 EXCEPTION(runtime_error,"GetValsOnTimeStamp - There is no an algorithm for representation of the field !!!");
406 string GenerateFieldName(const VISU::PFieldImpl theField,
407 const VISU::PValForTimeImpl theValForTime)
409 const VISU::TTime& aTime = theValForTime->myTime;
410 string aFieldName = theField->myMeshName + ", " + theField->myName + ": " +
411 VISU_Convertor::GenerateName(aTime);
415 void GetTimeStamp(VISU::TVTKSource& theStorage,
416 const VISU::PFieldImpl theField,
417 const VISU::PValForTimeImpl theValForTime)
419 int aNumberOfTuples = theField->myDataSize/theField->myNbComp;
420 string aFieldName = GenerateFieldName(theField,theValForTime);
421 MSG(MYDEBUG,"GetTimeStamp(TVTKSource) - aFieldName = "<<aFieldName<<
422 "; aNumberOfTuples = "<<aNumberOfTuples);
424 vtkDataSetAttributes* aDataSetAttributes;
425 switch(theField->myEntity){
426 case VISU::NODE_ENTITY :
427 aDataSetAttributes = theStorage->GetPointData();
430 aDataSetAttributes = theStorage->GetCellData();
433 vtkFloatArray *aFloatArray = vtkFloatArray::New();
434 switch(theField->myNbComp) {
436 aFloatArray->SetNumberOfComponents(1);
437 aDataSetAttributes->SetScalars(aFloatArray);
440 aFloatArray->SetNumberOfComponents(3);
441 aDataSetAttributes->SetVectors(aFloatArray);
444 GetValsOnTimeStamp(aFloatArray,aNumberOfTuples,aFieldName,theField,theValForTime);
447 void GetTimeStamp(VISU::TVTKAttribyteFilter& theAttribyteFilter,
448 VISU::TVTKMergetFilter& theMergeFilter,
449 VISU::TVTKExtractFilter& theExtractFilter,
450 const VISU::PFieldImpl theField,
451 const VISU::PValForTimeImpl theValForTime)
453 int aNumberOfTuples = theField->myDataSize/theField->myNbComp;
454 string aFieldName = GenerateFieldName(theField,theValForTime);
455 MSG(MYDEBUG,"GetTimeStamp(TVTKAttribyteFilter) - aFieldName = "<<aFieldName<<
456 "; aNumberOfTuples = "<<aNumberOfTuples);
458 vtkDataObject* aDataObject = vtkDataObject::New();
459 theMergeFilter->SetDataObject(aDataObject);
460 aDataObject->Delete();
462 theMergeFilter->SetInput(theExtractFilter->GetOutput());
463 theAttribyteFilter->SetInput(theMergeFilter->GetOutput());
465 switch(theField->myEntity){
466 case VISU::NODE_ENTITY :
467 theMergeFilter->SetOutputFieldToPointDataField();
468 theAttribyteFilter->SetInputFieldToPointDataField();
469 theAttribyteFilter->SetOutputAttributeDataToPointData();
472 theMergeFilter->SetOutputFieldToCellDataField();
473 theAttribyteFilter->SetInputFieldToCellDataField();
474 theAttribyteFilter->SetOutputAttributeDataToCellData();
477 vtkFloatArray *aFloatArray = vtkFloatArray::New();
478 switch(theField->myNbComp) {
480 aFloatArray->SetNumberOfComponents(1);
481 theAttribyteFilter->SetScalarComponent(0,aFieldName.c_str(),0);
484 aFloatArray->SetNumberOfComponents(3);
485 theAttribyteFilter->SetVectorComponent(0,aFieldName.c_str(),0);
486 theAttribyteFilter->SetVectorComponent(1,aFieldName.c_str(),1);
487 theAttribyteFilter->SetVectorComponent(2,aFieldName.c_str(),2);
490 vtkFieldData* aFieldData = aDataObject->GetFieldData();
491 aFieldData->AddArray(aFloatArray);
492 aFloatArray->Delete();
494 GetValsOnTimeStamp(aFloatArray,aNumberOfTuples,aFieldName,theField,theValForTime);
498 VISU_Convertor_impl::VISU_Convertor_impl() {
502 VISU_Convertor_impl::~VISU_Convertor_impl() {}
504 VISU_Convertor::TOutput*
505 VISU_Convertor_impl::GetMeshOnEntity(const string& theMeshName,
506 const VISU::TEntity& theEntity,
507 const string& theFamilyName)
509 MSG(MYDEBUG,"GetMeshOnEntity - theMeshName = '"<<theMeshName<<
510 "'; theEntity = "<<theEntity<<"; theFamilyName = '"<<theFamilyName<<"'");
511 //Cheching possibility do the query
512 TFindMeshOnEntity aFindMeshOnEntity =
513 FindMeshOnEntity(theMeshName,theEntity,theFamilyName);
514 VISU::TVTKSource* pSource;
515 PMeshImpl aMesh = boost::get<0>(aFindMeshOnEntity);;
516 PFamilyImpl aFamily = boost::get<2>(aFindMeshOnEntity);
517 PMeshOnEntityImpl aMeshOnEntity = boost::get<1>(aFindMeshOnEntity);
519 pSource = &(aFamily->myStorage);
521 pSource = &(aMeshOnEntity->myStorage);
522 VISU::TVTKSource& aSource = *pSource;
525 if(aSource.GetPointer() == NULL){
526 aSource = TOutput::New();
528 if(MYVTKDEBUG) aSource->DebugOn();
529 LoadMeshOnEntity(aMeshOnEntity,theFamilyName);
530 GetPoints(aSource,aMesh);
531 GetCellsOnEntity(aSource,aMeshOnEntity,theFamilyName);
532 if(MYDEBUGWITHFILES){
533 string aMeshName = QString(theMeshName.c_str()).simplifyWhiteSpace().latin1();
534 string aFamilyName = QString(theFamilyName.c_str()).simplifyWhiteSpace().latin1();
535 string aFileName = string("/users/")+getenv("USER")+"/"+getenv("USER")+"-";
536 aFileName += aMeshName + dtos("-%d-",int(theEntity)) + aFamilyName + "-Conv.vtk";
537 VISU::WriteToFile(aSource.GetPointer(),aFileName);
541 GetMeshOnEntitySize(theMeshName,theEntity,theFamilyName);
542 vtkUnstructuredGrid* aDataSet = aSource.GetPointer();
544 MSG(MYVTKDEBUG,"GetMeshOnEntity - GetPoints() = "<<float(aDataSet->GetPoints()->GetActualMemorySize()*1000));
545 MSG(MYVTKDEBUG,"GetMeshOnEntity - GetCells() = "<<float(aDataSet->GetCells()->GetActualMemorySize()*1000));
546 MSG(MYVTKDEBUG,"GetMeshOnEntity - GetCellTypesArray() = "<<float(aDataSet->GetCellTypesArray()->GetActualMemorySize()*1000));
547 MSG(MYVTKDEBUG,"GetMeshOnEntity - GetCellLocationsArray() = "<<float(aDataSet->GetCellLocationsArray()->GetActualMemorySize()*1000));
548 aDataSet->BuildLinks();
549 MSG(MYVTKDEBUG,"GetMeshOnEntity - GetCellLinks() = "<<float(aDataSet->GetCellLinks()->GetActualMemorySize()*1000));
550 MSG(MYVTKDEBUG,"GetMeshOnEntity - GetActualMemorySize() = "<<float(aDataSet->GetActualMemorySize()*1000));
553 aSource = vtkSmartPointerBase();
556 return aSource.GetPointer();
559 VISU_Convertor::TOutput*
560 VISU_Convertor_impl::GetMeshOnGroup(const string& theMeshName,
561 const string& theGroupName)
563 MSG(MYDEBUG,"GetMeshOnGroup - theMeshName = '"<<theMeshName<<
564 "'; theGroupName = '"<<theGroupName<<"'");
565 //Cheching possibility do the query
566 TFindMeshOnGroup aFindMeshOnGroup = FindMeshOnGroup(theMeshName,theGroupName);
567 PMeshImpl aMesh = boost::get<0>(aFindMeshOnGroup);
568 PGroupImpl aGroup = boost::get<1>(aFindMeshOnGroup);
569 const VISU::TFamilyAndEntitySet& aFamilyAndEntitySet = aGroup->myFamilyAndEntitySet;
570 VISU::TVTKSource& aSource = aGroup->myStorage;
573 if(aSource.GetPointer() == NULL){
574 aSource = TOutput::New();
576 LoadMeshOnGroup(aMesh,aFamilyAndEntitySet);
577 GetPoints(aSource,aMesh);
578 GetCellsOnGroup(aSource,aMesh,aFamilyAndEntitySet);
579 if(MYDEBUGWITHFILES){
580 string aMeshName = QString(theMeshName.c_str()).simplifyWhiteSpace().latin1();
581 string aGroupName = QString(theGroupName.c_str()).simplifyWhiteSpace().latin1();
582 string aFileName = string("/users/")+getenv("USER")+"/"+getenv("USER")+"-";
583 aFileName += aMeshName + "-" + aGroupName + "-Conv.vtk";
584 VISU::WriteToFile(aSource.GetPointer(),aFileName);
588 aSource = vtkSmartPointerBase();
591 return aSource.GetPointer();
594 VISU_Convertor::TOutput*
595 VISU_Convertor_impl::GetTimeStampOnMesh(const string& theMeshName,
596 const VISU::TEntity& theEntity,
597 const string& theFieldName,
600 MSG(MYDEBUG,"GetTimeStampOnMesh - theMeshName = '"<<theMeshName<<"; theEntity = "<<theEntity);
601 MSG(MYDEBUG,"GetTimeStampOnMesh - theFieldName = '"<<theFieldName<<"'; theStampsNum = "<<theStampsNum);
603 //Cheching possibility do the query
604 TFindTimeStamp aFindTimeStamp =
605 FindTimeStamp(theMeshName,theEntity,theFieldName,theStampsNum);
606 PMeshImpl aMesh = boost::get<0>(aFindTimeStamp);
607 PMeshOnEntityImpl aMeshOnEntity = boost::get<1>(aFindTimeStamp);
608 PMeshOnEntityImpl aVTKMeshOnEntity = boost::get<2>(aFindTimeStamp);
609 PValForTimeImpl aValForTime = boost::get<4>(aFindTimeStamp);
610 PFieldImpl aField = boost::get<3>(aFindTimeStamp);
612 VISU::TVTKAttribyteFilter& anAttribyteFilter = aValForTime->myAttribyteFilter;
613 VISU::TVTKSource& aSource = aValForTime->myStorage;
614 TOutput* anOutput = NULL;
617 if(aSource.GetPointer())
618 return aSource.GetPointer();
619 else if(anAttribyteFilter.GetPointer())
620 return anAttribyteFilter->GetUnstructuredGridOutput();
622 LoadFieldOnMesh(aMesh,aMeshOnEntity,aField,aValForTime);
624 VISU::TVTKExtractFilter& anExtractFilter = aField->myExtractFilter;
625 if(anExtractFilter.GetPointer() == NULL){
626 anExtractFilter = VISU_ExtractUnstructuredGrid::New();
627 anExtractFilter->Delete();
628 //anExtractFilter->DebugOn();
630 LoadMeshOnEntity(aVTKMeshOnEntity);
631 }catch(std::exception& exc){
632 aVTKMeshOnEntity = aMeshOnEntity;
633 MSG(MYDEBUG,"Follow exception was occured :\n"<<exc.what());
635 aVTKMeshOnEntity = aMeshOnEntity;
636 MSG(MYDEBUG,"Unknown exception was occured!");
638 GetMeshOnEntity(aVTKMeshOnEntity->myMeshName,aVTKMeshOnEntity->myEntity);
640 anExtractFilter->SetInput(aVTKMeshOnEntity->myStorage.GetPointer());
641 ::InitProfile(anExtractFilter,aMeshOnEntity,aValForTime);
643 if(!anExtractFilter->IsRemoving()){
644 aSource = TOutput::New();
646 aSource->ShallowCopy(aVTKMeshOnEntity->myStorage.GetPointer());
647 ::GetTimeStamp(aSource,aField,aValForTime);
648 anOutput = aSource.GetPointer();
650 anAttribyteFilter = vtkFieldDataToAttributeDataFilter::New();
651 anAttribyteFilter->Delete();
652 //anAttribyteFilter->DebugOn();
654 VISU::TVTKMergetFilter& aMergeFilter = aValForTime->myMergeFilter;
655 aMergeFilter = vtkMergeDataObjectFilter::New();
656 aMergeFilter->Delete();
657 //aMergeFilter->DebugOn();
659 ::GetTimeStamp(anAttribyteFilter,aMergeFilter,anExtractFilter,
661 anOutput = anAttribyteFilter->GetUnstructuredGridOutput();
663 if(MYDEBUGWITHFILES){
664 string aMeshName = QString(theMeshName.c_str()).simplifyWhiteSpace().latin1();
665 string aFieldName = QString(theFieldName.c_str()).simplifyWhiteSpace().latin1();
666 string aPrefix = string("/users/")+getenv("USER")+"/"+getenv("USER")+"-";
667 string aFileName = aPrefix + aMeshName + dtos("-%d-",int(theEntity)) +
668 aFieldName + dtos("-%d",theStampsNum) + "-Conv.vtk";
669 VISU::WriteToFile(anOutput,aFileName);
672 GetTimeStampSize(theMeshName,theEntity,theFieldName,theStampsNum);
673 vtkUnstructuredGrid *aDataSet = anAttribyteFilter->GetUnstructuredGridOutput();
675 if(theEntity == VISU::NODE_ENTITY)
676 MSG(MYVTKDEBUG,"GetTimeStampOnMesh - GetData() = "<<float(aDataSet->GetPointData()->GetActualMemorySize()*1000));
678 MSG(MYVTKDEBUG,"GetMeshOnEntity - GetData() = "<<float(aDataSet->GetCellData()->GetActualMemorySize()*1000));
679 MSG(MYVTKDEBUG,"GetTimeStampOnMesh - GetActualMemorySize() = "<<float(aDataSet->GetActualMemorySize()*1000));
683 aSource = vtkSmartPointerBase();
684 anAttribyteFilter = vtkSmartPointerBase();
691 VISU_Convertor_impl::FindMesh(const string& theMeshName)
694 TMeshMap::iterator aMeshMapIter = myMeshMap.find(theMeshName);
695 if(aMeshMapIter == myMeshMap.end())
696 EXCEPTION(runtime_error,"FindMesh >> There is no mesh with the name - '"<<theMeshName<<"'!!!");
698 PMeshImpl aMesh = aMeshMapIter->second;
703 VISU_Convertor_impl::TFindMeshOnEntity
704 VISU_Convertor_impl::FindMeshOnEntity(const string& theMeshName,
705 const VISU::TEntity& theEntity,
706 const string& theFamilyName)
708 PMeshImpl aMesh = FindMesh(theMeshName);
709 VISU::TMeshOnEntityMap& aMeshOnEntityMap = aMesh->myMeshOnEntityMap;
710 VISU::TMeshOnEntityMap::const_iterator aMeshOnEntityMapIter = aMeshOnEntityMap.find(theEntity);
711 if(aMeshOnEntityMapIter == aMeshOnEntityMap.end())
712 EXCEPTION(runtime_error,"FindMeshOnEntity >> There is no mesh on the entity - "<<theEntity<<"!!!");
714 PMeshOnEntityImpl aMeshOnEntity = aMeshOnEntityMapIter->second;
716 return TFindMeshOnEntity(aMesh,
717 aMeshOnEntityMap[theEntity],
718 GetFamily(aMeshOnEntity,theFamilyName));
722 float VISU_Convertor_impl::GetSize() {
724 const VISU::TMeshMap& aMeshMap = GetMeshMap();
725 VISU::TMeshMap::const_iterator aMeshMapIter = aMeshMap.begin();
726 for(; aMeshMapIter != aMeshMap.end(); aMeshMapIter++){
727 const string& aMeshName = aMeshMapIter->first;
728 const VISU::PMesh aMesh = aMeshMapIter->second;
729 const VISU::TMeshOnEntityMap& aMeshOnEntityMap = aMesh->myMeshOnEntityMap;
730 VISU::TMeshOnEntityMap::const_iterator aMeshOnEntityMapIter;
732 aMeshOnEntityMapIter = aMeshOnEntityMap.begin();
733 for(; aMeshOnEntityMapIter != aMeshOnEntityMap.end(); aMeshOnEntityMapIter++){
734 const VISU::TEntity& anEntity = aMeshOnEntityMapIter->first;
735 const VISU::PMeshOnEntity aMeshOnEntity = aMeshOnEntityMapIter->second;
736 const VISU::TFieldMap& aFieldMap = aMeshOnEntity->myFieldMap;
737 VISU::TFieldMap::const_iterator aFieldMapIter = aFieldMap.begin();
738 for(; aFieldMapIter != aFieldMap.end(); aFieldMapIter++){
739 const string& aFieldName = aFieldMapIter->first;
740 const VISU::PField aField = aFieldMapIter->second;
741 const VISU::TValField& aValField = aField->myValField;
742 VISU::TValField::const_iterator aValFieldIter = aValField.begin();
743 for(; aValFieldIter != aValField.end(); aValFieldIter++){
744 int aTimeStamp = aValFieldIter->first;
745 aResult += GetTimeStampSize(aMeshName,anEntity,aFieldName,aTimeStamp);
749 const VISU::TGroupMap& aGroupMap = aMesh->myGroupMap;
750 VISU::TGroupMap::const_iterator aGroupMapIter = aGroupMap.begin();
751 for(; aGroupMapIter != aGroupMap.end(); aGroupMapIter++){
752 const string& aGroupName = aGroupMapIter->first;
753 aResult += GetMeshOnGroupSize(aMeshName,aGroupName);
756 const VISU::TFamilyMap& aFamilyMap = aMeshOnEntity->myFamilyMap;
757 VISU::TFamilyMap::const_iterator aFamilyMapIter = aFamilyMap.begin();
758 for(; aFamilyMapIter != aFamilyMap.end(); aFamilyMapIter++){
759 const string& aFamilyName = aFamilyMapIter->first;
760 aResult += GetMeshOnEntitySize(aMeshName,anEntity,aFamilyName);
762 //Import mesh on entity
763 aResult += GetMeshOnEntitySize(aMeshName,anEntity);
766 MSG(MYDEBUG,"GetSize - aResult = "<<float(aResult));
771 float VISU_Convertor_impl::GetMeshOnEntitySize(const std::string& theMeshName,
772 const VISU::TEntity& theEntity,
773 const std::string& theFamilyName)
775 TFindMeshOnEntity aFindMeshOnEntity =
776 FindMeshOnEntity(theMeshName,theEntity,theFamilyName);
777 PMeshImpl aMesh = boost::get<0>(aFindMeshOnEntity);
778 PFamilyImpl aFamily = boost::get<2>(aFindMeshOnEntity);
779 PMeshOnEntityImpl aMeshOnEntity = boost::get<1>(aFindMeshOnEntity);
781 vtkIdType aPointsSize = 3*aMesh->myNbPoints*sizeof(VISU::TCoord);
782 vtkIdType aNbCells, aCellsSize;
785 aNbCells = aMeshOnEntity->myNbCells;
786 aCellsSize = aMeshOnEntity->myCellsSize;
788 aNbCells = aFamily->myNbCells;
789 aCellsSize = aFamily->myCellsSize;
792 vtkIdType aConnectivitySize = aCellsSize*sizeof(vtkIdType);
793 vtkIdType aTypesSize = aNbCells*sizeof(char);
794 vtkIdType aLocationsSize = aNbCells*sizeof(int);
795 float aNbCellsPerPoint = aCellsSize / aNbCells - 1;
796 vtkIdType aLinksSize = aMesh->myNbPoints *
797 (vtkIdType(sizeof(vtkIdType)*aNbCellsPerPoint) + sizeof(vtkCellLinks::Link));
799 vtkIdType aResult = aPointsSize + aConnectivitySize + aTypesSize + aLocationsSize + aLinksSize;
801 MSG(MYVTKDEBUG,"GetMeshOnEntitySize - aPointsSize = "<<float(aPointsSize));
802 MSG(MYVTKDEBUG,"GetMeshOnEntitySize - aConnectivitySize = "<<float(aConnectivitySize));
803 MSG(MYVTKDEBUG,"GetMeshOnEntitySize - aTypesSize = "<<float(aTypesSize));
804 MSG(MYVTKDEBUG,"GetMeshOnEntitySize - aLocationsSize = "<<float(aLocationsSize));
805 MSG(MYVTKDEBUG,"GetMeshOnEntitySize - aLinksSize = "<<float(aLinksSize));
807 MSG(MYDEBUG,"GetMeshOnEntitySize - aResult = "<<float(aResult)<<"; theMeshName = '"<<theMeshName<<
808 "'; theEntity = "<<theEntity<<"; theFamilyName = '"<<theFamilyName<<"'");
810 aResult = vtkIdType(aResult*ERR_SIZE_CALC);
815 VISU_Convertor_impl::TFindMeshOnGroup
816 VISU_Convertor_impl::FindMeshOnGroup(const std::string& theMeshName,
817 const std::string& theGroupName)
819 PMeshImpl aMesh = FindMesh(theMeshName);
820 VISU::TGroupMap& aGroupMap = aMesh->myGroupMap;
821 VISU::TGroupMap::iterator aGroupMapIter = aGroupMap.find(theGroupName);
822 if(aGroupMapIter == aGroupMap.end())
823 EXCEPTION(runtime_error,"FindMesh >> There is no the group in the mesh!!! - '"<<theGroupName<<"'");
825 VISU::PGroupImpl aGroup = aGroupMapIter->second;
826 return TFindMeshOnGroup(aMesh,aGroup);
830 float VISU_Convertor_impl::GetMeshOnGroupSize(const std::string& theMeshName,
831 const std::string& theGroupName)
833 TFindMeshOnGroup aFindMeshOnGroup = FindMeshOnGroup(theMeshName,theGroupName);
834 PMeshImpl aMesh = boost::get<0>(aFindMeshOnGroup);
835 PGroupImpl aGroup = boost::get<1>(aFindMeshOnGroup);
837 vtkIdType aPointsSize = 3*aMesh->myNbPoints*sizeof(VISU::TCoord);
838 vtkIdType aNbCells = aGroup->myNbCells, aCellsSize = aGroup->myCellsSize;
839 vtkIdType aConnectivityAndTypesSize = aCellsSize*sizeof(vtkIdType);
840 vtkIdType aLocationsSize = aNbCells*sizeof(int);
841 float aNbCellsPerPoint = aCellsSize / aNbCells - 1;
842 vtkIdType aLinksSize = aMesh->myNbPoints *
843 (vtkIdType(sizeof(vtkIdType)*aNbCellsPerPoint) + sizeof(short));
845 vtkIdType aResult = aPointsSize + aConnectivityAndTypesSize + aLocationsSize + aLinksSize;
847 MSG(MYVTKDEBUG,"GetMeshOnGroupSize - aPointsSize = "<<float(aPointsSize));
848 MSG(MYVTKDEBUG,"GetMeshOnGroupSize - aConnectivityAndTypesSize = "<<float(aConnectivityAndTypesSize));
849 MSG(MYVTKDEBUG,"GetMeshOnGroupSize - aLocationsSize = "<<float(aLocationsSize));
850 MSG(MYVTKDEBUG,"GetMeshOnGroupSize - aLinksSize = "<<float(aLinksSize));
852 MSG(MYDEBUG,"GetMeshOnGroupSize - aResult = "<<float(aResult)<<"; theMeshName = '"
853 <<theMeshName<<"'; theGroupName = '"<<theGroupName<<"'");
855 aResult = vtkIdType(aResult*ERR_SIZE_CALC);
859 VISU_Convertor_impl::TFindField
860 VISU_Convertor_impl::FindField(const string& theMeshName,
861 const VISU::TEntity& theEntity,
862 const string& theFieldName)
864 TFindMeshOnEntity aFindMeshOnEntity = FindMeshOnEntity(theMeshName,theEntity,"");
865 PMeshImpl aMesh = boost::get<0>(aFindMeshOnEntity);;
866 PFamilyImpl aFamily = boost::get<2>(aFindMeshOnEntity);
867 PMeshOnEntityImpl aMeshOnEntity = boost::get<1>(aFindMeshOnEntity);
869 VISU::TMeshOnEntityMap& aMeshOnEntityMap = aMesh->myMeshOnEntityMap;
870 PMeshOnEntityImpl aVTKMeshOnEntity;
871 if(theEntity == VISU::NODE_ENTITY){
872 if(aMeshOnEntityMap.find(VISU::CELL_ENTITY) != aMeshOnEntityMap.end())
873 aVTKMeshOnEntity = aMeshOnEntityMap[VISU::CELL_ENTITY];
874 else if(aMeshOnEntityMap.find(VISU::FACE_ENTITY) != aMeshOnEntityMap.end())
875 aVTKMeshOnEntity = aMeshOnEntityMap[VISU::FACE_ENTITY];
876 else if(aMeshOnEntityMap.find(VISU::NODE_ENTITY) != aMeshOnEntityMap.end())
877 aVTKMeshOnEntity = aMeshOnEntityMap[VISU::EDGE_ENTITY];
879 aVTKMeshOnEntity = aMeshOnEntity;
881 VISU::TFieldMap& aFieldMap = aMeshOnEntity->myFieldMap;
882 VISU::TFieldMap::const_iterator aFieldIter= aFieldMap.find(theFieldName);
883 if(aFieldIter == aFieldMap.end())
884 EXCEPTION(runtime_error,"FindField >> There is no field on the mesh!!!");
886 PFieldImpl aField = aFieldIter->second;
888 return TFindField(aMesh,
895 float VISU_Convertor_impl::GetFieldOnMeshSize(const std::string& theMeshName,
896 const VISU::TEntity& theEntity,
897 const std::string& theFieldName)
899 TFindField aFindField = FindField(theMeshName,theEntity,theFieldName);
900 PMeshOnEntityImpl aVTKMeshOnEntity = boost::get<2>(aFindField);
901 PField aField = boost::get<3>(aFindField);
903 float aMeshSize = GetMeshOnEntitySize(theMeshName,aVTKMeshOnEntity->myEntity);
904 float aFieldOnMeshSize = float(aField->myDataSize*sizeof(float)*aField->myValField.size()*ERR_SIZE_CALC);
905 float aResult = aMeshSize + aFieldOnMeshSize;
907 MSG(MYVTKDEBUG,"GetFieldOnMeshSize - aFieldOnMeshSize = "<<float(aFieldOnMeshSize));
908 MSG(MYDEBUG,"GetFieldOnMeshSize - aResult = "<<float(aResult)<<"; theMeshName = '"<<theMeshName<<
909 "'; theEntity = "<<theEntity<<"; theFieldName = '"<<theFieldName<<"'");
915 VISU_Convertor_impl::TFindTimeStamp
916 VISU_Convertor_impl::FindTimeStamp(const std::string& theMeshName,
917 const VISU::TEntity& theEntity,
918 const std::string& theFieldName,
921 TFindField aFindField = FindField(theMeshName,theEntity,theFieldName);
922 PField aField = boost::get<3>(aFindField);
924 VISU::TValField& aValField = aField->myValField;
925 VISU::TValField::const_iterator aValFieldIter= aValField.find(theStampsNum);
926 if(aValFieldIter == aValField.end())
927 EXCEPTION(runtime_error,"FindTimeStamp >> There is no field with the timestamp!!!");
929 PMeshImpl aMesh = boost::get<0>(aFindField);
930 PMeshOnEntityImpl aMeshOnEntity = boost::get<1>(aFindField);
931 PMeshOnEntityImpl aVTKMeshOnEntity = boost::get<2>(aFindField);
932 PValForTimeImpl aValForTime = aValFieldIter->second;
934 return TFindTimeStamp(aMesh,
942 float VISU_Convertor_impl::GetTimeStampSize(const std::string& theMeshName,
943 const VISU::TEntity& theEntity,
944 const std::string& theFieldName,
947 TFindTimeStamp aFindTimeStamp =
948 FindTimeStamp(theMeshName,theEntity,theFieldName,theStampsNum);
949 PMeshOnEntityImpl aVTKMeshOnEntity = boost::get<2>(aFindTimeStamp);
950 PField aField = boost::get<3>(aFindTimeStamp);
952 float aMeshSize = GetMeshOnEntitySize(theMeshName,aVTKMeshOnEntity->myEntity);
953 float aTimeStampSize = float(aField->myDataSize*sizeof(float) * ERR_SIZE_CALC);
954 float aResult = aMeshSize + aTimeStampSize;
956 MSG(MYDEBUG && MYVTKDEBUG,"GetTimeStampSize - aTimeStampSize = "<<float(aTimeStampSize));
957 MSG(MYDEBUG,"GetTimeStampSize - aResult = "<<float(aResult)<<
958 "; theMeshName = '"<<theMeshName<<"'; theEntity = "<<theEntity<<
959 "; theFieldName = '"<<theFieldName<<"'; theStampsNum = "<<theStampsNum);
966 VISU_Convertor_impl::GetField(const string& theMeshName,
967 VISU::TEntity theEntity,
968 const string& theFieldName)
970 TFindField aFindField = FindField(theMeshName,theEntity,theFieldName);
971 PField aField = boost::get<3>(aFindField);
976 const VISU::PValForTime
977 VISU_Convertor_impl::GetTimeStamp(const std::string& theMeshName,
978 const VISU::TEntity& theEntity,
979 const std::string& theFieldName,
982 TFindTimeStamp aFindTimeStamp =
983 FindTimeStamp(theMeshName,theEntity,theFieldName,theStampsNum);
984 PValForTime aValForTime = boost::get<4>(aFindTimeStamp);