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;
211 theStorage->SetPoints(aPoints);
215 inline void PrintCells(int& theStartId,
216 vtkCellArray* theConnectivity,
217 const VISU::TMeshOnEntityImpl::TConnect& theVector)
219 vtkIdList *anIdList = vtkIdList::New();
220 int kEnd = theVector.size();
221 anIdList->SetNumberOfIds(kEnd);
222 for(int k = 0; k < kEnd; k++)
223 anIdList->SetId(k,theVector[k]);
224 theConnectivity->InsertNextCell(anIdList);
228 void GetCellsOnEntity(VISU::TVTKSource& theStorage,
229 const VISU::PMeshOnEntityImpl theMeshOnEntity,
230 const string& theFamilyName)
232 //Check on existing family
233 PFamilyImpl aFamily = GetFamily(theMeshOnEntity,theFamilyName);
235 pair<int,int> aCellsDim = theMeshOnEntity->GetCellsDims(theFamilyName);
236 int aNbCells = aCellsDim.first, aCellsSize = aCellsDim.second;
237 vtkCellArray* aConnectivity = vtkCellArray::New();
238 aConnectivity->Allocate(aCellsSize,0);
239 vtkUnsignedCharArray* aCellTypesArray = vtkUnsignedCharArray::New();
240 aCellTypesArray->SetNumberOfComponents(1);
241 aCellTypesArray->SetNumberOfTuples(aNbCells);
242 MSG(MYDEBUG,"GetCellsOnEntity - isFamilyPresent = "<<bool(aFamily));
243 const VISU::TMeshOnEntityImpl::TCellsConn &aCellsConn = theMeshOnEntity->myCellsConn;
244 VISU::TMeshOnEntityImpl::TCellsConn::const_iterator aCellsConnIter = aCellsConn.begin();
245 for(int i = 0, j = 0; aCellsConnIter != aCellsConn.end(); aCellsConnIter++){
246 const VISU::TMeshOnEntityImpl::TConnForCellType& anArray = aCellsConnIter->second;
247 int aVtkType = aCellsConnIter->first;
248 MSG(MYDEBUG,"GetCellsOnEntity - aVtkType = "<<aVtkType<<"; anArray.size() = "<<anArray.size());
250 for(int k = 0, kEnd = anArray.size(); k < kEnd; k++, i++){
251 PrintCells(i,aConnectivity,anArray[k]);
252 aCellTypesArray->SetValue(j++,(unsigned char)aVtkType);
255 const VISU::TFamilyImpl::TSubMesh& aSubMesh = aFamily->mySubMesh;
257 EXCEPTION(runtime_error,"GetCells >> There is no elements on the family !!!");
258 VISU::TFamilyImpl::TSubMesh::const_iterator aSubMeshIter = aSubMesh.find(aVtkType);
259 if(aSubMeshIter == aSubMesh.end()) continue;
260 const VISU::TFamilyImpl::TSubMeshOnCellType& aSubMeshOnCellType = aSubMeshIter->second;
261 MSG(MYDEBUG,"GetCellsOnEntity - aSubMeshOnCellType.size() = "<<aSubMeshOnCellType.size());
262 VISU::TFamilyImpl::TSubMeshOnCellType::const_iterator aSubMeshOnCellTypeIter = aSubMeshOnCellType.begin();
263 for(; aSubMeshOnCellTypeIter != aSubMeshOnCellType.end(); aSubMeshOnCellTypeIter++, i++){
264 PrintCells(i,aConnectivity,anArray[*aSubMeshOnCellTypeIter]);
265 aCellTypesArray->SetValue(j++,(unsigned char)aVtkType);
269 vtkIdType *pts = 0, npts = 0;
270 vtkIntArray* aCellLocationsArray = vtkIntArray::New();
271 aCellLocationsArray->SetNumberOfComponents(1);
272 aCellLocationsArray->SetNumberOfTuples(aNbCells);
273 aConnectivity->InitTraversal();
274 for(int i=0; aConnectivity->GetNextCell(npts,pts); i++)
275 aCellLocationsArray->SetValue(i,aConnectivity->GetTraversalLocation(npts));
276 theStorage->SetCells(aCellTypesArray,aCellLocationsArray,aConnectivity);
277 if(MYVTKDEBUG) aConnectivity->DebugOn();
278 aCellLocationsArray->Delete();
279 aCellTypesArray->Delete();
280 aConnectivity->Delete();
284 void GetCellsOnGroup(VISU::TVTKSource& theStorage,
285 VISU::PMeshImpl theMesh,
286 const VISU::TFamilyAndEntitySet& theFamilyAndEntitySet)
288 //Calculate dimentions of the group
289 int aNbCells = 0, aCellsSize = 0;
290 VISU::TFamilyAndEntitySet::const_iterator aFamilyAndEntitySetIter = theFamilyAndEntitySet.begin();
291 for(; aFamilyAndEntitySetIter != theFamilyAndEntitySet.end(); aFamilyAndEntitySetIter++){
292 const VISU::TFamilyAndEntity& aFamilyAndEntity = *aFamilyAndEntitySetIter;
293 const string& aFamilyName = aFamilyAndEntity.first;
294 const VISU::TEntity& anEntity = aFamilyAndEntity.second;
295 const VISU::PMeshOnEntity aMeshOnEntity = theMesh->myMeshOnEntityMap[anEntity];
296 pair<int,int> aCellsDim = aMeshOnEntity->GetCellsDims(aFamilyName);
297 aNbCells += aCellsDim.first;
298 aCellsSize += aCellsDim.second;
300 vtkCellArray* aConnectivity = vtkCellArray::New();
301 aConnectivity->Allocate(aCellsSize,0);
302 vtkUnsignedCharArray* aCellTypesArray = vtkUnsignedCharArray::New();
303 aCellTypesArray->SetNumberOfComponents(1);
304 aCellTypesArray->SetNumberOfTuples(aNbCells);
305 aFamilyAndEntitySetIter = theFamilyAndEntitySet.begin();
306 for(int i = 0, j = 0; aFamilyAndEntitySetIter != theFamilyAndEntitySet.end(); aFamilyAndEntitySetIter++){
307 const VISU::TFamilyAndEntity& aFamilyAndEntity = *aFamilyAndEntitySetIter;
308 const string& aFamilyName = aFamilyAndEntity.first;
309 const VISU::TEntity& anEntity = aFamilyAndEntity.second;
310 PMeshOnEntityImpl aMeshOnEntity = theMesh->myMeshOnEntityMap[anEntity];
311 PFamilyImpl aFamily = GetFamily(aMeshOnEntity,aFamilyName);
312 const VISU::TMeshOnEntityImpl::TCellsConn &aCellsConn = aMeshOnEntity->myCellsConn;
313 VISU::TMeshOnEntityImpl::TCellsConn::const_iterator aCellsConnIter = aCellsConn.begin();
314 for(; aCellsConnIter != aCellsConn.end(); aCellsConnIter++){
315 const VISU::TMeshOnEntityImpl::TConnForCellType& anArray = aCellsConnIter->second;
316 int aVtkType = aCellsConnIter->first;
317 MSG(MYDEBUG,"GetCellsOnGroup - aVtkType = "<<aVtkType<<"; anArray.size() = "<<anArray.size());
318 const VISU::TFamilyImpl::TSubMesh& aSubMesh = aFamily->mySubMesh;
320 EXCEPTION(runtime_error,"GetCells >> There is no elements on the family !!!");
321 VISU::TFamilyImpl::TSubMesh::const_iterator aSubMeshIter = aSubMesh.find(aVtkType);
322 if(aSubMeshIter == aSubMesh.end()) continue;
323 const VISU::TFamilyImpl::TSubMeshOnCellType& aSubMeshOnCellType = aSubMeshIter->second;
324 MSG(MYDEBUG,"GetCellsOnGroup - aSubMeshOnCellType.size() = "<<aSubMeshOnCellType.size());
325 VISU::TFamilyImpl::TSubMeshOnCellType::const_iterator aSubMeshOnCellTypeIter = aSubMeshOnCellType.begin();
326 for(; aSubMeshOnCellTypeIter != aSubMeshOnCellType.end(); aSubMeshOnCellTypeIter++, i++){
327 PrintCells(i,aConnectivity,anArray[*aSubMeshOnCellTypeIter]);
328 aCellTypesArray->SetValue(j++,(unsigned char)aVtkType);
332 vtkIdType *pts = 0, npts = 0;
333 vtkIntArray* aCellLocationsArray = vtkIntArray::New();
334 aCellLocationsArray->SetNumberOfComponents(1);
335 aCellLocationsArray->SetNumberOfTuples(aNbCells);
336 aConnectivity->InitTraversal();
337 for(int i = 0; aConnectivity->GetNextCell(npts,pts); i++)
338 aCellLocationsArray->SetValue(i,aConnectivity->GetTraversalLocation(npts));
339 theStorage->SetCells(aCellTypesArray,aCellLocationsArray,aConnectivity);
340 aCellLocationsArray->Delete();
341 aCellTypesArray->Delete();
342 aConnectivity->Delete();
346 void InitProfile(VISU::TVTKExtractFilter& theFilter,
347 PMeshOnEntityImpl theMeshOnEntity,
348 PValForTimeImpl theValForTime)
350 const VISU::TValForTimeImpl::TValForCells& aValForCells = theValForTime->myValForCells;
351 const VISU::TMeshOnEntityImpl::TCellsConn &aCellsConn = theMeshOnEntity->myCellsConn;
352 VISU::TMeshOnEntityImpl::TCellsConn::const_iterator aCellsConnIter = aCellsConn.begin();
353 for(; aCellsConnIter != aCellsConn.end(); aCellsConnIter++){
354 const vtkIdType& aCellType = aCellsConnIter->first;
355 if(aValForCells.find(aCellType) == aValForCells.end())
356 theFilter->RemoveCellsWithType(aCellType);
361 void GetValsOnTimeStamp(vtkFloatArray *theFloatArray,
362 const vtkIdType& theNumberOfTuples,
363 const std::string& theFieldName,
364 VISU::PFieldImpl theField,
365 VISU::PValForTimeImpl theValForTime)
367 //theFloatArray->DebugOn();
368 theFloatArray->SetNumberOfTuples(theNumberOfTuples);
369 theFloatArray->SetName(theFieldName.c_str());
370 MSG(MYDEBUG,"GetValsOnTimeStamp - theNumberOfTuples = "<<theNumberOfTuples);
371 const VISU::TValForTimeImpl::TValForCells& aValForCells = theValForTime->myValForCells;
372 VISU::TValForTimeImpl::TValForCells::const_iterator aValForCellsIter = aValForCells.begin();
373 for(int k = 0; aValForCellsIter != aValForCells.end(); aValForCellsIter++) {
374 const VISU::TValForTimeImpl::TValForCellsWithType& anArray = aValForCellsIter->second;
375 int iEnd = anArray.size()/theField->myNbComp;
376 int aVtkType = aValForCellsIter->first;
377 MSG(MYDEBUG,"GetValsOnTimeStamp - iEnd = "<<iEnd<<"; aVtkType = "<<aVtkType);
378 switch(theField->myNbComp) {
380 for (int i = 0; i < iEnd; i++)
381 theFloatArray->SetTuple1(k++,anArray[i]);
384 for (int i = 0, ji = 0; i < iEnd; ++i, ji = i*2)
385 theFloatArray->SetTuple3(k++,anArray[ji],anArray[ji+1],0.0);
388 for (int i = 0, ji = 0; i < iEnd; ++i, ji = i*3)
389 theFloatArray->SetTuple3(k++,anArray[ji],anArray[ji+1],anArray[ji+2]);
392 for (int i = 0, ji = 0; i < iEnd; ++i, ji = i*4)
393 theFloatArray->SetTuple3(k++,anArray[ji],anArray[ji+2],0.0);
396 for (int i = 0, ji = 0; i < iEnd; ++i, ji = i*4)
397 theFloatArray->SetTuple3(k++,anArray[ji],anArray[ji+2],anArray[ji+5]);
400 EXCEPTION(runtime_error,"GetValsOnTimeStamp - There is no an algorithm for representation of the field !!!");
405 string GenerateFieldName(const VISU::PFieldImpl theField,
406 const VISU::PValForTimeImpl theValForTime)
408 const VISU::TTime& aTime = theValForTime->myTime;
409 string aFieldName = theField->myMeshName + ", " + theField->myName + ": " +
410 VISU_Convertor::GenerateName(aTime);
414 void GetTimeStamp(VISU::TVTKSource& theStorage,
415 const VISU::PFieldImpl theField,
416 const VISU::PValForTimeImpl theValForTime)
418 int aNumberOfTuples = theField->myDataSize/theField->myNbComp;
419 string aFieldName = GenerateFieldName(theField,theValForTime);
420 MSG(MYDEBUG,"GetTimeStamp(TVTKSource) - aFieldName = "<<aFieldName<<
421 "; aNumberOfTuples = "<<aNumberOfTuples);
423 vtkDataSetAttributes* aDataSetAttributes;
424 switch(theField->myEntity){
425 case VISU::NODE_ENTITY :
426 aDataSetAttributes = theStorage->GetPointData();
429 aDataSetAttributes = theStorage->GetCellData();
432 vtkFloatArray *aFloatArray = vtkFloatArray::New();
433 switch(theField->myNbComp) {
435 aFloatArray->SetNumberOfComponents(1);
436 aDataSetAttributes->SetScalars(aFloatArray);
439 aFloatArray->SetNumberOfComponents(3);
440 aDataSetAttributes->SetVectors(aFloatArray);
443 GetValsOnTimeStamp(aFloatArray,aNumberOfTuples,aFieldName,theField,theValForTime);
446 void GetTimeStamp(VISU::TVTKAttribyteFilter& theAttribyteFilter,
447 VISU::TVTKMergetFilter& theMergeFilter,
448 VISU::TVTKExtractFilter& theExtractFilter,
449 const VISU::PFieldImpl theField,
450 const VISU::PValForTimeImpl theValForTime)
452 int aNumberOfTuples = theField->myDataSize/theField->myNbComp;
453 string aFieldName = GenerateFieldName(theField,theValForTime);
454 MSG(MYDEBUG,"GetTimeStamp(TVTKAttribyteFilter) - aFieldName = "<<aFieldName<<
455 "; aNumberOfTuples = "<<aNumberOfTuples);
457 vtkDataObject* aDataObject = vtkDataObject::New();
458 theMergeFilter->SetDataObject(aDataObject);
459 aDataObject->Delete();
461 theMergeFilter->SetInput(theExtractFilter->GetOutput());
462 theAttribyteFilter->SetInput(theMergeFilter->GetOutput());
464 switch(theField->myEntity){
465 case VISU::NODE_ENTITY :
466 theMergeFilter->SetOutputFieldToPointDataField();
467 theAttribyteFilter->SetInputFieldToPointDataField();
468 theAttribyteFilter->SetOutputAttributeDataToPointData();
471 theMergeFilter->SetOutputFieldToCellDataField();
472 theAttribyteFilter->SetInputFieldToCellDataField();
473 theAttribyteFilter->SetOutputAttributeDataToCellData();
476 vtkFloatArray *aFloatArray = vtkFloatArray::New();
477 switch(theField->myNbComp) {
479 aFloatArray->SetNumberOfComponents(1);
480 theAttribyteFilter->SetScalarComponent(0,aFieldName.c_str(),0);
483 aFloatArray->SetNumberOfComponents(3);
484 theAttribyteFilter->SetVectorComponent(0,aFieldName.c_str(),0);
485 theAttribyteFilter->SetVectorComponent(1,aFieldName.c_str(),1);
486 theAttribyteFilter->SetVectorComponent(2,aFieldName.c_str(),2);
489 vtkFieldData* aFieldData = aDataObject->GetFieldData();
490 aFieldData->AddArray(aFloatArray);
491 aFloatArray->Delete();
493 GetValsOnTimeStamp(aFloatArray,aNumberOfTuples,aFieldName,theField,theValForTime);
497 VISU_Convertor_impl::VISU_Convertor_impl() {
501 VISU_Convertor_impl::~VISU_Convertor_impl() {}
503 VISU_Convertor::TOutput*
504 VISU_Convertor_impl::GetMeshOnEntity(const string& theMeshName,
505 const VISU::TEntity& theEntity,
506 const string& theFamilyName)
508 MSG(MYDEBUG,"GetMeshOnEntity - theMeshName = '"<<theMeshName<<
509 "'; theEntity = "<<theEntity<<"; theFamilyName = '"<<theFamilyName<<"'");
510 //Cheching possibility do the query
511 TFindMeshOnEntity aFindMeshOnEntity =
512 FindMeshOnEntity(theMeshName,theEntity,theFamilyName);
513 VISU::TVTKSource* pSource;
514 PMeshImpl aMesh = boost::get<0>(aFindMeshOnEntity);;
515 PFamilyImpl aFamily = boost::get<2>(aFindMeshOnEntity);
516 PMeshOnEntityImpl aMeshOnEntity = boost::get<1>(aFindMeshOnEntity);
518 pSource = &(aFamily->myStorage);
520 pSource = &(aMeshOnEntity->myStorage);
521 VISU::TVTKSource& aSource = *pSource;
524 if(aSource.GetPointer() == NULL){
525 aSource = TOutput::New();
527 if(MYVTKDEBUG) aSource->DebugOn();
528 LoadMeshOnEntity(aMeshOnEntity,theFamilyName);
529 GetPoints(aSource,aMesh);
530 GetCellsOnEntity(aSource,aMeshOnEntity,theFamilyName);
531 if(MYDEBUGWITHFILES){
532 string aMeshName = QString(theMeshName.c_str()).simplifyWhiteSpace().latin1();
533 string aFamilyName = QString(theFamilyName.c_str()).simplifyWhiteSpace().latin1();
534 string aFileName = string("/users/")+getenv("USER")+"/"+getenv("USER")+"-";
535 aFileName += aMeshName + dtos("-%d-",int(theEntity)) + aFamilyName + "-Conv.vtk";
536 VISU::WriteToFile(aSource.GetPointer(),aFileName);
540 GetMeshOnEntitySize(theMeshName,theEntity,theFamilyName);
541 vtkUnstructuredGrid* aDataSet = aSource.GetPointer();
543 MSG(MYVTKDEBUG,"GetMeshOnEntity - GetPoints() = "<<float(aDataSet->GetPoints()->GetActualMemorySize()*1000));
544 MSG(MYVTKDEBUG,"GetMeshOnEntity - GetCells() = "<<float(aDataSet->GetCells()->GetActualMemorySize()*1000));
545 MSG(MYVTKDEBUG,"GetMeshOnEntity - GetCellTypesArray() = "<<float(aDataSet->GetCellTypesArray()->GetActualMemorySize()*1000));
546 MSG(MYVTKDEBUG,"GetMeshOnEntity - GetCellLocationsArray() = "<<float(aDataSet->GetCellLocationsArray()->GetActualMemorySize()*1000));
547 aDataSet->BuildLinks();
548 MSG(MYVTKDEBUG,"GetMeshOnEntity - GetCellLinks() = "<<float(aDataSet->GetCellLinks()->GetActualMemorySize()*1000));
549 MSG(MYVTKDEBUG,"GetMeshOnEntity - GetActualMemorySize() = "<<float(aDataSet->GetActualMemorySize()*1000));
552 aSource = vtkSmartPointerBase();
555 return aSource.GetPointer();
558 VISU_Convertor::TOutput*
559 VISU_Convertor_impl::GetMeshOnGroup(const string& theMeshName,
560 const string& theGroupName)
562 MSG(MYDEBUG,"GetMeshOnGroup - theMeshName = '"<<theMeshName<<
563 "'; theGroupName = '"<<theGroupName<<"'");
564 //Cheching possibility do the query
565 TFindMeshOnGroup aFindMeshOnGroup = FindMeshOnGroup(theMeshName,theGroupName);
566 PMeshImpl aMesh = boost::get<0>(aFindMeshOnGroup);
567 PGroupImpl aGroup = boost::get<1>(aFindMeshOnGroup);
568 const VISU::TFamilyAndEntitySet& aFamilyAndEntitySet = aGroup->myFamilyAndEntitySet;
569 VISU::TVTKSource& aSource = aGroup->myStorage;
572 if(aSource.GetPointer() == NULL){
573 aSource = TOutput::New();
575 LoadMeshOnGroup(aMesh,aFamilyAndEntitySet);
576 GetPoints(aSource,aMesh);
577 GetCellsOnGroup(aSource,aMesh,aFamilyAndEntitySet);
578 if(MYDEBUGWITHFILES){
579 string aMeshName = QString(theMeshName.c_str()).simplifyWhiteSpace().latin1();
580 string aGroupName = QString(theGroupName.c_str()).simplifyWhiteSpace().latin1();
581 string aFileName = string("/users/")+getenv("USER")+"/"+getenv("USER")+"-";
582 aFileName += aMeshName + "-" + aGroupName + "-Conv.vtk";
583 VISU::WriteToFile(aSource.GetPointer(),aFileName);
587 aSource = vtkSmartPointerBase();
590 return aSource.GetPointer();
593 VISU_Convertor::TOutput*
594 VISU_Convertor_impl::GetTimeStampOnMesh(const string& theMeshName,
595 const VISU::TEntity& theEntity,
596 const string& theFieldName,
599 MSG(MYDEBUG,"GetTimeStampOnMesh - theMeshName = '"<<theMeshName<<"; theEntity = "<<theEntity);
600 MSG(MYDEBUG,"GetTimeStampOnMesh - theFieldName = '"<<theFieldName<<"'; theStampsNum = "<<theStampsNum);
602 //Cheching possibility do the query
603 TFindTimeStamp aFindTimeStamp =
604 FindTimeStamp(theMeshName,theEntity,theFieldName,theStampsNum);
605 PMeshImpl aMesh = boost::get<0>(aFindTimeStamp);
606 PMeshOnEntityImpl aMeshOnEntity = boost::get<1>(aFindTimeStamp);
607 PMeshOnEntityImpl aVTKMeshOnEntity = boost::get<2>(aFindTimeStamp);
608 PValForTimeImpl aValForTime = boost::get<4>(aFindTimeStamp);
609 PFieldImpl aField = boost::get<3>(aFindTimeStamp);
611 VISU::TVTKAttribyteFilter& anAttribyteFilter = aValForTime->myAttribyteFilter;
612 VISU::TVTKSource& aSource = aValForTime->myStorage;
613 TOutput* anOutput = NULL;
616 if(aSource.GetPointer())
617 return aSource.GetPointer();
618 else if(anAttribyteFilter.GetPointer())
619 return anAttribyteFilter->GetUnstructuredGridOutput();
621 LoadFieldOnMesh(aMesh,aMeshOnEntity,aField,aValForTime);
623 VISU::TVTKExtractFilter& anExtractFilter = aField->myExtractFilter;
624 if(anExtractFilter.GetPointer() == NULL){
625 anExtractFilter = VISU_ExtractUnstructuredGrid::New();
626 anExtractFilter->Delete();
627 //anExtractFilter->DebugOn();
629 LoadMeshOnEntity(aVTKMeshOnEntity);
630 }catch(std::exception& exc){
631 aVTKMeshOnEntity = aMeshOnEntity;
632 MSG(MYDEBUG,"Follow exception was occured :\n"<<exc.what());
634 aVTKMeshOnEntity = aMeshOnEntity;
635 MSG(MYDEBUG,"Unknown exception was occured!");
637 GetMeshOnEntity(aVTKMeshOnEntity->myMeshName,aVTKMeshOnEntity->myEntity);
639 anExtractFilter->SetInput(aVTKMeshOnEntity->myStorage.GetPointer());
640 ::InitProfile(anExtractFilter,aMeshOnEntity,aValForTime);
642 if(!anExtractFilter->IsRemoving()){
643 aSource = TOutput::New();
645 aSource->ShallowCopy(aVTKMeshOnEntity->myStorage.GetPointer());
646 ::GetTimeStamp(aSource,aField,aValForTime);
647 anOutput = aSource.GetPointer();
649 anAttribyteFilter = vtkFieldDataToAttributeDataFilter::New();
650 anAttribyteFilter->Delete();
651 //anAttribyteFilter->DebugOn();
653 VISU::TVTKMergetFilter& aMergeFilter = aValForTime->myMergeFilter;
654 aMergeFilter = vtkMergeDataObjectFilter::New();
655 aMergeFilter->Delete();
656 //aMergeFilter->DebugOn();
658 ::GetTimeStamp(anAttribyteFilter,aMergeFilter,anExtractFilter,
660 anOutput = anAttribyteFilter->GetUnstructuredGridOutput();
662 if(MYDEBUGWITHFILES){
663 string aMeshName = QString(theMeshName.c_str()).simplifyWhiteSpace().latin1();
664 string aFieldName = QString(theFieldName.c_str()).simplifyWhiteSpace().latin1();
665 string aPrefix = string("/users/")+getenv("USER")+"/"+getenv("USER")+"-";
666 string aFileName = aPrefix + aMeshName + dtos("-%d-",int(theEntity)) +
667 aFieldName + dtos("-%d",theStampsNum) + "-Conv.vtk";
668 VISU::WriteToFile(anOutput,aFileName);
671 GetTimeStampSize(theMeshName,theEntity,theFieldName,theStampsNum);
672 vtkUnstructuredGrid *aDataSet = anAttribyteFilter->GetUnstructuredGridOutput();
674 if(theEntity == VISU::NODE_ENTITY)
675 MSG(MYVTKDEBUG,"GetTimeStampOnMesh - GetData() = "<<float(aDataSet->GetPointData()->GetActualMemorySize()*1000));
677 MSG(MYVTKDEBUG,"GetMeshOnEntity - GetData() = "<<float(aDataSet->GetCellData()->GetActualMemorySize()*1000));
678 MSG(MYVTKDEBUG,"GetTimeStampOnMesh - GetActualMemorySize() = "<<float(aDataSet->GetActualMemorySize()*1000));
682 aSource = vtkSmartPointerBase();
683 anAttribyteFilter = vtkSmartPointerBase();
690 VISU_Convertor_impl::FindMesh(const string& theMeshName)
693 TMeshMap::iterator aMeshMapIter = myMeshMap.find(theMeshName);
694 if(aMeshMapIter == myMeshMap.end())
695 EXCEPTION(runtime_error,"FindMesh >> There is no mesh with the name - '"<<theMeshName<<"'!!!");
697 PMeshImpl aMesh = aMeshMapIter->second;
702 VISU_Convertor_impl::TFindMeshOnEntity
703 VISU_Convertor_impl::FindMeshOnEntity(const string& theMeshName,
704 const VISU::TEntity& theEntity,
705 const string& theFamilyName)
707 PMeshImpl aMesh = FindMesh(theMeshName);
708 VISU::TMeshOnEntityMap& aMeshOnEntityMap = aMesh->myMeshOnEntityMap;
709 VISU::TMeshOnEntityMap::const_iterator aMeshOnEntityMapIter = aMeshOnEntityMap.find(theEntity);
710 if(aMeshOnEntityMapIter == aMeshOnEntityMap.end())
711 EXCEPTION(runtime_error,"FindMeshOnEntity >> There is no mesh on the entity - "<<theEntity<<"!!!");
713 PMeshOnEntityImpl aMeshOnEntity = aMeshOnEntityMapIter->second;
715 return TFindMeshOnEntity(aMesh,
716 aMeshOnEntityMap[theEntity],
717 GetFamily(aMeshOnEntity,theFamilyName));
721 float VISU_Convertor_impl::GetSize() {
723 const VISU::TMeshMap& aMeshMap = GetMeshMap();
724 VISU::TMeshMap::const_iterator aMeshMapIter = aMeshMap.begin();
725 for(; aMeshMapIter != aMeshMap.end(); aMeshMapIter++){
726 const string& aMeshName = aMeshMapIter->first;
727 const VISU::PMesh aMesh = aMeshMapIter->second;
728 const VISU::TMeshOnEntityMap& aMeshOnEntityMap = aMesh->myMeshOnEntityMap;
729 VISU::TMeshOnEntityMap::const_iterator aMeshOnEntityMapIter;
731 aMeshOnEntityMapIter = aMeshOnEntityMap.begin();
732 for(; aMeshOnEntityMapIter != aMeshOnEntityMap.end(); aMeshOnEntityMapIter++){
733 const VISU::TEntity& anEntity = aMeshOnEntityMapIter->first;
734 const VISU::PMeshOnEntity aMeshOnEntity = aMeshOnEntityMapIter->second;
735 const VISU::TFieldMap& aFieldMap = aMeshOnEntity->myFieldMap;
736 VISU::TFieldMap::const_iterator aFieldMapIter = aFieldMap.begin();
737 for(; aFieldMapIter != aFieldMap.end(); aFieldMapIter++){
738 const string& aFieldName = aFieldMapIter->first;
739 const VISU::PField aField = aFieldMapIter->second;
740 const VISU::TValField& aValField = aField->myValField;
741 VISU::TValField::const_iterator aValFieldIter = aValField.begin();
742 for(; aValFieldIter != aValField.end(); aValFieldIter++){
743 int aTimeStamp = aValFieldIter->first;
744 aResult += GetTimeStampSize(aMeshName,anEntity,aFieldName,aTimeStamp);
748 const VISU::TGroupMap& aGroupMap = aMesh->myGroupMap;
749 VISU::TGroupMap::const_iterator aGroupMapIter = aGroupMap.begin();
750 for(; aGroupMapIter != aGroupMap.end(); aGroupMapIter++){
751 const string& aGroupName = aGroupMapIter->first;
752 aResult += GetMeshOnGroupSize(aMeshName,aGroupName);
755 const VISU::TFamilyMap& aFamilyMap = aMeshOnEntity->myFamilyMap;
756 VISU::TFamilyMap::const_iterator aFamilyMapIter = aFamilyMap.begin();
757 for(; aFamilyMapIter != aFamilyMap.end(); aFamilyMapIter++){
758 const string& aFamilyName = aFamilyMapIter->first;
759 aResult += GetMeshOnEntitySize(aMeshName,anEntity,aFamilyName);
761 //Import mesh on entity
762 aResult += GetMeshOnEntitySize(aMeshName,anEntity);
765 MSG(MYDEBUG,"GetSize - aResult = "<<float(aResult));
770 float VISU_Convertor_impl::GetMeshOnEntitySize(const std::string& theMeshName,
771 const VISU::TEntity& theEntity,
772 const std::string& theFamilyName)
774 TFindMeshOnEntity aFindMeshOnEntity =
775 FindMeshOnEntity(theMeshName,theEntity,theFamilyName);
776 PMeshImpl aMesh = boost::get<0>(aFindMeshOnEntity);
777 PFamilyImpl aFamily = boost::get<2>(aFindMeshOnEntity);
778 PMeshOnEntityImpl aMeshOnEntity = boost::get<1>(aFindMeshOnEntity);
780 vtkIdType aPointsSize = 3*aMesh->myNbPoints*sizeof(VISU::TCoord);
781 vtkIdType aNbCells, aCellsSize;
784 aNbCells = aMeshOnEntity->myNbCells;
785 aCellsSize = aMeshOnEntity->myCellsSize;
787 aNbCells = aFamily->myNbCells;
788 aCellsSize = aFamily->myCellsSize;
791 vtkIdType aConnectivitySize = aCellsSize*sizeof(vtkIdType);
792 vtkIdType aTypesSize = aNbCells*sizeof(char);
793 vtkIdType aLocationsSize = aNbCells*sizeof(int);
794 float aNbCellsPerPoint = aCellsSize / aNbCells - 1;
795 vtkIdType aLinksSize = aMesh->myNbPoints *
796 (vtkIdType(sizeof(vtkIdType)*aNbCellsPerPoint) + sizeof(vtkCellLinks::Link));
798 vtkIdType aResult = aPointsSize + aConnectivitySize + aTypesSize + aLocationsSize + aLinksSize;
800 MSG(MYVTKDEBUG,"GetMeshOnEntitySize - aPointsSize = "<<float(aPointsSize));
801 MSG(MYVTKDEBUG,"GetMeshOnEntitySize - aConnectivitySize = "<<float(aConnectivitySize));
802 MSG(MYVTKDEBUG,"GetMeshOnEntitySize - aTypesSize = "<<float(aTypesSize));
803 MSG(MYVTKDEBUG,"GetMeshOnEntitySize - aLocationsSize = "<<float(aLocationsSize));
804 MSG(MYVTKDEBUG,"GetMeshOnEntitySize - aLinksSize = "<<float(aLinksSize));
806 MSG(MYDEBUG,"GetMeshOnEntitySize - aResult = "<<float(aResult)<<"; theMeshName = '"<<theMeshName<<
807 "'; theEntity = "<<theEntity<<"; theFamilyName = '"<<theFamilyName<<"'");
809 aResult = vtkIdType(aResult*ERR_SIZE_CALC);
814 VISU_Convertor_impl::TFindMeshOnGroup
815 VISU_Convertor_impl::FindMeshOnGroup(const std::string& theMeshName,
816 const std::string& theGroupName)
818 PMeshImpl aMesh = FindMesh(theMeshName);
819 VISU::TGroupMap& aGroupMap = aMesh->myGroupMap;
820 VISU::TGroupMap::iterator aGroupMapIter = aGroupMap.find(theGroupName);
821 if(aGroupMapIter == aGroupMap.end())
822 EXCEPTION(runtime_error,"FindMesh >> There is no the group in the mesh!!! - '"<<theGroupName<<"'");
824 VISU::PGroupImpl aGroup = aGroupMapIter->second;
825 return TFindMeshOnGroup(aMesh,aGroup);
829 float VISU_Convertor_impl::GetMeshOnGroupSize(const std::string& theMeshName,
830 const std::string& theGroupName)
832 TFindMeshOnGroup aFindMeshOnGroup = FindMeshOnGroup(theMeshName,theGroupName);
833 PMeshImpl aMesh = boost::get<0>(aFindMeshOnGroup);
834 PGroupImpl aGroup = boost::get<1>(aFindMeshOnGroup);
836 vtkIdType aPointsSize = 3*aMesh->myNbPoints*sizeof(VISU::TCoord);
837 vtkIdType aNbCells = aGroup->myNbCells, aCellsSize = aGroup->myCellsSize;
838 vtkIdType aConnectivityAndTypesSize = aCellsSize*sizeof(vtkIdType);
839 vtkIdType aLocationsSize = aNbCells*sizeof(int);
840 float aNbCellsPerPoint = aCellsSize / aNbCells - 1;
841 vtkIdType aLinksSize = aMesh->myNbPoints *
842 (vtkIdType(sizeof(vtkIdType)*aNbCellsPerPoint) + sizeof(short));
844 vtkIdType aResult = aPointsSize + aConnectivityAndTypesSize + aLocationsSize + aLinksSize;
846 MSG(MYVTKDEBUG,"GetMeshOnGroupSize - aPointsSize = "<<float(aPointsSize));
847 MSG(MYVTKDEBUG,"GetMeshOnGroupSize - aConnectivityAndTypesSize = "<<float(aConnectivityAndTypesSize));
848 MSG(MYVTKDEBUG,"GetMeshOnGroupSize - aLocationsSize = "<<float(aLocationsSize));
849 MSG(MYVTKDEBUG,"GetMeshOnGroupSize - aLinksSize = "<<float(aLinksSize));
851 MSG(MYDEBUG,"GetMeshOnGroupSize - aResult = "<<float(aResult)<<"; theMeshName = '"
852 <<theMeshName<<"'; theGroupName = '"<<theGroupName<<"'");
854 aResult = vtkIdType(aResult*ERR_SIZE_CALC);
858 VISU_Convertor_impl::TFindField
859 VISU_Convertor_impl::FindField(const string& theMeshName,
860 const VISU::TEntity& theEntity,
861 const string& theFieldName)
863 TFindMeshOnEntity aFindMeshOnEntity = FindMeshOnEntity(theMeshName,theEntity,"");
864 PMeshImpl aMesh = boost::get<0>(aFindMeshOnEntity);;
865 PFamilyImpl aFamily = boost::get<2>(aFindMeshOnEntity);
866 PMeshOnEntityImpl aMeshOnEntity = boost::get<1>(aFindMeshOnEntity);
868 VISU::TMeshOnEntityMap& aMeshOnEntityMap = aMesh->myMeshOnEntityMap;
869 PMeshOnEntityImpl aVTKMeshOnEntity;
870 if(theEntity == VISU::NODE_ENTITY){
871 if(aMeshOnEntityMap.find(VISU::CELL_ENTITY) != aMeshOnEntityMap.end())
872 aVTKMeshOnEntity = aMeshOnEntityMap[VISU::CELL_ENTITY];
873 else if(aMeshOnEntityMap.find(VISU::FACE_ENTITY) != aMeshOnEntityMap.end())
874 aVTKMeshOnEntity = aMeshOnEntityMap[VISU::FACE_ENTITY];
875 else if(aMeshOnEntityMap.find(VISU::NODE_ENTITY) != aMeshOnEntityMap.end())
876 aVTKMeshOnEntity = aMeshOnEntityMap[VISU::EDGE_ENTITY];
878 aVTKMeshOnEntity = aMeshOnEntity;
880 VISU::TFieldMap& aFieldMap = aMeshOnEntity->myFieldMap;
881 VISU::TFieldMap::const_iterator aFieldIter= aFieldMap.find(theFieldName);
882 if(aFieldIter == aFieldMap.end())
883 EXCEPTION(runtime_error,"FindField >> There is no field on the mesh!!!");
885 PFieldImpl aField = aFieldIter->second;
887 return TFindField(aMesh,
894 float VISU_Convertor_impl::GetFieldOnMeshSize(const std::string& theMeshName,
895 const VISU::TEntity& theEntity,
896 const std::string& theFieldName)
898 TFindField aFindField = FindField(theMeshName,theEntity,theFieldName);
899 PMeshOnEntityImpl aVTKMeshOnEntity = boost::get<2>(aFindField);
900 PField aField = boost::get<3>(aFindField);
902 float aMeshSize = GetMeshOnEntitySize(theMeshName,aVTKMeshOnEntity->myEntity);
903 float aFieldOnMeshSize = float(aField->myDataSize*sizeof(float)*aField->myValField.size()*ERR_SIZE_CALC);
904 float aResult = aMeshSize + aFieldOnMeshSize;
906 MSG(MYVTKDEBUG,"GetFieldOnMeshSize - aFieldOnMeshSize = "<<float(aFieldOnMeshSize));
907 MSG(MYDEBUG,"GetFieldOnMeshSize - aResult = "<<float(aResult)<<"; theMeshName = '"<<theMeshName<<
908 "'; theEntity = "<<theEntity<<"; theFieldName = '"<<theFieldName<<"'");
914 VISU_Convertor_impl::TFindTimeStamp
915 VISU_Convertor_impl::FindTimeStamp(const std::string& theMeshName,
916 const VISU::TEntity& theEntity,
917 const std::string& theFieldName,
920 TFindField aFindField = FindField(theMeshName,theEntity,theFieldName);
921 PField aField = boost::get<3>(aFindField);
923 VISU::TValField& aValField = aField->myValField;
924 VISU::TValField::const_iterator aValFieldIter= aValField.find(theStampsNum);
925 if(aValFieldIter == aValField.end())
926 EXCEPTION(runtime_error,"FindTimeStamp >> There is no field with the timestamp!!!");
928 PMeshImpl aMesh = boost::get<0>(aFindField);
929 PMeshOnEntityImpl aMeshOnEntity = boost::get<1>(aFindField);
930 PMeshOnEntityImpl aVTKMeshOnEntity = boost::get<2>(aFindField);
931 PValForTimeImpl aValForTime = aValFieldIter->second;
933 return TFindTimeStamp(aMesh,
941 float VISU_Convertor_impl::GetTimeStampSize(const std::string& theMeshName,
942 const VISU::TEntity& theEntity,
943 const std::string& theFieldName,
946 TFindTimeStamp aFindTimeStamp =
947 FindTimeStamp(theMeshName,theEntity,theFieldName,theStampsNum);
948 PMeshOnEntityImpl aVTKMeshOnEntity = boost::get<2>(aFindTimeStamp);
949 PField aField = boost::get<3>(aFindTimeStamp);
951 float aMeshSize = GetMeshOnEntitySize(theMeshName,aVTKMeshOnEntity->myEntity);
952 float aTimeStampSize = float(aField->myDataSize*sizeof(float) * ERR_SIZE_CALC);
953 float aResult = aMeshSize + aTimeStampSize;
955 MSG(MYDEBUG && MYVTKDEBUG,"GetTimeStampSize - aTimeStampSize = "<<float(aTimeStampSize));
956 MSG(MYDEBUG,"GetTimeStampSize - aResult = "<<float(aResult)<<
957 "; theMeshName = '"<<theMeshName<<"'; theEntity = "<<theEntity<<
958 "; theFieldName = '"<<theFieldName<<"'; theStampsNum = "<<theStampsNum);
965 VISU_Convertor_impl::GetField(const string& theMeshName,
966 VISU::TEntity theEntity,
967 const string& theFieldName)
969 TFindField aFindField = FindField(theMeshName,theEntity,theFieldName);
970 PField aField = boost::get<3>(aFindField);
975 const VISU::PValForTime
976 VISU_Convertor_impl::GetTimeStamp(const std::string& theMeshName,
977 const VISU::TEntity& theEntity,
978 const std::string& theFieldName,
981 TFindTimeStamp aFindTimeStamp =
982 FindTimeStamp(theMeshName,theEntity,theFieldName,theStampsNum);
983 PValForTime aValForTime = boost::get<4>(aFindTimeStamp);