1 // Copyright (C) 2007-2008 CEA/DEN, EDF R&D, OPEN CASCADE
3 // Copyright (C) 2003-2007 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
22 // VISU OBJECT : interactive object for VISU entities implementation
27 #include "VISU_MeshValue.hxx"
28 #include "VISU_ElnoMeshValue.hxx"
29 #include "VISU_Structures_impl.hxx"
30 #include "VISU_ConvertorUtils.hxx"
32 #include "VISU_PointCoords.hxx"
33 #include "VISU_VTKTypeList.hxx"
35 #include <vtkUnstructuredGrid.h>
36 #include <vtkPolyData.h>
38 #include <vtkPointData.h>
39 #include <vtkCellData.h>
41 #include <vtkCharArray.h>
42 #include <vtkUnsignedCharArray.h>
43 #include <vtkShortArray.h>
44 #include <vtkUnsignedShortArray.h>
45 #include <vtkIntArray.h>
46 #include <vtkUnsignedIntArray.h>
47 #include <vtkLongArray.h>
48 #include <vtkUnsignedLongArray.h>
49 #include <vtkFloatArray.h>
50 #include <vtkDoubleArray.h>
56 static int MYDEBUG = 0;
58 static int MYDEBUG = 0;
63 //---------------------------------------------------------------
65 GenerateFieldName(const PFieldImpl& theField,
66 const PValForTimeImpl& theValForTime)
68 const VISU::TTime& aTime = theValForTime->myTime;
69 std::string aFieldName = theField->myMeshName + ", " + theField->myName + ": " +
70 VISU_Convertor::GenerateName(aTime);
75 //---------------------------------------------------------------
78 ::Init(vtkIdType theNbElem,
83 myNbGauss = theNbGauss;
85 myStep = theNbComp*theNbGauss;
113 return myNbElem * myStep;
117 //----------------------------------------------------------------------------
118 template<int EDataType>
120 InitTimeStampOnProfile(const PUnstructuredGrid& theSource,
121 const PFieldImpl& theField,
122 const PValForTimeImpl& theValForTime,
123 const VISU::TEntity& theEntity);
126 //----------------------------------------------------------------------------
128 GetTimeStampOnProfile(const PUnstructuredGrid& theSource,
129 const PFieldImpl& theField,
130 const PValForTimeImpl& theValForTime,
131 const VISU::TEntity& theEntity)
133 vtkIdType aDataType = theField->GetDataType();
136 InitTimeStampOnProfile<VTK_DOUBLE>(theSource, theField, theValForTime, theEntity);
139 InitTimeStampOnProfile<VTK_FLOAT>(theSource, theField, theValForTime, theEntity);
142 InitTimeStampOnProfile<VTK_INT>(theSource, theField, theValForTime, theEntity);
145 InitTimeStampOnProfile<VTK_LONG>(theSource, theField, theValForTime, theEntity);
148 EXCEPTION(std::runtime_error,
149 "GetTimeStampOnProfile - handling unsupported data type - "<<aDataType);
154 //----------------------------------------------------------------------------
155 template<int EDataType>
156 struct TDataArrayHolder
158 typedef typename TL::TEnum2VTKArrayType<EDataType>::TResult TVTKDataArray;
159 typedef typename TL::TEnum2VTKBasicType<EDataType>::TResult TVTKBasicType;
160 TVTKDataArray* myDataArray;
162 TDataArrayHolder(TVTKDataArray* theDataArray):
163 myDataArray(theDataArray)
167 WritePointer(TVTKDataArray* theDataArray,
168 vtkIdType theTupleId,
169 TVTKBasicType* thePointer)
171 vtkIdType aNumberOfComponents = theDataArray->GetNumberOfComponents();
172 vtkIdType aPosition = theTupleId * aNumberOfComponents;
173 TVTKBasicType *aPtr = theDataArray->WritePointer(aPosition, aNumberOfComponents);
174 for(vtkIdType anId = 0; anId < aNumberOfComponents; anId++)
175 *aPtr++ = *thePointer++;
180 SetTuple(vtkIdType theTupleId,
181 TVTKBasicType* thePointer)
183 this->WritePointer(myDataArray, theTupleId, thePointer);
188 //----------------------------------------------------------------------------
189 template<int EDataType>
190 struct TDataArrayHolder2: TDataArrayHolder<EDataType>
192 typedef TDataArrayHolder<EDataType> TSuperClass;
193 typedef typename TSuperClass::TVTKDataArray TVTKDataArray;
194 typedef typename TSuperClass::TVTKBasicType TVTKBasicType;
195 TVTKDataArray* myDataArray2;
197 TDataArrayHolder2(TVTKDataArray* theDataArray,
198 TVTKDataArray* theDataArray2):
199 TSuperClass(theDataArray),
200 myDataArray2(theDataArray2)
205 SetTuple(vtkIdType theTupleId,
206 TVTKBasicType* thePointer)
208 this->WritePointer(this->myDataArray, theTupleId, thePointer);
209 this->WritePointer(this->myDataArray2, theTupleId, thePointer);
214 //----------------------------------------------------------------------------
215 template<int EDataType>
216 struct TTimeStampOnProfileInitArray
218 typedef typename TL::TEnum2VTKArrayType<EDataType>::TResult TVTKDataArray;
219 typedef typename TL::TEnum2VTKBasicType<EDataType>::TResult TVTKBasicType;
220 typedef TTMeshValue<TVTKBasicType> TMeshValue;
221 typedef MED::SharedPtr<TMeshValue> TMeshValuePtr;
223 typedef TDataArrayHolder<EDataType> TTDataArrayHolder;
224 typedef MED::SharedPtr<TTDataArrayHolder> PDataArrayHolder;
225 PDataArrayHolder myDataArrayHolder;
227 TTimeStampOnProfileInitArray(const PDataArrayHolder& theDataArrayHolder):
228 myDataArrayHolder(theDataArrayHolder)
232 Execute(const PFieldImpl& theField,
233 const PValForTimeImpl& theValForTime,
234 const TGaussMetric theGaussMetric = VISU::AVERAGE_METRIC)
236 vtkIdType aNbComp = theField->myNbComp;
237 vtkIdType aSize = std::max(vtkIdType(3), aNbComp);
238 TVector<TVTKBasicType> aDataValues (aSize);
240 const TGeom2MeshValue& aGeom2MeshValue = theValForTime->GetGeom2MeshValue();
241 TGeom2MeshValue::const_iterator anIter = aGeom2MeshValue.begin();
242 for (int aTupleId = 0; anIter != aGeom2MeshValue.end(); anIter++) {
243 EGeometry aEGeom = anIter->first;
244 const TMeshValuePtr aMeshValue = anIter->second;
246 vtkIdType aNbElem = aMeshValue->GetNbElem();
247 vtkIdType aNbGauss = aMeshValue->GetNbGauss();
250 "- aEGeom = "<<aEGeom<<
251 "; aNbElem = "<<aNbElem<<
252 "; aNbGauss = "<<aNbGauss<<
255 for (vtkIdType iElem = 0; iElem < aNbElem; iElem++, aTupleId++) {
256 typename TMeshValue::TCValueSliceArr aValueSliceArr = aMeshValue->GetCompValueSliceArr(iElem);
258 for (vtkIdType iComp = 0; iComp < aNbComp; iComp++) {
259 const typename TMeshValue::TCValueSlice& aValueSlice = aValueSliceArr[iComp];
260 //jfa//aDataValues[iComp] = TVTKBasicType();
261 for (vtkIdType iGauss = 0; iGauss < aNbGauss; iGauss++) {
262 TVTKBasicType aValue = aValueSlice[iGauss];
263 //jfa//if (iGauss == 0 && theGaussMetric != VISU::AVERAGE_METRIC) {
265 aDataValues[iComp] = aValue; // init
268 switch (theGaussMetric) {
269 case VISU::AVERAGE_METRIC:
270 aDataValues[iComp] += aValue; break;
271 case VISU::MINIMUM_METRIC:
272 aDataValues[iComp] = std::min( aValue, aDataValues[iComp] ); break;
273 case VISU::MAXIMUM_METRIC:
274 aDataValues[iComp] = std::max( aValue, aDataValues[iComp] ); break;
278 if (theGaussMetric == VISU::AVERAGE_METRIC)
279 aDataValues[iComp] /= aNbGauss;
282 this->myDataArrayHolder->SetTuple(aTupleId, &aDataValues[0]);
288 template<int EDataType>
289 struct TTimeStampOnProfileInitModulus
291 typedef typename TL::TEnum2VTKArrayType<EDataType>::TResult TVTKDataArray;
292 typedef typename TL::TEnum2VTKBasicType<EDataType>::TResult TVTKBasicType;
293 typedef TTMeshValue<TVTKBasicType> TMeshValue;
294 typedef MED::SharedPtr<TMeshValue> TMeshValuePtr;
296 typedef TDataArrayHolder<EDataType> TTDataArrayHolder;
297 typedef MED::SharedPtr<TTDataArrayHolder> PDataArrayHolder;
298 PDataArrayHolder myDataArrayHolder;
300 TTimeStampOnProfileInitModulus(const PDataArrayHolder& theDataArrayHolder):
301 myDataArrayHolder(theDataArrayHolder)
305 Execute(const PFieldImpl& theField, const PValForTimeImpl& theValForTime)
307 vtkIdType aNbComp = theField->myNbComp;
308 vtkIdType aSize = vtkIdType(3); // Minimum, Maximum and Average modulus
309 TVector<TVTKBasicType> aDataValues (aSize);
311 const TGeom2MeshValue& aGeom2MeshValue = theValForTime->GetGeom2MeshValue();
312 TGeom2MeshValue::const_iterator anIter = aGeom2MeshValue.begin();
313 for (int aTupleId = 0; anIter != aGeom2MeshValue.end(); anIter++) {
314 EGeometry aEGeom = anIter->first;
315 const TMeshValuePtr aMeshValue = anIter->second;
317 vtkIdType aNbElem = aMeshValue->GetNbElem();
318 vtkIdType aNbGauss = aMeshValue->GetNbGauss();
321 "- aEGeom = "<<aEGeom<<
322 "; aNbElem = "<<aNbElem<<
323 "; aNbGauss = "<<aNbGauss<<
326 for (vtkIdType iElem = 0; iElem < aNbElem; iElem++, aTupleId++) {
327 typename TMeshValue::TCValueSliceArr aValueSliceArr = aMeshValue->GetCompValueSliceArr(iElem);
329 // modules of all gauss points
330 TVector<TVTKBasicType> aModules (aNbGauss);
332 for (vtkIdType iComp = 0; iComp < aNbComp; iComp++) {
333 const typename TMeshValue::TCValueSlice& aValueSlice = aValueSliceArr[iComp];
334 for (vtkIdType iGauss = 0; iGauss < aNbGauss; iGauss++) {
335 TVTKBasicType aValue = aValueSlice[iGauss];
337 // modules of all gauss points
339 aModules[iGauss] = aValue * aValue; // init
341 aModules[iGauss] += aValue * aValue;
345 TVTKBasicType aModule = (TVTKBasicType)sqrt(aModules[0]);
346 aDataValues[0] = aModule; // init Min
347 aDataValues[1] = aModule; // init Max
348 aDataValues[2] = aModule; // init Average
350 for (vtkIdType ig = 0; ig < aNbGauss; ig++) {
351 aModule = (TVTKBasicType)sqrt(aModules[ig]);
353 aDataValues[0] = std::min(TVTKBasicType(aModule), aDataValues[0]); // Min
354 aDataValues[1] = std::max(TVTKBasicType(aModule), aDataValues[1]); // Max
355 aDataValues[2] += aModule; // Average
357 aDataValues[2] /= aNbGauss; // Average
359 this->myDataArrayHolder->SetTuple(aTupleId, &aDataValues[0]);
366 //----------------------------------------------------------------------------
367 template<int EDataType>
369 InitTimeStampOnProfile(const PUnstructuredGrid& theSource,
370 const PFieldImpl& theField,
371 const PValForTimeImpl& theValForTime,
372 const VISU::TEntity& theEntity)
374 vtkIdType aNbTuples = theField->myDataSize / theField->myNbComp;
375 std::string aFieldName = VISU::GenerateFieldName(theField, theValForTime);
377 vtkDataSetAttributes* aDataSetAttributes;
378 switch ( theEntity ) {
379 case VISU::NODE_ENTITY :
380 aDataSetAttributes = theSource->GetPointData();
383 aDataSetAttributes = theSource->GetCellData();
386 typedef typename TL::TEnum2VTKArrayType<EDataType>::TResult TVTKDataArray;
387 TVTKDataArray *aSelectedDataArray = TVTKDataArray::New();
388 vtkIdType aNbComp = theField->myNbComp;
392 aSelectedDataArray->SetNumberOfComponents( 1 );
393 aDataSetAttributes->SetScalars( aSelectedDataArray );
396 aSelectedDataArray->SetNumberOfComponents( 3 );
397 aDataSetAttributes->SetVectors( aSelectedDataArray );
399 aSelectedDataArray->SetNumberOfTuples( aNbTuples );
400 aSelectedDataArray->SetName( aFieldName.c_str() );
402 TVTKDataArray *aFullDataArray = TVTKDataArray::New();
403 aFullDataArray->SetNumberOfComponents( aNbComp );
404 aFullDataArray->SetNumberOfTuples( aNbTuples );
405 aFullDataArray->SetName( "VISU_FIELD" );
406 aDataSetAttributes->AddArray( aFullDataArray );
408 INITMSG(MYDEBUG,"InitTimeStampOnProfile "<<
409 "- theEntity = "<<theEntity<<
410 "; aNbTuples = "<<aNbTuples<<
411 "; aNbComp = "<<aNbComp<<
414 TTimerLog aTimerLog(MYDEBUG,"InitTimeStampOnProfile");
416 const TGeom2MeshValue& aGeom2MeshValue = theValForTime->GetGeom2MeshValue();
417 typedef typename TL::TEnum2VTKBasicType< EDataType >::TResult TVTKBasicType;
418 typedef TTMeshValue< TVTKBasicType > TMeshValue;
419 typedef MED::SharedPtr< TMeshValue > TMeshValuePtr;
421 typedef TDataArrayHolder< EDataType > TTDataArrayHolder;
422 typedef MED::SharedPtr< TTDataArrayHolder > PDataArrayHolder;
424 TMeshValuePtr aMeshValue = theValForTime->GetFirstMeshValue();
425 if ( aGeom2MeshValue.size() == 1 && aMeshValue->GetNbGauss() == 1 ) {
426 aFullDataArray->SetVoidArray(aMeshValue->GetPointer(),
429 INITMSG(MYDEBUG,"InitTimeStampOnProfile - aFullDataArray->SetVoidArray()"<<std::endl);
430 if ( aNbComp == 1 ) {
431 aSelectedDataArray->SetVoidArray( aMeshValue->GetPointer(),
434 INITMSG(MYDEBUG,"InitTimeStampOnProfile - aSelectedDataArray->SetVoidArray()"<<std::endl);
436 PDataArrayHolder aDataArrayHolder(new TTDataArrayHolder(aSelectedDataArray));
437 TTimeStampOnProfileInitArray<EDataType>(aDataArrayHolder).Execute(theField, theValForTime);
441 typedef TDataArrayHolder2<EDataType> TTDataArrayHolder2;
442 PDataArrayHolder aDataArrayHolder(new TTDataArrayHolder2(aSelectedDataArray, aFullDataArray));
443 TTimeStampOnProfileInitArray<EDataType>(aDataArrayHolder).Execute(theField, theValForTime);
445 if ( theValForTime->GetMaxNbGauss() > 1 ) { // at least one of geometry contains multiple gauss points
446 TVTKDataArray *aGaussMinDataArray = TVTKDataArray::New();
447 aGaussMinDataArray->SetNumberOfComponents( aNbComp );
448 aGaussMinDataArray->SetNumberOfTuples( aNbTuples );
449 aGaussMinDataArray->SetName( "VISU_FIELD_GAUSS_MIN" );
450 aDataSetAttributes->AddArray( aGaussMinDataArray );
452 PDataArrayHolder aGaussMinDataArrayHolder(new TTDataArrayHolder(aGaussMinDataArray));
453 TTimeStampOnProfileInitArray<EDataType>(aGaussMinDataArrayHolder).Execute
454 (theField, theValForTime, VISU::MINIMUM_METRIC);
455 aGaussMinDataArray->Delete();
457 TVTKDataArray *aGaussMaxDataArray = TVTKDataArray::New();
458 aGaussMaxDataArray->SetNumberOfComponents( aNbComp );
459 aGaussMaxDataArray->SetNumberOfTuples( aNbTuples );
460 aGaussMaxDataArray->SetName( "VISU_FIELD_GAUSS_MAX" );
461 aDataSetAttributes->AddArray( aGaussMaxDataArray );
463 PDataArrayHolder aGaussMaxDataArrayHolder(new TTDataArrayHolder(aGaussMaxDataArray));
464 TTimeStampOnProfileInitArray<EDataType>(aGaussMaxDataArrayHolder).Execute
465 (theField, theValForTime, VISU::MAXIMUM_METRIC);
466 aGaussMaxDataArray->Delete();
468 TVTKDataArray *aGaussModulusDataArray = TVTKDataArray::New();
469 aGaussModulusDataArray->SetNumberOfComponents( 3 ); // Min, Max and Average
470 aGaussModulusDataArray->SetNumberOfTuples( aNbTuples );
471 aGaussModulusDataArray->SetName( "VISU_FIELD_GAUSS_MOD" );
472 aDataSetAttributes->AddArray( aGaussModulusDataArray );
474 PDataArrayHolder aGaussModulusDataArrayHolder(new TTDataArrayHolder(aGaussModulusDataArray));
475 TTimeStampOnProfileInitModulus<EDataType>(aGaussModulusDataArrayHolder).Execute
476 (theField, theValForTime);
477 aGaussModulusDataArray->Delete();
481 aSelectedDataArray->Delete();
482 aFullDataArray->Delete();
484 // Process the case for ELNO data
485 //-------------------------------
486 if ( theField->myIsELNO ) {
487 // To calculate effective number of components for the VTK compatibel ELNO data representation
488 vtkIdType aEffectNbTuples = 0;
489 TGeom2MeshValue::const_iterator anIter = aGeom2MeshValue.begin();
490 for ( ; anIter != aGeom2MeshValue.end(); anIter++ ) {
491 const PMeshValue& aMeshValue = anIter->second;
492 aEffectNbTuples += aMeshValue->GetNbElem() * aMeshValue->GetNbGauss();
495 vtkIdType anEffectNbComp = ( aEffectNbTuples * aNbComp ) / aNbTuples + 1;
497 // To create corresponding VTK representation for the ELNO data
498 TSetElnoNodeData< EDataType > aSetElnoNodeData( anEffectNbComp,
502 "ELNO_COMPONENT_MAPPER" );
504 std::vector< TVTKBasicType > aDataValues( aNbComp ); // To reserve a temproary value holder
506 // To initilize these VTK representation for the ELNO data from the MED
507 anIter = aGeom2MeshValue.begin();
508 for ( ; anIter != aGeom2MeshValue.end(); anIter++ ) {
509 EGeometry aEGeom = anIter->first;
510 const TMeshValuePtr aMeshValue = anIter->second;
512 vtkIdType aNbElem = aMeshValue->GetNbElem();
513 vtkIdType aNbGauss = aMeshValue->GetNbGauss();
516 "- aEGeom = "<<aEGeom<<
517 "; aNbElem = "<<aNbElem<<
518 "; aNbGauss = "<<aNbGauss<<
520 std::vector<int> med2visu(aNbGauss);
521 InitMed2VisuArray(med2visu,aEGeom);
522 for ( vtkIdType iElem = 0; iElem < aNbElem; iElem++ ) {
523 const typename TMeshValue::TValueSliceArr& aValueSliceArr =
524 aMeshValue->GetGaussValueSliceArr( iElem );
526 for( vtkIdType iGauss = 0; iGauss < aNbGauss; iGauss++ ) {
527 const typename TMeshValue::TCValueSlice& aValueSlice = aValueSliceArr[ med2visu[iGauss] ];
529 for( vtkIdType iComp = 0; iComp < aNbComp; iComp++ ) {
530 aDataValues[ iComp ] = aValueSlice[ iComp ];
533 aSetElnoNodeData.AddNextPointData( &aDataValues[ 0 ] );
536 aSetElnoNodeData.InsertNextCellData();
540 // Assign the ELNO data on the corresponding VTK data set attribute
541 aSetElnoNodeData.AddData( aDataSetAttributes );
543 //-------------------------------
547 //----------------------------------------------------------------------------
548 template<int EDataType>
550 InitTimeStampOnGaussMesh(const PPolyData& theSource,
551 const PFieldImpl& theField,
552 const PValForTimeImpl& theValForTime);
555 GetTimeStampOnGaussMesh(const PPolyData& theSource,
556 const PFieldImpl& theField,
557 const PValForTimeImpl& theValForTime)
559 vtkIdType aDataType = theField->GetDataType();
562 InitTimeStampOnGaussMesh<VTK_DOUBLE>(theSource, theField, theValForTime);
565 InitTimeStampOnGaussMesh<VTK_FLOAT>(theSource, theField, theValForTime);
568 InitTimeStampOnGaussMesh<VTK_INT>(theSource, theField, theValForTime);
571 InitTimeStampOnGaussMesh<VTK_LONG>(theSource, theField, theValForTime);
574 EXCEPTION(std::runtime_error,
575 "GetTimeStampOnGaussMesh - handling unsupported data type - "<<aDataType);
579 //----------------------------------------------------------------------------
580 template<int EDataType>
581 struct TTimeStampOnGaussMeshInitArray
583 typedef typename TL::TEnum2VTKArrayType<EDataType>::TResult TVTKDataArray;
584 typedef typename TL::TEnum2VTKBasicType<EDataType>::TResult TVTKBasicType;
585 typedef TTMeshValue<TVTKBasicType> TMeshValue;
586 typedef MED::SharedPtr<TMeshValue> TMeshValuePtr;
588 typedef TDataArrayHolder<EDataType> TTDataArrayHolder;
589 typedef MED::SharedPtr<TTDataArrayHolder> PDataArrayHolder;
590 PDataArrayHolder myDataArrayHolder;
592 TTimeStampOnGaussMeshInitArray(const PDataArrayHolder& theDataArrayHolder):
593 myDataArrayHolder(theDataArrayHolder)
597 Execute(const PFieldImpl& theField,
598 const PValForTimeImpl& theValForTime)
600 vtkIdType aNbComp = theField->myNbComp;
601 vtkIdType aSize = std::max(vtkIdType(3), aNbComp);
602 TVector<TVTKBasicType> aDataValues(aSize);
604 const TGeom2MeshValue& aGeom2MeshValue = theValForTime->GetGeom2MeshValue();
606 PGaussMeshImpl aGaussMesh = theValForTime->myGaussMesh;
607 const TGeom2GaussSubMesh& aGeom2GaussSubMesh = aGaussMesh->myGeom2GaussSubMesh;
608 TGeom2GaussSubMesh::const_iterator anIter = aGeom2GaussSubMesh.begin();
609 for(int aTupleId = 0; anIter != aGeom2GaussSubMesh.end(); anIter++){
610 EGeometry aEGeom = anIter->first;
612 PGaussSubMeshImpl aGaussSubMesh = anIter->second;
613 if(!aGaussSubMesh->myIsDone)
616 TGeom2MeshValue::const_iterator anIter2 = aGeom2MeshValue.find(aEGeom);
617 if(anIter2 == aGeom2MeshValue.end()){
618 EXCEPTION(std::runtime_error,
619 "TTimeStampOnGaussMeshInitArray >> Can't find values for corresponding Gauss Points SubMesh");
621 TMeshValuePtr aMeshValue = anIter2->second;
622 vtkIdType aNbGauss = aMeshValue->GetNbGauss();
623 vtkIdType aNbElem = aMeshValue->GetNbElem();
628 const TPointCoords& aCoords = aGaussSubMesh->myPointCoords;
631 "- aEGeom = "<<aEGeom<<
632 "; aNbElem = "<<aNbElem<<
633 "; aNbGauss = "<<aNbGauss<<
634 "; aCoords.GetNbPoints() = "<<aCoords.GetNbPoints()<<
637 if(aCoords.GetNbPoints() == aNbElem*aNbGauss){
638 for(int iElem = 0; iElem < aNbElem; iElem++){
639 typename TMeshValue::TCValueSliceArr aValueSliceArr = aMeshValue->GetGaussValueSliceArr(iElem);
640 for(int iGauss = 0; iGauss < aNbGauss; iGauss++, aTupleId++){
641 const typename TMeshValue::TCValueSlice& aValueSlice = aValueSliceArr[iGauss];
642 for(int iComp = 0; iComp < aNbComp; iComp++){
643 aDataValues[iComp] = aValueSlice[iComp];
645 this->myDataArrayHolder->SetTuple(aTupleId, &aDataValues[0]);
649 for(int iElem = 0; iElem < aNbElem; iElem++, aTupleId++){
650 typename TMeshValue::TCValueSliceArr aValueSliceArr = aMeshValue->GetCompValueSliceArr(iElem);
651 for(int iComp = 0; iComp < aNbComp; iComp++){
652 const typename TMeshValue::TCValueSlice& aValueSlice = aValueSliceArr[iComp];
653 aDataValues[iComp] = TVTKBasicType();
654 for(int iGauss = 0; iGauss < aNbGauss; iGauss++){
655 aDataValues[iComp] += aValueSlice[iGauss];
657 aDataValues[iComp] /= aNbGauss;
659 this->myDataArrayHolder->SetTuple(aTupleId, &aDataValues[0]);
667 template<int EDataType>
669 InitTimeStampOnGaussMesh(const PPolyData& theSource,
670 const PFieldImpl& theField,
671 const PValForTimeImpl& theValForTime)
673 vtkIdType aNbTuples = theSource->GetNumberOfPoints();
674 std::string aFieldName = VISU::GenerateFieldName(theField, theValForTime);
676 vtkDataSetAttributes* aDataSetAttributes = theSource->GetPointData();
678 typedef typename TL::TEnum2VTKArrayType<EDataType>::TResult TVTKDataArray;
679 TVTKDataArray *aSelectedDataArray = TVTKDataArray::New();
680 vtkIdType aNbComp = theField->myNbComp;
683 aSelectedDataArray->SetNumberOfComponents(1);
684 aDataSetAttributes->SetScalars(aSelectedDataArray);
687 aSelectedDataArray->SetNumberOfComponents(3);
688 aDataSetAttributes->SetVectors(aSelectedDataArray);
690 aSelectedDataArray->SetNumberOfTuples(aNbTuples);
691 aSelectedDataArray->SetName(aFieldName.c_str());
693 TVTKDataArray *aFullDataArray = TVTKDataArray::New();
694 aFullDataArray->SetNumberOfComponents(aNbComp);
695 aFullDataArray->SetNumberOfTuples(aNbTuples);
696 aFullDataArray->SetName("VISU_FIELD");
697 aDataSetAttributes->AddArray(aFullDataArray);
699 INITMSG(MYDEBUG,"InitTimeStampOnGaussMesh "<<
700 "- aNbTuples = "<<aNbTuples<<
701 "; aNbComp = "<<aNbComp<<
703 TTimerLog aTimerLog(MYDEBUG,"InitTimeStampOnGaussMesh");
705 const TGeom2MeshValue& aGeom2MeshValue = theValForTime->GetGeom2MeshValue();
706 typedef typename TL::TEnum2VTKBasicType<EDataType>::TResult TVTKBasicType;
707 typedef TTMeshValue<TVTKBasicType> TMeshValue;
708 typedef MED::SharedPtr<TMeshValue> TMeshValuePtr;
710 typedef TDataArrayHolder<EDataType> TTDataArrayHolder;
711 typedef MED::SharedPtr<TTDataArrayHolder> PDataArrayHolder;
713 TMeshValuePtr aMeshValue = theValForTime->GetFirstMeshValue();
714 if(aGeom2MeshValue.size() == 1){
715 aFullDataArray->SetVoidArray(aMeshValue->GetPointer(),
718 INITMSG(MYDEBUG,"InitTimeStampOnGaussMesh - aFullDataArray->SetVoidArray()"<<std::endl);
719 if(aNbComp == 1 || aNbComp == 3){
720 aSelectedDataArray->SetVoidArray(aMeshValue->GetPointer(),
723 INITMSG(MYDEBUG,"InitTimeStampOnGaussMesh - aSelectedDataArray->SetVoidArray()"<<std::endl);
725 PDataArrayHolder aDataArrayHolder(new TTDataArrayHolder(aSelectedDataArray));
726 TTimeStampOnGaussMeshInitArray<EDataType>(aDataArrayHolder).Execute(theField, theValForTime);
729 typedef TDataArrayHolder2<EDataType> TTDataArrayHolder2;
730 PDataArrayHolder aDataArrayHolder(new TTDataArrayHolder2(aSelectedDataArray, aFullDataArray));
731 TTimeStampOnGaussMeshInitArray<EDataType>(aDataArrayHolder).Execute(theField, theValForTime);
734 aSelectedDataArray->Delete();
735 aFullDataArray->Delete();
739 //----------------------------------------------------------------------------
741 void InitMed2VisuArray(std::vector<int>& anArray, EGeometry aEGeom){
743 #if !(defined(VTK_QUADRATIC_EDGE) && defined(VISU_USE_VTK_QUADRATIC)) && defined(VISU_ENABLE_QUADRATIC)
751 #if !(defined(VTK_QUADRATIC_TRIANGLE) && defined(VISU_USE_VTK_QUADRATIC)) && defined(VISU_ENABLE_QUADRATIC)
763 #if !(defined(VTK_QUADRATIC_QUAD) && defined(VISU_USE_VTK_QUADRATIC)) && defined(VISU_ENABLE_QUADRATIC)
789 #if (defined(VTK_QUADRATIC_TETRA) && defined(VISU_USE_VTK_QUADRATIC)) && defined(VISU_ENABLE_QUADRATIC)
806 #if (defined(VTK_QUADRATIC_PYRAMID) && defined(VISU_USE_VTK_QUADRATIC)) && defined(VISU_ENABLE_QUADRATIC)
826 for(int i=0;i<anArray.size();i++){