1 // Copyright (C) 2010-2012 CEA/DEN, EDF R&D
3 // This library is free software; you can redistribute it and/or
4 // modify it under the terms of the GNU Lesser General Public
5 // License as published by the Free Software Foundation; either
6 // version 2.1 of the License.
8 // This library is distributed in the hope that it will be useful,
9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
11 // Lesser General Public License for more details.
13 // You should have received a copy of the GNU Lesser General Public
14 // License along with this library; if not, write to the Free Software
15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
20 #include "vtkMedReader.h"
22 #include "vtkMedFile.h"
23 #include "vtkMedDriver.h"
24 #include "vtkMedUtilities.h"
25 #include "vtkMedFamily.h"
26 #include "vtkMedGroup.h"
27 #include "vtkMedMesh.h"
28 #include "vtkMedUnstructuredGrid.h"
29 #include "vtkMedCurvilinearGrid.h"
30 #include "vtkMedRegularGrid.h"
31 #include "vtkMedEntityArray.h"
32 #include "vtkMedField.h"
33 #include "vtkMedFieldStep.h"
34 #include "vtkMedFieldOverEntity.h"
35 #include "vtkMedProfile.h"
36 #include "vtkMedIntArray.h"
37 #include "vtkMedLocalization.h"
38 #include "vtkMedFamilyOnEntity.h"
39 #include "vtkMedFieldOnProfile.h"
40 #include "vtkMedSelection.h"
41 #include "vtkMedFamilyOnEntityOnProfile.h"
42 #include "vtkMedLink.h"
43 #include "vtkMedInterpolation.h"
44 #include "vtkMedStructElement.h"
45 #include "vtkMedConstantAttribute.h"
47 #include "vtkObjectFactory.h"
48 #include "vtkMutableDirectedGraph.h"
49 #include "vtkStringArray.h"
50 #include "vtkDataSetAttributes.h"
51 #include "vtkUnsignedCharArray.h"
52 #include "vtkInformation.h"
53 #include "vtkInformationVector.h"
54 #include "vtkSmartPointer.h"
55 #include "vtkVariantArray.h"
56 #include "vtkMultiBlockDataSet.h"
57 #include "vtkExecutive.h"
58 #include "vtkStreamingDemandDrivenPipeline.h"
59 #include "vtkUnstructuredGrid.h"
61 #include "vtkPointData.h"
62 #include "vtkCellData.h"
63 #include "vtkFieldData.h"
64 #include "vtkInformationQuadratureSchemeDefinitionVectorKey.h"
65 #include "vtkQuadratureSchemeDefinition.h"
66 #include "vtkCellType.h"
67 #include "vtkCellArray.h"
68 #include "vtkDoubleArray.h"
69 #include "vtkConfigure.h"
70 #include "vtkMultiProcessController.h"
71 #include "vtkCommunicator.h"
73 #include "vtkSMDoubleVectorProperty.h"
74 #include "vtkInformationDataObjectKey.h"
88 vtkSmartPointer<vtkDataArray> DataArray;
89 vtkSmartPointer<vtkDataArray> Vectors;
90 vtkSmartPointer<vtkIdTypeArray> QuadratureIndexArray;
93 class vtkMedListOfFieldSteps : public std::list<vtkMedFieldStep*>
96 typedef int LocalizationKey;
98 class vtkMedReader::vtkMedReaderInternal
102 int CurrentPieceNumber;
104 double UpdateTimeStep;
105 vtkTimeStamp FileNameMTime;
106 vtkTimeStamp MetaDataMTime;
107 vtkTimeStamp GroupSelectionMTime;
108 vtkTimeStamp FamilySelectionMTime;
110 int RealAnimationMode;
111 vtkMedSelection* Families;
113 // this stores the aggregation of all compute steps from
114 // both meshes and fields.
115 std::map<med_float, std::set<med_int> > GlobalComputeStep;
117 // Store the vtkMutableDirectedGraph that represents links between family, groups and cell types
118 vtkMutableDirectedGraph* SIL;
120 // this map is used to keep clean data sets in the cache, without any field.
121 // for each support, store the vtkDataSet
122 map<vtkMedFamilyOnEntityOnProfile*, vtkSmartPointer<vtkDataSet> > DataSetCache;
124 // this is the current dataset for the given support.
125 map<vtkMedFamilyOnEntityOnProfile*, vtkDataSet*> CurrentDataSet;
127 // for each support, cache the VTK arrays that correspond to a given field at the given step.
128 map<vtkMedFamilyOnEntityOnProfile*, map<vtkMedFieldOnProfile*, VTKField> > FieldCache;
129 //map<vtkMedFamilyOnEntity*, map<vtkMedFieldOnProfile*, bool> > FieldMatchCache;
131 // This list keep tracks of all the currently selected supports
132 set<vtkMedFamilyOnEntityOnProfile*> UsedSupports;
134 // this map keeps for each support, the quadrature offset array so that
135 // different fields on the same support can use
136 // the same offset array, provided they use the same gauss points definitions
137 map<vtkMedFamilyOnEntityOnProfile*,
138 map<LocalizationKey, vtkSmartPointer<vtkIdTypeArray> > >
139 QuadratureOffsetCache;
141 map<vtkMedFamilyOnEntityOnProfile*,
142 map<vtkMedFieldOnProfile*, LocalizationKey> > QuadOffsetKey;
144 std::map<std::string, vtkSmartPointer<vtkMedFile> > MedFiles;
146 vtkMedReaderInternal()
148 this->SIL=vtkMutableDirectedGraph::New();
149 this->SILUpdateStamp=-1;
150 this->RealAnimationMode=vtkMedReader::PhysicalTime;
151 this->Families=vtkMedSelection::New();
152 this->FamilySelectionMTime.Modified();
153 this->GroupSelectionMTime.Modified();
155 ~vtkMedReaderInternal()
158 this->Families->Delete();
161 void ClearSupportCache(vtkMedFamilyOnEntityOnProfile* foep)
163 //this->Med2VTKPointIndex.erase(foep);
164 this->QuadratureOffsetCache.erase(foep);
165 //this->FieldMatchCache.erase(foe);
166 this->FieldCache.erase(foep);
167 this->CurrentDataSet.erase(foep);
168 this->DataSetCache.erase(foep);
171 vtkIdType GetChild(vtkIdType parent, const vtkStdString childName)
173 vtkStringArray* names=vtkStringArray::SafeDownCast(
174 this->SIL->GetVertexData()->GetArray("Names"));
177 vtkIdType nedges=this->SIL->GetOutDegree(parent);
178 for(vtkIdType id=0; id<nedges; id++)
180 vtkOutEdgeType edge=this->SIL->GetOutEdge(parent, id);
181 if(names->GetValue(edge.Target)==childName)
187 vtkIdType GetGroupId(const char* key)
189 vtkstd::string meshname, celltypename, groupname;
190 vtkMedUtilities::SplitGroupKey(key, meshname, celltypename, groupname);
191 vtkIdType root=GetChild(0, "SIL");
194 vtkIdType mesh=GetChild(root, meshname);
197 vtkIdType type=GetChild(mesh, celltypename);
200 return GetChild(type, groupname);
206 vtkCxxRevisionMacro(vtkMedReader, "$Revision$");
207 vtkStandardNewMacro(vtkMedReader);
209 vtkMedReader::vtkMedReader()
212 this->SetNumberOfInputPorts(0);
213 this->PointFields=vtkMedSelection::New();
214 this->CellFields=vtkMedSelection::New();
215 this->QuadratureFields=vtkMedSelection::New();
216 this->ElnoFields=vtkMedSelection::New();
217 this->Entities=vtkMedSelection::New();
218 this->Groups=vtkMedSelection::New();
219 this->Frequencies=vtkMedSelection::New();
220 this->AnimationMode=Default;
221 this->TimeIndexForIterations=0;
222 this->CacheStrategy=CacheGeometry;
223 this->Internal=new vtkMedReaderInternal;
224 this->TimePrecision=0.00001;
225 this->AvailableTimes=vtkDoubleArray::New();
226 this->GenerateVectors = 0;
229 vtkMedReader::~vtkMedReader()
231 this->SetFileName(NULL);
232 this->PointFields->Delete();
233 this->CellFields->Delete();
234 this->QuadratureFields->Delete();
235 this->ElnoFields->Delete();
236 this->Entities->Delete();
237 this->Groups->Delete();
238 this->Frequencies->Delete();
239 delete this->Internal;
240 this->AvailableTimes->Delete();
243 int vtkMedReader::GetSILUpdateStamp()
245 return this->Internal->SILUpdateStamp;
248 void vtkMedReader::SetFileName(const char* fname)
251 if(fname==this->FileName)
253 if(fname&&this->FileName&&!strcmp(fname, this->FileName))
257 delete[] this->FileName;
260 size_t fnl=strlen(fname)+1;
261 char* dst=new char[fnl];
262 const char* src=fname;
277 this->Internal->MedFiles.clear();
278 this->Internal->FileNameMTime.Modified();
282 int vtkMedReader::CanReadFile(const char* fname)
284 // the factory give a driver only when it can read the file version,
285 // or it returns a NULL pointer.
286 vtkSmartPointer<vtkMedFile> file=vtkSmartPointer<vtkMedFile>::New();
287 file->SetFileName(fname);
289 if(!file->CreateDriver())
295 int vtkMedReader::RequestInformation(vtkInformation *request,
296 vtkInformationVector **inputVector, vtkInformationVector *outputVector)
298 vtkInformation* outInfo = outputVector->GetInformationObject(0);
299 outInfo->Set(vtkStreamingDemandDrivenPipeline::MAXIMUM_NUMBER_OF_PIECES(),-1);
301 if(this->Internal->MetaDataMTime <= this->Internal->FileNameMTime)
303 this->ClearCaches(Initialize);
305 vtkMedFile* file=vtkMedFile::New();
306 file->SetFileName(this->FileName);
307 this->Internal->MedFiles[this->FileName] = file;
310 std::list<vtkMedFile*> file_stack;
311 file_stack.push_back(file);
313 // if the file cannot create a driver, this means that the filename is not
314 // valid, or that I do not know how to read this file.
315 while (file_stack.size() > 0)
317 vtkMedFile* file = file_stack.front();
318 file_stack.pop_front();
319 // This reads information from disk
320 file->ReadInformation();
322 // add all files linked to in the current file to the files to read.
323 for (int linkid=0; linkid<file->GetNumberOfLink(); linkid++)
325 vtkMedLink* link = file->GetLink(linkid);
326 const char* filename = link->GetFullLink(file->GetFileName());
327 if(this->Internal->MedFiles.find(filename) == this->Internal->MedFiles.end())
329 vtkMedFile* newfile = vtkMedFile::New();
330 newfile->SetFileName(filename);
331 this->Internal->MedFiles[filename] = newfile;
332 file_stack.push_back(newfile);
338 // This computes some meta information, like which field use which
339 // support, but do not read large data from disk.
342 // This computes the initial global id of each cell type.
343 this->InitializeCellGlobalIds();
345 this->ClearSelections();
347 this->BuildSIL(this->Internal->SIL);
348 this->Internal->SILUpdateStamp++;
350 this->GatherComputeSteps();
352 this->Internal->MetaDataMTime.Modified();
355 outInfo->Set(vtkDataObject::SIL(), this->Internal->SIL);
356 request->AppendUnique(vtkExecutive::KEYS_TO_COPY(),
357 vtkDataObject::SIL());
358 request->AppendUnique(vtkExecutive::KEYS_TO_COPY(),
359 vtkMedUtilities::BLOCK_NAME());
360 this->AdvertiseTime(outInfo);
364 int vtkMedReader::RequestData(vtkInformation *request,
365 vtkInformationVector **inputVector, vtkInformationVector *outputVector)
367 if(this->FileName==NULL)
369 vtkWarningMacro( << "FileName must be set and meta data loaded");
373 vtkInformation *info=outputVector->GetInformationObject(0);
375 vtkMultiBlockDataSet *output=vtkMultiBlockDataSet::SafeDownCast(info->Get(
376 vtkDataObject::DATA_OBJECT()));
378 if (info->Has(vtkStreamingDemandDrivenPipeline::UPDATE_NUMBER_OF_PIECES()))
380 this->Internal->NumberOfPieces=info->Get(
381 vtkStreamingDemandDrivenPipeline::UPDATE_NUMBER_OF_PIECES());
385 vtkMultiProcessController* controller =
386 vtkMultiProcessController::GetGlobalController();
389 this->Internal->NumberOfPieces=controller->GetNumberOfProcesses();
393 this->Internal->NumberOfPieces = 1;
396 if (info->Has(vtkStreamingDemandDrivenPipeline::UPDATE_PIECE_NUMBER()))
398 this->Internal->CurrentPieceNumber=info->Get(
399 vtkStreamingDemandDrivenPipeline::UPDATE_PIECE_NUMBER());
403 this->Internal->CurrentPieceNumber=0;
404 vtkMultiProcessController* controller =
405 vtkMultiProcessController::GetGlobalController();
408 this->Internal->CurrentPieceNumber= controller->GetLocalProcessId();
412 this->Internal->CurrentPieceNumber=0;
417 vtkStreamingDemandDrivenPipeline::UPDATE_NUMBER_OF_GHOST_LEVELS()))
419 this->Internal->GhostLevel=info->Get(
420 vtkStreamingDemandDrivenPipeline::UPDATE_NUMBER_OF_GHOST_LEVELS());
424 this->Internal->GhostLevel=0;
427 if (info->Has(vtkStreamingDemandDrivenPipeline::UPDATE_TIME_STEPS()))
429 this->Internal->UpdateTimeStep=info->Get(
430 vtkStreamingDemandDrivenPipeline::UPDATE_TIME_STEPS())[0];
434 this->Internal->UpdateTimeStep=0;
437 this->InitializeParallelRead();
438 output->Initialize();
440 this->ChooseRealAnimationMode();
442 std::list<vtkMedDriver::FileOpen> openlist;
443 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
444 fileit = this->Internal->MedFiles.begin();
445 while(fileit != this->Internal->MedFiles.end())
447 vtkMedFile* file = fileit->second;
448 openlist.push_back(vtkMedDriver::FileOpen(file->GetMedDriver()));
452 // clear the dataset cache of unneeded geometry
453 this->ClearCaches(StartRequest);
455 // This call create the vtkMedSupports, but do not create the corresponding vtkDataSet;
456 this->CreateMedSupports();
457 this->ClearCaches(AfterCreateMedSupports);
458 // This call creates the actual vtkDataSet that corresponds to each support
461 int maxprogress=2*this->Internal->UsedSupports.size();
464 for(set<vtkMedFamilyOnEntityOnProfile*>::iterator it=
465 this->Internal->UsedSupports.begin(); it
466 !=this->Internal->UsedSupports.end(); it++)
469 vtkMedFamilyOnEntityOnProfile* foep = *it;
470 sstr<<"Support : "<<vtkMedUtilities::SimplifyName(
471 foep->GetFamilyOnEntity()->GetFamily()->GetName());
472 this->SetProgressText(sstr.str().c_str());
473 int doBuildSupportField = 1;
475 this->BuildVTKSupport(foep, doBuildSupportField);
476 this->UpdateProgress((float) progress/((float) maxprogress-1));
481 this->ClearCaches(EndBuildVTKSupports);
482 // This call maps the fields to the supports
483 for(set<vtkMedFamilyOnEntityOnProfile*>::iterator it=
484 this->Internal->UsedSupports.begin(); it
485 !=this->Internal->UsedSupports.end(); it++)
487 vtkMedFamilyOnEntityOnProfile* foep = *it;
488 if((foep->GetValid() == 0) && (this->Internal->NumberOfPieces == 1))
491 sstr<<"Loading fields on "<<vtkMedUtilities::SimplifyName(
492 foep->GetFamilyOnEntity()->GetFamily()->GetName());
493 this->SetProgressText(sstr.str().c_str());
495 this->MapFieldsOnSupport(*it, doMapField);
496 this->UpdateProgress((float) progress/((float) maxprogress-1));
501 // This call clean up caches (what is actually done depends of the CacheStrategy)
502 this->ClearCaches(EndRequest);
506 void vtkMedReader::InitializeCellGlobalIds()
508 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
509 fileit = this->Internal->MedFiles.begin();
510 while(fileit != this->Internal->MedFiles.end())
512 vtkMedFile* file = fileit->second;
514 for(int m=0; m<file->GetNumberOfMesh(); m++)
516 vtkMedMesh* mesh = file->GetMesh(m);
517 med_int nstep = mesh->GetNumberOfGridStep();
518 for(med_int stepid=0; stepid<nstep; stepid++)
520 vtkMedGrid* grid = mesh->GetGridStep(stepid);
521 grid->InitializeCellGlobalIds();
527 // Method to create the filters for the MED parallel read functions
528 // It is defined here as we have all information for initialization
529 void vtkMedReader::InitializeParallelRead()
531 // If there is only one process for reading no need to enter here
532 if (this->Internal->NumberOfPieces <= 1)
537 // FIRST: Generate filters for the cells
539 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
540 meshfit = this->Internal->MedFiles.begin();
541 while(meshfit != this->Internal->MedFiles.end())
543 vtkMedFile* meshfile = meshfit->second;
545 med_idt pFileID = meshfile->GetMedDriver()->GetParallelFileId();
547 for(int mid=0; mid<meshfile->GetNumberOfMesh(); mid++)
549 vtkMedMesh* mesh = meshfile->GetMesh(mid);
550 for(int gid=0; gid<mesh->GetNumberOfGridStep(); gid++)
552 vtkMedGrid* grid = mesh->GetGridStep(gid);
553 // read point family data and create EntityArrays
555 for(int eid=0; eid < grid->GetNumberOfEntityArray(); eid++)
557 vtkMedEntityArray* array = grid->GetEntityArray(eid);
559 // Next continue is to avoid to create filters for the
560 // points, at the moment we charge the points in all nodes
561 if (array->GetEntity().GeometryType == MED_POINT1) // !MED_NODE
564 med_int nbofconstituentpervalue = vtkMedUtilities::GetNumberOfNodes(
565 array->GetEntity().GeometryType);
566 if (nbofconstituentpervalue == -1)
567 vtkErrorMacro("Still not implemented for MED_POLYGON and MED_POLYHEDRON"); // à gerer
569 // Calculating block sizes
570 int nEntity = array->GetNumberOfEntity();
571 int block_size = ( nEntity / this->Internal->NumberOfPieces );
572 med_size start = block_size * this->Internal->CurrentPieceNumber + 1;
573 med_size stride = block_size;
575 med_size blocksize = block_size;
576 med_size lastblocksize = (nEntity % this->Internal->NumberOfPieces);
577 if ((this->Internal->NumberOfPieces ==
578 this->Internal->CurrentPieceNumber+1) && (lastblocksize != 0))
580 blocksize += lastblocksize;
581 stride += lastblocksize;
585 vtkMedFilter *filter = vtkMedFilter::New();
586 filter->SetFilterSizes( start, stride, count, blocksize, lastblocksize );
587 array->SetFilter(filter);
593 // SECOND: Filters for the Fields
595 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
597 // link the FieldOnProfile with the profiles
598 fieldfileit = this->Internal->MedFiles.begin();
599 while(fieldfileit != this->Internal->MedFiles.end())
601 vtkMedFile* fieldfile = fieldfileit->second;
603 med_idt pFileID = fieldfile->GetMedDriver()->GetParallelFileId();
605 for(int fid=0; fid<fieldfile->GetNumberOfField(); fid++)
607 vtkMedField* field = fieldfile->GetField(fid);
609 if (field->GetFieldType() == vtkMedField::CellField)
611 for(int sid = 0; sid< field->GetNumberOfFieldStep(); sid++)
613 vtkMedFieldStep* step = field->GetFieldStep(sid);
615 for(int foeid = 0; foeid < step->GetNumberOfFieldOverEntity(); foeid++)
616 // TODO : seul le premier pas de temps est dispo au debut
618 vtkMedFieldOverEntity* fieldOverEntity = step->GetFieldOverEntity(foeid);
620 for(int fopid = 0; fopid < fieldOverEntity->GetNumberOfFieldOnProfile(); fopid++)
622 vtkMedFieldOnProfile* fop = fieldOverEntity->GetFieldOnProfile(fopid);
623 // Here implement the filters as before:
624 // 1- Modify vtkMedFieldOnProfile to contain a filter
625 // 2- Create the filters here only if they are on CELLs (use GetFieldType)
626 med_int nbofconstituentpervalue = field->GetNumberOfComponent();
628 int nVectors = fop->GetNumberOfValues();
630 int block_size = ( nVectors / this->Internal->NumberOfPieces );
631 int start = block_size * this->Internal->CurrentPieceNumber + 1;
632 int stride = block_size;
634 int blocksize = block_size;
635 int lastblocksize = (nVectors % this->Internal->NumberOfPieces);
636 if ((this->Internal->NumberOfPieces ==
637 this->Internal->CurrentPieceNumber+1) && (lastblocksize != 0))
639 blocksize += lastblocksize;
640 stride += lastblocksize;
644 vtkMedFilter *filter = vtkMedFilter::New();
645 filter->SetFilterSizes( start, stride, count, blocksize, lastblocksize );
646 fop->SetFilter(filter);
656 void vtkMedReader::LinkMedInfo()
658 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
660 // link the FieldOnProfile with the profiles
661 fieldfileit = this->Internal->MedFiles.begin();
662 while(fieldfileit != this->Internal->MedFiles.end())
664 vtkMedFile* fieldfile = fieldfileit->second;
667 for(int fid=0; fid<fieldfile->GetNumberOfField(); fid++)
669 vtkMedField* field = fieldfile->GetField(fid);
671 for(int sid = 0; sid< field->GetNumberOfFieldStep(); sid++)
673 vtkMedFieldStep* step = field->GetFieldStep(sid);
675 for(int foeid = 0; foeid < step->GetNumberOfFieldOverEntity(); foeid++)
677 vtkMedFieldOverEntity* fieldOverEntity = step->GetFieldOverEntity(foeid);
679 for(int fopid = 0; fopid < fieldOverEntity->GetNumberOfFieldOnProfile(); fopid++)
681 vtkMedFieldOnProfile* fop = fieldOverEntity->GetFieldOnProfile(fopid);
683 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
684 profilefileit = this->Internal->MedFiles.begin();
685 while(profilefileit != this->Internal->MedFiles.end() && fop->GetProfile() == NULL)
687 vtkMedFile* profilefile = profilefileit->second;
690 for(int pid = 0; pid < profilefile->GetNumberOfProfile(); pid++)
692 vtkMedProfile *profile = profilefile->GetProfile(pid);
693 if(strcmp(profile->GetName(), fop->GetProfileName()) == 0)
695 fop->SetProfile(profile);
706 // first, add a familyOnEntityOnProfile to all FamilyOnEntity with a NULL
707 // profile. This is used if no field is mapped to this support.
708 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
709 meshfit = this->Internal->MedFiles.begin();
710 while(meshfit != this->Internal->MedFiles.end())
712 vtkMedFile* meshfile = meshfit->second;
715 for(int mid=0; mid<meshfile->GetNumberOfMesh(); mid++)
717 vtkMedMesh* mesh = meshfile->GetMesh(mid);
719 for(int gid=0; gid<mesh->GetNumberOfGridStep(); gid++)
721 vtkMedGrid* grid = mesh->GetGridStep(gid);
722 // read point family data and create EntityArrays
724 for(int eid=0; eid < grid->GetNumberOfEntityArray(); eid++)
726 vtkMedEntityArray* array = grid->GetEntityArray(eid);
728 for(int fid=0; fid < array->GetNumberOfFamilyOnEntity(); fid++)
730 vtkMedFamilyOnEntity* foe = array->GetFamilyOnEntity(fid);
731 if(foe->GetFamilyOnEntityOnProfile((vtkMedProfile*)NULL) == NULL)
733 vtkMedFamilyOnEntityOnProfile* foep =
734 vtkMedFamilyOnEntityOnProfile::New();
735 foep->SetFamilyOnEntity(foe);
736 foep->SetProfile(NULL);
737 foe->AddFamilyOnEntityOnProfile(foep);
746 fieldfileit = this->Internal->MedFiles.begin();
747 while(fieldfileit != this->Internal->MedFiles.end())
749 vtkMedFile* fieldfile = fieldfileit->second;
752 for(int fieldid=0; fieldid < fieldfile->GetNumberOfField(); fieldid++)
754 vtkMedField* field = fieldfile->GetField(fieldid);
756 for(int fstepid=0; fstepid < field->GetNumberOfFieldStep(); fstepid++)
758 vtkMedFieldStep* step = field->GetFieldStep(fstepid);
760 vtkMedComputeStep meshcs = step->GetMeshComputeStep();
762 for(int foeid=0; foeid<step->GetNumberOfFieldOverEntity() ;foeid++)
764 vtkMedFieldOverEntity* fieldOverEntity = step->GetFieldOverEntity(foeid);
765 const vtkMedEntity& fieldentity = fieldOverEntity->GetEntity();
768 fopid < fieldOverEntity->GetNumberOfFieldOnProfile(); fopid++)
770 vtkMedFieldOnProfile* fop =
771 fieldOverEntity->GetFieldOnProfile(fopid);
773 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
774 meshfileit = this->Internal->MedFiles.begin();
775 while(meshfileit != this->Internal->MedFiles.end())
777 vtkMedFile* meshfile = meshfileit->second;
780 if(field->GetLocal() == 1 && (meshfile != fieldfile))
783 for(int mid=0; mid<meshfile->GetNumberOfMesh(); mid++)
785 vtkMedMesh* mesh = meshfile->GetMesh(mid);
787 // the field must be on this mesh.
788 if(strcmp(mesh->GetName(),
789 field->GetMeshName()) != 0)
792 vtkMedGrid* grid = mesh->GetGridStep(meshcs);
795 vtkErrorMacro("the field " << field->GetName()
796 << " at step iteration:"
797 << step->GetComputeStep().IterationIt
799 << step->GetComputeStep().TimeIt
800 << " uses mesh at step "
801 << meshcs.TimeIt << " " << meshcs.IterationIt
802 << "which does not exists!");
806 for(int eid=0; eid < grid->GetNumberOfEntityArray(); eid++)
808 vtkMedEntityArray* array = grid->GetEntityArray(eid);
810 // if the support is on points,
811 // the field must also be on points
812 if(array->GetEntity().EntityType == MED_NODE &&
813 fieldentity.EntityType != MED_NODE)
816 if(array->GetEntity().EntityType != MED_NODE &&
817 fieldentity.EntityType == MED_NODE)
820 // for fields not on points, the geometry type
821 // of the support must match
822 if(array->GetEntity().EntityType != MED_NODE &&
823 array->GetEntity().GeometryType != fieldentity.GeometryType)
826 for(int fid = 0; fid < array->GetNumberOfFamilyOnEntity(); fid++)
828 vtkMedFamilyOnEntity* foe = array->GetFamilyOnEntity(fid);
829 if(foe->GetFamilyOnEntityOnProfile(fop->GetProfile()) == NULL)
831 vtkMedFamilyOnEntityOnProfile* foep =
832 vtkMedFamilyOnEntityOnProfile::New();
833 foep->SetProfile(fop->GetProfile());
834 foep->SetFamilyOnEntity(foe);
835 foe->AddFamilyOnEntityOnProfile(foep);
838 // also add the family on entity with no profile.
839 if(foe->GetFamilyOnEntityOnProfile((vtkMedProfile*)NULL) == NULL)
841 vtkMedFamilyOnEntityOnProfile* foep =
842 vtkMedFamilyOnEntityOnProfile::New();
843 foep->SetProfile(NULL);
844 foep->SetFamilyOnEntity(foe);
845 foe->AddFamilyOnEntityOnProfile(foep);
858 // Now, link localizations and interpolations
859 fieldfileit = this->Internal->MedFiles.begin();
860 while(fieldfileit != this->Internal->MedFiles.end())
862 vtkMedFile* fieldfile = fieldfileit->second;
865 for(int locid = 0; locid < fieldfile->GetNumberOfLocalization(); locid++)
867 vtkMedLocalization* loc = fieldfile->GetLocalization(locid);
869 for(int fid = 0; fid < fieldfile->GetNumberOfField() &&
870 loc->GetInterpolation() == NULL; fid++)
872 vtkMedField* field = fieldfile->GetField(fid);
873 for(int interpid = 0; interpid < field->GetNumberOfInterpolation();
876 vtkMedInterpolation* interp = field->GetInterpolation(interpid);
877 if(strcmp(loc->GetInterpolationName(),
878 interp->GetName()) == 0)
880 loc->SetInterpolation(interp);
888 // now that the interpolation is set, build the shape functions.
889 fieldfileit = this->Internal->MedFiles.begin();
890 while(fieldfileit != this->Internal->MedFiles.end())
892 vtkMedFile* fieldfile = fieldfileit->second;
895 for(int locid = 0; locid < fieldfile->GetNumberOfLocalization(); locid++)
897 vtkMedLocalization* loc = fieldfile->GetLocalization(locid);
898 loc->BuildShapeFunction();
902 // set the supportmesh pointer in the structural element
903 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
904 fileit = this->Internal->MedFiles.begin();
905 while(fileit != this->Internal->MedFiles.end())
907 vtkMedFile* file = fileit->second;
910 for(int structelemit = 0;
911 structelemit < file->GetNumberOfStructElement();
914 vtkMedStructElement* structElem =
915 file->GetStructElement(structelemit);
917 for(int supportmeshit = 0;
918 supportmeshit < file->GetNumberOfSupportMesh();
921 vtkMedMesh* supportMesh =
922 file->GetSupportMesh(supportmeshit);
924 if(strcmp(supportMesh->GetName(), structElem->GetName()) == 0 )
926 structElem->SetSupportMesh(supportMesh);
933 // set the pointer to the profile used by the constant attributes
934 fileit = this->Internal->MedFiles.begin();
935 while(fileit != this->Internal->MedFiles.end())
937 vtkMedFile* file = fileit->second;
940 for(int structelemit = 0;
941 structelemit < file->GetNumberOfStructElement();
944 vtkMedStructElement* structElem =
945 file->GetStructElement(structelemit);
947 for(int cstattit = 0; cstattit < structElem->GetNumberOfConstantAttribute(); cstattit++)
949 vtkMedConstantAttribute* cstatt = structElem->GetConstantAttribute(cstattit);
952 profit < file->GetNumberOfProfile();
955 vtkMedProfile* profile =
956 file->GetProfile(profit);
958 if(strcmp(profile->GetName(), cstatt->GetProfileName()) == 0 )
960 cstatt->SetProfile(profile);
968 meshfit = this->Internal->MedFiles.begin();
969 while(meshfit != this->Internal->MedFiles.end())
971 vtkMedFile* meshfile = meshfit->second;
974 for(int mid=0; mid<meshfile->GetNumberOfMesh(); mid++)
976 vtkMedMesh* mesh = meshfile->GetMesh(mid);
978 for(int gid=0; gid<mesh->GetNumberOfGridStep(); gid++)
980 vtkMedGrid* grid = mesh->GetGridStep(gid);
981 // read point family data and create EntityArrays
983 for(int eid=0; eid < grid->GetNumberOfEntityArray(); eid++)
985 vtkMedEntityArray* array = grid->GetEntityArray(eid);
986 if(array->GetEntity().EntityType != MED_STRUCT_ELEMENT)
989 for(int structelemit = 0; structelemit < meshfile->GetNumberOfStructElement(); structelemit++)
991 vtkMedStructElement* structelem = meshfile->GetStructElement(structelemit);
992 if(structelem->GetGeometryType() == array->GetEntity().GeometryType)
994 array->SetStructElement(structelem);
1004 int vtkMedReader::GetFrequencyArrayStatus(const char* name)
1006 return this->Frequencies->GetKeyStatus(name);
1009 void vtkMedReader::SetFrequencyArrayStatus(const char* name, int status)
1011 if(this->Frequencies->GetKeyStatus(name) == status)
1016 this->Frequencies->SetKeyStatus(name, status);
1021 int vtkMedReader::GetNumberOfFrequencyArrays()
1023 return this->Frequencies->GetNumberOfKey();
1026 const char* vtkMedReader::GetFrequencyArrayName(int index)
1028 return this->Frequencies->GetKey(index);
1033 bool operator()(pair<double, int> i, pair<double, int> j)
1035 if(i.first!=j.first)
1036 return (i.first<j.first);
1037 return i.second<j.second;
1041 vtkDoubleArray* vtkMedReader::GetAvailableTimes()
1043 this->AvailableTimes->Initialize();
1044 this->AvailableTimes->SetNumberOfComponents(1);
1046 std::set<std::string> newFrequencies;
1049 std::map<med_float, std::set<med_int> >::iterator it =
1050 this->Internal->GlobalComputeStep.begin();
1051 while(it != this->Internal->GlobalComputeStep.end())
1053 double time = it->first;
1054 this->AvailableTimes->InsertNextValue(time);
1055 string name = vtkMedUtilities::GetModeKey(tid, time, this->Internal->GlobalComputeStep.size()-1);
1056 this->Frequencies->AddKey(name.c_str());
1057 newFrequencies.insert(name);
1062 // now check if old frequencies have been removed
1063 for(int f = 0; f < this->Frequencies->GetNumberOfKey(); f++)
1065 const char* name = this->Frequencies->GetKey(f);
1066 if(newFrequencies.find(name) == newFrequencies.end())
1068 this->Frequencies->RemoveKeyByIndex(f);
1073 return this->AvailableTimes;
1076 void vtkMedReader::ChooseRealAnimationMode()
1078 if(this->AnimationMode!=Default)
1080 this->Internal->RealAnimationMode=this->AnimationMode;
1084 // if there is exactly one physical time and more than one iteration
1085 // set the animation mode to iteration, else default to physical time.
1086 if (this->Internal->GlobalComputeStep.size() == 1 &&
1087 this->Internal->GlobalComputeStep[0].size() > 1)
1089 this->Internal->RealAnimationMode=Iteration;
1093 this->Internal->RealAnimationMode=PhysicalTime;
1096 int vtkMedReader::GetEntityStatus(const vtkMedEntity& entity)
1098 if (entity.EntityType==MED_NODE)
1100 if(entity.EntityType == MED_DESCENDING_FACE
1101 || entity.EntityType == MED_DESCENDING_EDGE)
1104 return this->Entities->GetKeyStatus(vtkMedUtilities::EntityKey(entity).c_str());
1107 int vtkMedReader::GetFamilyStatus(vtkMedMesh* mesh, vtkMedFamily* family)
1112 if(this->Internal->GroupSelectionMTime > this->Internal->FamilySelectionMTime)
1114 this->SelectFamiliesFromGroups();
1117 int status = this->Internal->Families->GetKeyStatus(vtkMedUtilities::FamilyKey(
1118 mesh->GetName(), family->GetPointOrCell(),
1119 family->GetName()).c_str());
1124 int vtkMedReader::IsMeshSelected(vtkMedMesh* mesh)
1126 for(int fam=0; fam<mesh->GetNumberOfPointFamily(); fam++)
1128 if(this->GetFamilyStatus(mesh, mesh->GetPointFamily(fam))!=0)
1132 for(int fam=0; fam<mesh->GetNumberOfCellFamily(); fam++)
1134 if(this->GetFamilyStatus(mesh, mesh->GetCellFamily(fam))!=0)
1140 void vtkMedReader::GatherComputeSteps()
1142 this->Internal->GlobalComputeStep.clear();
1144 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
1145 fieldfileit = this->Internal->MedFiles.begin();
1146 while(fieldfileit != this->Internal->MedFiles.end())
1148 vtkMedFile* file = fieldfileit->second;
1151 // first loop over all fields to gather their compute steps
1152 for(int fieldId=0; fieldId<file->GetNumberOfField(); fieldId++)
1154 vtkMedField* field=file->GetField(fieldId);
1156 for(int stepId=0; stepId<field->GetNumberOfFieldStep(); stepId++)
1158 vtkMedFieldStep* step=field->GetFieldStep(stepId);
1159 const vtkMedComputeStep& cs = step->GetComputeStep();
1160 this->Internal->GlobalComputeStep[cs.TimeOrFrequency].insert(cs.IterationIt);
1164 // then loop over all meshes to gather their grid steps too.
1165 // for meshes, do not add the MED_UNDEF_DT time
1166 for(int meshId=0; meshId<file->GetNumberOfMesh(); meshId++)
1168 vtkMedMesh* mesh=file->GetMesh(meshId);
1170 for(int stepId=0; stepId<mesh->GetNumberOfGridStep(); stepId++)
1172 vtkMedGrid* grid=mesh->GetGridStep(stepId);
1173 const vtkMedComputeStep& cs = grid->GetComputeStep();
1174 if(cs.TimeOrFrequency != MED_UNDEF_DT || cs.TimeIt != MED_NO_DT)
1176 this->Internal->GlobalComputeStep[cs.TimeOrFrequency].insert(cs.IterationIt);
1181 if(this->Internal->GlobalComputeStep.size() == 0)
1183 this->Internal->GlobalComputeStep[MED_UNDEF_DT].insert(MED_NO_IT);
1187 int vtkMedReader::IsFieldSelected(vtkMedField* field)
1189 return this->IsPointFieldSelected(field)||this->IsCellFieldSelected(field)
1190 ||this->IsQuadratureFieldSelected(field) || this->IsElnoFieldSelected(field);
1193 int vtkMedReader::IsPointFieldSelected(vtkMedField* field)
1195 return field->GetFieldType()==vtkMedField::PointField
1196 &&this->GetPointFieldArrayStatus(vtkMedUtilities::SimplifyName(
1197 field->GetName()).c_str());
1200 int vtkMedReader::IsCellFieldSelected(vtkMedField* field)
1202 return field->GetFieldType()==vtkMedField::CellField
1203 &&this->GetCellFieldArrayStatus(vtkMedUtilities::SimplifyName(
1204 field->GetName()).c_str());
1207 vtkMedProfile* vtkMedReader::GetProfile(const char* pname)
1209 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
1210 fileit = this->Internal->MedFiles.begin();
1211 while(fileit != this->Internal->MedFiles.end())
1213 vtkMedFile* file = fileit->second;
1215 vtkMedProfile* profile = file->GetProfile(pname);
1222 vtkMedLocalization* vtkMedReader::GetLocalization(const char* lname)
1224 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
1225 fileit = this->Internal->MedFiles.begin();
1226 while(fileit != this->Internal->MedFiles.end())
1228 vtkMedFile* file = fileit->second;
1230 vtkMedLocalization* loc = file->GetLocalization(lname);
1238 int vtkMedReader::GetLocalizationKey(vtkMedFieldOnProfile* fop)
1240 vtkMedLocalization* def=this->GetLocalization(fop->GetLocalizationName());
1242 // This is not a quadrature field with explicit definition.
1243 // There are two possible cases : either the intergration point is
1244 // at the center of the cell
1245 //1 quadrature point at the cell center
1247 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
1248 fileit = this->Internal->MedFiles.begin();
1249 while(fileit != this->Internal->MedFiles.end())
1251 vtkMedFile* file = fileit->second;
1254 if(def && def->GetParentFile() == file)
1255 return nloc + def->GetMedIterator() - 1;
1257 nloc += file->GetNumberOfLocalization();
1261 if(fop->GetNumberOfIntegrationPoint()==1)
1262 return nloc + 1 + fop->GetParentFieldOverEntity()->GetEntity().GeometryType;
1264 // or it is an elno field (field stored on nodes of the cells,
1265 // but with discontinuities at the vertices)
1266 return -fop->GetParentFieldOverEntity()->GetEntity().GeometryType;//ELNO
1269 int vtkMedReader::IsQuadratureFieldSelected(vtkMedField* field)
1271 return field->GetFieldType()==vtkMedField::QuadratureField
1272 &&this->GetQuadratureFieldArrayStatus(vtkMedUtilities::SimplifyName(
1273 field->GetName()).c_str());
1276 int vtkMedReader::IsElnoFieldSelected(vtkMedField* field)
1278 return field->GetFieldType()==vtkMedField::ElnoField
1279 &&this->GetElnoFieldArrayStatus(vtkMedUtilities::SimplifyName(
1280 field->GetName()).c_str());
1284 // Give the animation steps to the pipeline
1285 void vtkMedReader::AdvertiseTime(vtkInformation* info)
1287 this->ChooseRealAnimationMode();
1289 if(this->Internal->RealAnimationMode==PhysicalTime)
1291 // I advertise the union of all times available
1292 // in all selected fields and meshes
1293 set<double> timeset;
1295 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
1296 fieldfileit = this->Internal->MedFiles.begin();
1297 while(fieldfileit != this->Internal->MedFiles.end())
1299 vtkMedFile* file = fieldfileit->second;
1302 // first loop over all fields to gather their compute steps
1303 for(int fieldId=0; fieldId<file->GetNumberOfField(); fieldId++)
1305 vtkMedField* field=file->GetField(fieldId);
1307 if(!this->IsFieldSelected(field))
1310 field->GatherFieldTimes(timeset);
1313 // then loop over all meshes to gather their grid steps too.
1314 for(int meshId=0; meshId<file->GetNumberOfMesh(); meshId++)
1316 vtkMedMesh* mesh=file->GetMesh(meshId);
1318 if(!this->IsMeshSelected(mesh))
1321 mesh->GatherGridTimes(timeset);
1325 if(timeset.size() > 0)
1327 // remove MED_UNDEF_DT if there are other time step
1328 if(timeset.size() > 1)
1329 timeset.erase(MED_UNDEF_DT);
1331 vector<double> times;
1332 set<double>::iterator it = timeset.begin();
1333 while(it != timeset.end())
1335 times.push_back(*it);
1338 sort(times.begin(), times.end());
1340 info->Set(vtkStreamingDemandDrivenPipeline::TIME_STEPS(), ×[0],
1342 double timeRange[2];
1343 timeRange[0]=times[0];
1344 timeRange[1]=times[times.size()-1];
1345 info->Set(vtkStreamingDemandDrivenPipeline::TIME_RANGE(), &timeRange[0],
1350 info->Remove(vtkStreamingDemandDrivenPipeline::TIME_STEPS());
1351 info->Remove(vtkStreamingDemandDrivenPipeline::TIME_RANGE());
1354 else if(this->Internal->RealAnimationMode==Iteration)
1356 // I advertise the union of all iterations available at the given
1357 // Time for all selected fields.
1358 set<med_int> iterationsets;
1359 med_float time = MED_UNDEF_DT;
1360 if(this->TimeIndexForIterations >= 0 &&
1361 this->TimeIndexForIterations <
1362 this->AvailableTimes->GetNumberOfTuples())
1364 time = this->AvailableTimes->
1365 GetValue((vtkIdType)this->TimeIndexForIterations);
1368 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
1369 fieldfileit = this->Internal->MedFiles.begin();
1370 while(fieldfileit != this->Internal->MedFiles.end())
1372 vtkMedFile* file = fieldfileit->second;
1375 for(int fieldId=0; fieldId<file->GetNumberOfField(); fieldId++)
1377 vtkMedField* field=file->GetField(fieldId);
1378 if(!this->IsFieldSelected(field))
1381 field->GatherFieldIterations(time, iterationsets);
1383 // then loop over all meshes to gather their grid steps too.
1384 for(int meshId=0; meshId<file->GetNumberOfMesh(); meshId++)
1386 vtkMedMesh* mesh=file->GetMesh(meshId);
1388 if(!this->IsMeshSelected(mesh))
1391 mesh->GatherGridIterations(time, iterationsets);
1395 if(iterationsets.size()>0)
1397 // remove MED_NO_IT if there are other available iterations.
1398 if(iterationsets.size()>1)
1399 iterationsets.erase(MED_NO_IT);
1401 vector<double> iterations;
1402 set<med_int>::iterator it=iterationsets.begin();
1403 while(it!=iterationsets.end())
1405 iterations.push_back((double)*it);
1408 sort(iterations.begin(), iterations.end());
1409 info->Set(vtkStreamingDemandDrivenPipeline::TIME_STEPS(), &iterations[0],
1411 double timeRange[2];
1412 timeRange[0]=iterations[0];
1413 timeRange[1]=iterations[iterations.size()-1];
1414 info->Set(vtkStreamingDemandDrivenPipeline::TIME_RANGE(), &timeRange[0],
1419 info->Remove(vtkStreamingDemandDrivenPipeline::TIME_STEPS());
1420 info->Remove(vtkStreamingDemandDrivenPipeline::TIME_RANGE());
1425 info->Remove(vtkStreamingDemandDrivenPipeline::TIME_STEPS());
1426 info->Remove(vtkStreamingDemandDrivenPipeline::TIME_RANGE());
1430 vtkIdType vtkMedReader::GetFrequencyIndex(double freq)
1432 return this->AvailableTimes->LookupValue(freq);
1435 int vtkMedReader::RequestDataObject(vtkInformation* request,
1436 vtkInformationVector** vtkNotUsed(inputVector), vtkInformationVector* outputVector)
1438 vtkInformation *info = outputVector->GetInformationObject(0);
1439 if (vtkMultiBlockDataSet::SafeDownCast(
1440 info->Get(vtkDataObject::DATA_OBJECT())))
1442 // The output is already created
1447 vtkMultiBlockDataSet* output=vtkMultiBlockDataSet::New();
1448 this->GetExecutive()->SetOutputData(0, output);
1450 this->GetOutputPortInformation(0)->Set(vtkDataObject::DATA_EXTENT_TYPE(),
1451 output->GetExtentType());
1456 void vtkMedReader::ClearSelections()
1458 this->PointFields->Initialize();
1459 this->CellFields->Initialize();
1460 this->QuadratureFields->Initialize();
1461 this->ElnoFields->Initialize();
1463 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
1464 fieldfileit = this->Internal->MedFiles.begin();
1465 while(fieldfileit != this->Internal->MedFiles.end())
1467 vtkMedFile* file = fieldfileit->second;
1470 for(int index=0; index < file->GetNumberOfField(); index++)
1472 vtkMedField* field = file->GetField(index);
1473 switch(field->GetFieldType())
1475 case vtkMedField::PointField :
1476 this->PointFields->AddKey(vtkMedUtilities::SimplifyName(
1477 field->GetName()).c_str());
1479 case vtkMedField::CellField :
1480 this->CellFields->AddKey(vtkMedUtilities::SimplifyName(
1481 field->GetName()).c_str());
1483 case vtkMedField::QuadratureField :
1484 this->QuadratureFields->AddKey(vtkMedUtilities::SimplifyName(
1485 field->GetName()).c_str());
1487 case vtkMedField::ElnoField :
1488 this->ElnoFields->AddKey(vtkMedUtilities::SimplifyName(
1489 field->GetName()).c_str());
1494 this->Internal->Families->Initialize();
1495 this->Groups->Initialize();
1496 for(int meshIndex=0; meshIndex < file->GetNumberOfMesh(); meshIndex++)
1498 vtkMedMesh* mesh = file->GetMesh(meshIndex);
1499 for(int famIndex=0; famIndex<mesh->GetNumberOfPointFamily(); famIndex++)
1501 vtkMedFamily* fam=mesh->GetPointFamily(famIndex);
1503 int ng=fam->GetNumberOfGroup();
1504 for(int gindex=0; gindex<ng; gindex++)
1506 vtkMedGroup* group=fam->GetGroup(gindex);
1507 string gname=vtkMedUtilities::GroupKey(mesh->GetName(),
1508 fam->GetPointOrCell(), group->GetName());
1510 this->Groups->AddKey(gname.c_str());
1511 this->Groups->SetKeyStatus(gname.c_str(), 0);
1514 for(int famIndex=0; famIndex<mesh->GetNumberOfCellFamily(); famIndex++)
1516 vtkMedFamily* fam=mesh->GetCellFamily(famIndex);
1518 int ng=fam->GetNumberOfGroup();
1519 for(int gindex=0; gindex<ng; gindex++)
1521 vtkMedGroup* group=fam->GetGroup(gindex);
1522 string gname=vtkMedUtilities::GroupKey(mesh->GetName(),
1523 fam->GetPointOrCell(), group->GetName());
1525 this->Groups->AddKey(gname.c_str());
1526 this->Groups->SetKeyStatus(gname.c_str(), 1);
1530 this->Internal->GroupSelectionMTime.Modified();
1532 for(int meshIndex=0; meshIndex< file->GetNumberOfMesh(); meshIndex++)
1534 if(file->GetMesh(meshIndex)->GetNumberOfGridStep() == 0)
1537 vtkMedGrid* grid=file->GetMesh(meshIndex)->GetGridStep(0);
1539 for(int entityIndex=0; entityIndex<grid->GetNumberOfEntityArray();
1542 vtkMedEntityArray* array=grid->GetEntityArray(entityIndex);
1543 string name=vtkMedUtilities::EntityKey(array->GetEntity());
1544 this->Entities->AddKey(name.c_str());
1551 void vtkMedReader::SelectFamiliesFromGroups()
1553 this->Internal->Families->Initialize();
1555 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
1556 meshfileit = this->Internal->MedFiles.begin();
1557 while(meshfileit != this->Internal->MedFiles.end())
1559 vtkMedFile* file = meshfileit->second;
1562 for(int meshIndex=0; meshIndex < file->GetNumberOfMesh(); meshIndex++)
1564 vtkMedMesh* mesh = file->GetMesh(meshIndex);
1565 for(int famIndex=0; famIndex<mesh->GetNumberOfFamily(); famIndex++)
1567 vtkMedFamily* fam=mesh->GetFamily(famIndex);
1568 string name=vtkMedUtilities::FamilyKey(mesh->GetName(),
1569 fam->GetPointOrCell(), fam->GetName());
1571 this->Internal->Families->SetKeyStatus(name.c_str(), 0);
1573 for(int gindex=0; gindex<fam->GetNumberOfGroup(); gindex++)
1575 vtkMedGroup* group=fam->GetGroup(gindex);
1576 string gname=vtkMedUtilities::GroupKey(mesh->GetName(),
1577 fam->GetPointOrCell(), group->GetName());
1578 int state=this->Groups->GetKeyStatus(gname.c_str());
1582 this->SetFamilyStatus(name.c_str(), 1);
1589 this->Internal->FamilySelectionMTime.Modified();
1592 int vtkMedReader::GetNumberOfPointFieldArrays()
1594 return this->PointFields->GetNumberOfKey();
1598 vtkMedReader::GetPointFieldArrayName(int index)
1600 return this->PointFields->GetKey(index);
1603 int vtkMedReader::GetPointFieldArrayStatus(const char* name)
1605 return this->PointFields->GetKeyStatus(name);
1608 void vtkMedReader::SetPointFieldArrayStatus(const char* name, int status)
1610 if(this->PointFields->KeyExists(name)&&this->PointFields->GetKeyStatus(
1614 this->PointFields->SetKeyStatus(name, status);
1619 int vtkMedReader::GetNumberOfCellFieldArrays()
1621 return this->CellFields->GetNumberOfKey();
1625 vtkMedReader::GetCellFieldArrayName(int index)
1627 return this->CellFields->GetKey(index);
1630 int vtkMedReader::GetCellFieldArrayStatus(const char* name)
1632 return this->CellFields->GetKeyStatus(name);
1635 void vtkMedReader::SetCellFieldArrayStatus(const char* name, int status)
1637 if(this->CellFields->KeyExists(name)&&this->CellFields->GetKeyStatus(
1641 this->CellFields->SetKeyStatus(name, status);
1646 int vtkMedReader::GetNumberOfQuadratureFieldArrays()
1648 return this->QuadratureFields->GetNumberOfKey();
1651 const char* vtkMedReader::GetQuadratureFieldArrayName(int index)
1653 return this->QuadratureFields->GetKey(index);
1656 int vtkMedReader::GetQuadratureFieldArrayStatus(const char* name)
1658 return this->QuadratureFields->GetKeyStatus(name);
1661 void vtkMedReader::SetQuadratureFieldArrayStatus(const char* name, int status)
1663 if(this->QuadratureFields->KeyExists(name)
1664 &&this->QuadratureFields->GetKeyStatus(name)==status)
1667 this->QuadratureFields->SetKeyStatus(name, status);
1672 int vtkMedReader::GetNumberOfElnoFieldArrays()
1674 return this->ElnoFields->GetNumberOfKey();
1677 const char* vtkMedReader::GetElnoFieldArrayName(int index)
1679 return this->ElnoFields->GetKey(index);
1682 int vtkMedReader::GetElnoFieldArrayStatus(const char* name)
1684 return this->ElnoFields->GetKeyStatus(name);
1687 void vtkMedReader::SetElnoFieldArrayStatus(const char* name, int status)
1689 if(this->ElnoFields->KeyExists(name)
1690 &&this->ElnoFields->GetKeyStatus(name)==status)
1693 this->ElnoFields->SetKeyStatus(name, status);
1698 void vtkMedReader::SetEntityStatus(const char* name, int status)
1700 if(this->Entities->KeyExists(name)&&this->Entities->GetKeyStatus(name)
1704 this->Entities->SetKeyStatus(name, status);
1709 void vtkMedReader::SetFamilyStatus(const char* name, int status)
1711 if(this->Internal->Families->KeyExists(name)
1712 &&this->Internal->Families->GetKeyStatus(name)==status)
1715 this->Internal->Families->SetKeyStatus(name, status);
1718 void vtkMedReader::SetGroupStatus(const char* name, int status)
1721 if(this->Groups->KeyExists(name)&&this->Groups->GetKeyStatus(name)
1725 this->Groups->SetKeyStatus(name, status);
1729 this->Internal->GroupSelectionMTime.Modified();
1732 int vtkMedReader::GetGroupStatus(const char* key)
1734 return this->Groups->GetKeyStatus(key);
1737 void vtkMedReader::AddQuadratureSchemeDefinition(vtkInformation* info,
1738 vtkMedLocalization* loc)
1740 if(info==NULL||loc==NULL)
1743 vtkInformationQuadratureSchemeDefinitionVectorKey *key=
1744 vtkQuadratureSchemeDefinition::DICTIONARY();
1746 vtkQuadratureSchemeDefinition* def=vtkQuadratureSchemeDefinition::New();
1747 int cellType=vtkMedUtilities::GetVTKCellType(loc->GetGeometryType());
1748 def->Initialize(cellType, vtkMedUtilities::GetNumberOfPoint(
1749 loc->GetGeometryType()), loc->GetNumberOfQuadraturePoint(),
1750 (double*)loc->GetShapeFunction()->GetVoidPointer(0),
1751 (double*)loc->GetWeights()->GetVoidPointer(0));
1752 key->Set(info, def, cellType);
1757 void vtkMedReader::LoadConnectivity(vtkMedEntityArray* array)
1759 vtkMedGrid* grid = array->GetParentGrid();
1760 array->LoadConnectivity();
1761 if (array->GetConnectivity()==MED_NODAL||vtkMedUtilities::GetDimension(
1762 array->GetEntity().GeometryType)<2
1763 || grid->GetParentMesh()->GetMeshType() == MED_STRUCTURED_MESH)
1766 vtkMedEntity subentity;
1767 subentity.EntityType = vtkMedUtilities::GetSubType(array->GetEntity().EntityType);
1769 vtkMedUnstructuredGrid* ugrid = vtkMedUnstructuredGrid::SafeDownCast(grid);
1773 for(int index=0; index<vtkMedUtilities::GetNumberOfSubEntity(
1774 array->GetEntity().GeometryType); index++)
1776 subentity.GeometryType = vtkMedUtilities::GetSubGeometry(
1777 array->GetEntity().GeometryType, index);
1778 vtkMedEntityArray* subarray=ugrid->GetEntityArray(subentity);
1782 vtkErrorMacro("DESC connectivity used, but sub types do not exist in file.");
1785 subarray->LoadConnectivity();
1789 vtkDataSet* vtkMedReader::CreateUnstructuredGridForPointSupport(
1790 vtkMedFamilyOnEntityOnProfile* foep)
1792 vtkUnstructuredGrid* vtkgrid = vtkUnstructuredGrid::New();
1793 foep->ComputeIntersection(NULL);
1794 vtkMedFamilyOnEntity* foe = foep->GetFamilyOnEntity();
1795 vtkMedGrid* medgrid=foe->GetParentGrid();
1796 vtkMedUnstructuredGrid* medugrid=vtkMedUnstructuredGrid::SafeDownCast(
1798 vtkMedCurvilinearGrid* medcgrid=vtkMedCurvilinearGrid::SafeDownCast(
1801 medgrid->LoadCoordinates();
1803 vtkIdType npts=medgrid->GetNumberOfPoints();
1805 bool shallowCopy= (medugrid != NULL || medcgrid!=NULL);
1806 if(medgrid->GetParentMesh()->GetSpaceDimension()!=3)
1812 shallowCopy = foep->CanShallowCopyPointField(NULL);
1815 vtkDataArray* coords = NULL;
1817 if(medugrid != NULL)
1818 coords = medugrid->GetCoordinates();
1819 if(medcgrid != NULL)
1820 coords = medcgrid->GetCoordinates();
1823 vtkIdType numberOfPoints;
1824 vtkPoints* points=vtkPoints::New(coords->GetDataType());
1825 vtkgrid->SetPoints(points);
1828 vtkIdTypeArray* pointGlobalIds=vtkIdTypeArray::New();
1829 pointGlobalIds->SetName("MED_POINT_ID");
1830 pointGlobalIds->SetNumberOfComponents(1);
1831 vtkgrid->GetPointData()->SetGlobalIds(pointGlobalIds);
1832 pointGlobalIds->Delete();
1836 vtkgrid->GetPoints()->SetData(coords);
1837 numberOfPoints=npts;
1839 pointGlobalIds->SetNumberOfTuples(numberOfPoints);
1840 vtkIdType* ptr=pointGlobalIds->GetPointer(0);
1841 for(int pid=0; pid<numberOfPoints; pid++)
1846 vtkIdType currentIndex=0;
1848 for(vtkIdType index=0; index<medgrid->GetNumberOfPoints(); index++)
1850 if (!foep->KeepPoint(index))
1855 double coord[3]={0.0, 0.0, 0.0};
1856 double * tuple=medgrid->GetCoordTuple(index);
1857 for(int dim=0; dim<medgrid->GetParentMesh()->GetSpaceDimension()&&dim<3; dim++)
1859 coord[dim]=tuple[dim];
1861 vtkgrid->GetPoints()->InsertPoint(currentIndex, coord);
1862 pointGlobalIds->InsertNextValue(index+1);
1865 vtkgrid->GetPoints()->Squeeze();
1866 pointGlobalIds->Squeeze();
1867 numberOfPoints=currentIndex;
1870 // now create the VTK_VERTEX cells
1871 for(vtkIdType id=0; id<numberOfPoints; id++)
1873 vtkgrid->InsertNextCell(VTK_VERTEX, 1, &id);
1877 // in this particular case, the global ids of the cells is the same as the global ids of the points.
1878 vtkgrid->GetCellData()->SetGlobalIds(vtkgrid->GetPointData()->GetGlobalIds());
1883 vtkMedGrid* vtkMedReader::FindGridStep(vtkMedMesh* mesh)
1885 if(this->Internal->RealAnimationMode == vtkMedReader::PhysicalTime)
1887 vtkMedComputeStep cs;
1888 cs.TimeOrFrequency = this->Internal->UpdateTimeStep;
1889 return mesh->FindGridStep(cs, vtkMedReader::PhysicalTime);
1891 else if(this->Internal->RealAnimationMode == vtkMedReader::Modes)
1893 vtkMedComputeStep cs;
1894 cs.IterationIt = MED_NO_IT;
1895 cs.TimeIt = MED_NO_DT;
1896 cs.TimeOrFrequency = MED_NO_DT;
1897 return mesh->FindGridStep(cs, vtkMedReader::Modes);
1901 vtkMedComputeStep cs;
1902 // the time is set by choosing its index in the global
1903 // array giving the available times : this->TimeIndexForIterations
1904 cs.TimeOrFrequency = (med_int)this->AvailableTimes->GetValue(
1905 (vtkIdType)this->TimeIndexForIterations);
1906 // the iteration is asked by the pipeline
1907 cs.IterationIt = (med_int)this->Internal->UpdateTimeStep;
1908 return mesh->FindGridStep(cs, vtkMedReader::Iteration);
1913 void vtkMedReader::CreateMedSupports()
1915 this->Internal->UsedSupports.clear();
1917 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
1918 meshfileit = this->Internal->MedFiles.begin();
1919 while(meshfileit != this->Internal->MedFiles.end())
1921 vtkMedFile* file = meshfileit->second;
1924 for(int meshIndex=0; meshIndex<file->GetNumberOfMesh(); meshIndex++)
1926 vtkMedMesh* mesh = file->GetMesh(meshIndex);
1927 vtkMedGrid* grid = this->FindGridStep(mesh);
1931 for(int entityIndex=0; entityIndex<grid->GetNumberOfEntityArray();
1934 vtkMedEntityArray* array=grid->GetEntityArray(entityIndex);
1935 if(this->GetEntityStatus(array->GetEntity())==0)
1940 file->GetMedDriver()->LoadFamilyIds(array);
1941 for(int foeIndex=0; foeIndex<array->GetNumberOfFamilyOnEntity();
1944 vtkMedFamilyOnEntity* foe=array->GetFamilyOnEntity(foeIndex);
1945 vtkMedFamily* family=foe->GetFamily();
1946 if(this->GetFamilyStatus(mesh, family)==0)
1949 // now, I look over all non-point fields to see which profiles
1950 // have to be used on points.
1951 bool selectedSupport = false;
1953 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
1954 fieldfileit = this->Internal->MedFiles.begin();
1955 while(fieldfileit != this->Internal->MedFiles.end())
1957 vtkMedFile* fieldfile = fieldfileit->second;
1960 for(int fieldId=0; fieldId<fieldfile->GetNumberOfField(); fieldId++)
1962 vtkMedField* field=fieldfile->GetField(fieldId);
1964 if (!this->IsFieldSelected(field))
1967 vtkMedListOfFieldSteps steps;
1969 this->GatherFieldSteps(field, steps);
1971 vtkMedListOfFieldSteps::iterator it=steps.begin();
1972 while(it!=steps.end())
1974 vtkMedFieldStep *step = *it;
1975 step->LoadInformation();
1978 for(int eid=0; eid<step->GetNumberOfFieldOverEntity(); eid++)
1980 vtkMedFieldOverEntity* fieldOverEntity =
1981 step->GetFieldOverEntity(eid);
1983 for(int pid = 0; pid < fieldOverEntity->GetNumberOfFieldOnProfile(); pid++)
1985 vtkMedFieldOnProfile* fop =
1986 fieldOverEntity->GetFieldOnProfile(pid);
1987 vtkMedProfile* prof = fop->GetProfile();
1988 vtkMedFamilyOnEntityOnProfile* foep =
1989 foe->GetFamilyOnEntityOnProfile(prof);
1992 this->Internal->UsedSupports.insert(foep);
1993 selectedSupport = true;
2000 // If no field use this family on entity, I nevertheless create the
2001 // support, with an empty profile.
2002 if(!selectedSupport)
2004 vtkMedFamilyOnEntityOnProfile* foep =
2005 foe->GetFamilyOnEntityOnProfile((vtkMedProfile*)NULL);
2008 foep = vtkMedFamilyOnEntityOnProfile::New();
2009 foep->SetFamilyOnEntity(foe);
2010 foep->SetProfile(NULL);
2011 foe->AddFamilyOnEntityOnProfile(foep);
2014 this->Internal->UsedSupports.insert(foep);
2022 bool vtkMedReader::BuildVTKSupport(
2023 vtkMedFamilyOnEntityOnProfile* foep,
2027 vtkMedFamilyOnEntity* foe = foep->GetFamilyOnEntity();
2030 vtkMultiProcessController* controller = vtkMultiProcessController::GetGlobalController();
2031 if (controller != NULL)
2033 numProc = controller->GetNumberOfProcesses();
2036 if ((foep->GetValid() == 0) && numProc == 1)
2041 vtkMedGrid* grid=foe->GetParentGrid();
2043 vtkMedEntityArray* array=foe->GetEntityArray();
2044 vtkMedMesh* mesh=grid->GetParentMesh();
2045 vtkSmartPointer<vtkStringArray> path = vtkSmartPointer<vtkStringArray>::New();
2046 string meshName=vtkMedUtilities::SimplifyName(mesh->GetName());
2047 path->InsertNextValue(meshName);
2048 std::string finalName;
2050 if(foe->GetPointOrCell()==vtkMedUtilities::OnPoint)
2052 path->InsertNextValue(vtkMedUtilities::OnPointName);
2053 finalName=vtkMedUtilities::SimplifyName(foe->GetFamily()->GetName());
2057 path->InsertNextValue(vtkMedUtilities::OnCellName);
2058 path->InsertNextValue(vtkMedUtilities::SimplifyName(foe->GetFamily()->GetName()));
2059 finalName=vtkMedUtilities::EntityKey(array->GetEntity());
2062 if(foep->GetProfile() != NULL)
2064 path->InsertNextValue(finalName);
2065 finalName = foep->GetProfile()->GetName();
2068 ostringstream progressBarTxt;
2069 for(int depth=0; depth<path->GetNumberOfValues(); depth++)
2071 progressBarTxt<<path->GetValue(depth)<<" ";
2073 progressBarTxt<<finalName;
2074 SetProgressText(progressBarTxt.str().c_str());
2076 vtkDataSet* cachedDataSet = NULL;
2077 if(this->Internal->DataSetCache.find(foep)!=this->Internal->DataSetCache.end())
2079 cachedDataSet = this->Internal->DataSetCache[foep];
2083 vtkDataSet* dataset = NULL;
2086 if(foe->GetPointOrCell()==vtkMedUtilities::OnPoint)
2088 dataset = this->CreateUnstructuredGridForPointSupport(foep);
2092 dataset = foep->GetFamilyOnEntity()->GetParentGrid()->
2093 CreateVTKDataSet(foep);
2102 this->Internal->DataSetCache[foep]=dataset;
2103 cachedDataSet = dataset;
2108 vtkMultiBlockDataSet* root=vtkMedUtilities::GetParent(this->GetOutput(), path);
2109 int nb=root->GetNumberOfBlocks();
2111 if(cachedDataSet != NULL)
2113 vtkDataSet* realDataSet=cachedDataSet->NewInstance();
2114 root->SetBlock(nb, realDataSet);
2115 realDataSet->Delete();
2117 root->GetMetaData(nb)->Set(vtkCompositeDataSet::NAME(), finalName.c_str());
2118 realDataSet->ShallowCopy(cachedDataSet);
2120 this->Internal->DataSetCache[foep]=cachedDataSet;
2121 this->Internal->CurrentDataSet[foep]=realDataSet;
2123 path->InsertNextValue(finalName);
2124 path->SetName("BLOCK_NAME");
2125 realDataSet->GetFieldData()->AddArray(path);
2126 realDataSet->GetInformation()->Remove(vtkMedUtilities::BLOCK_NAME());
2127 for(int depth=0; depth<path->GetNumberOfValues(); depth++)
2129 realDataSet->GetInformation()->Set(vtkMedUtilities::BLOCK_NAME(),
2130 path->GetValue(depth), depth);
2135 void vtkMedReader::MapFieldOnSupport(vtkMedFieldOnProfile* fop,
2136 vtkMedFamilyOnEntityOnProfile* foep,
2139 bool cached = false;
2141 if(this->Internal->FieldCache.find(foep)
2142 !=this->Internal->FieldCache.end())
2144 map<vtkMedFieldOnProfile*, VTKField>& fieldCache =
2145 this->Internal->FieldCache[foep];
2146 if(fieldCache.find(fop)!=fieldCache.end())
2154 this->CreateVTKFieldOnSupport(fop, foep, doCreateField);
2156 this->SetVTKFieldOnSupport(fop, foep);
2159 void vtkMedReader::MapFieldsOnSupport(vtkMedFamilyOnEntityOnProfile* foep,
2162 // now loop over all fields to map it to the created grids
2163 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
2164 fieldfileit = this->Internal->MedFiles.begin();
2165 while(fieldfileit != this->Internal->MedFiles.end())
2167 vtkMedFile* file = fieldfileit->second;
2170 for(int fieldId=0; fieldId<file->GetNumberOfField(); fieldId++)
2172 vtkMedField* field=file->GetField(fieldId);
2174 if(strcmp(foep->GetFamilyOnEntity()->
2175 GetParentGrid()->GetParentMesh()->GetName(),
2176 field->GetMeshName()) != 0)
2179 if(strcmp(foep->GetFamilyOnEntity()->
2180 GetParentGrid()->GetParentMesh()->GetName(),
2181 field->GetMeshName()) != 0)
2184 if (!this->IsFieldSelected(field))
2187 vtkMedListOfFieldSteps steps;
2189 this->GatherFieldSteps(field, steps);
2191 vtkMedListOfFieldSteps::iterator it=steps.begin();
2192 while(it!=steps.end())
2194 vtkMedFieldStep *step = *it;
2195 step->LoadInformation();
2198 for(int eid=0; eid<step->GetNumberOfFieldOverEntity(); eid++)
2200 vtkMedFieldOverEntity* fieldOverEntity = step->GetFieldOverEntity(eid);
2201 for(int pid = 0; pid < fieldOverEntity->GetNumberOfFieldOnProfile(); pid++)
2203 vtkMedFieldOnProfile* fop = fieldOverEntity->GetFieldOnProfile(pid);
2204 if(foep->CanMapField(fop))
2206 this->MapFieldOnSupport(fop, foep, doCreateField);
2215 void vtkMedReader::GatherFieldSteps(vtkMedField* field,
2216 vtkMedListOfFieldSteps& steps)
2218 if(this->Internal->RealAnimationMode == vtkMedReader::PhysicalTime)
2220 vtkMedComputeStep cs;
2221 cs.TimeOrFrequency = this->Internal->UpdateTimeStep;
2222 vtkMedFieldStep* fs =
2223 field->FindFieldStep(cs, vtkMedReader::PhysicalTime);
2225 steps.push_back(fs);
2227 else if(this->Internal->RealAnimationMode == vtkMedReader::Iteration)
2229 vtkMedComputeStep cs;
2230 cs.IterationIt = (med_int)this->Internal->UpdateTimeStep;
2231 cs.TimeOrFrequency = (med_int)this->AvailableTimes->GetValue(
2232 (vtkIdType)this->TimeIndexForIterations);
2233 vtkMedFieldStep* fs =
2234 field->FindFieldStep(cs, vtkMedReader::Iteration);
2236 steps.push_back(fs);
2240 for(int modeid = 0; modeid < this->Frequencies->GetNumberOfKey(); modeid++)
2242 if(this->Frequencies->GetKeyStatus(
2243 this->Frequencies->GetKey(modeid)) == 0)
2248 vtkMedComputeStep cs;
2249 cs.TimeOrFrequency = this->AvailableTimes->GetValue(modeid);
2250 vtkMedFieldStep* fs =
2251 field->FindFieldStep(cs, vtkMedReader::PhysicalTime);
2253 steps.push_back(fs);
2258 void vtkMedReader::SetVTKFieldOnSupport(vtkMedFieldOnProfile* fop,
2259 vtkMedFamilyOnEntityOnProfile* foep)
2261 vtkMedFamilyOnEntity* foe = foep->GetFamilyOnEntity();
2262 vtkMedFieldOverEntity* fieldOverEntity = fop->GetParentFieldOverEntity();
2263 vtkMedFieldStep* step = fieldOverEntity->GetParentStep();
2264 vtkMedField* field = step->GetParentField();
2266 const vtkMedComputeStep& cs = step->GetComputeStep();
2268 vtkDataSet* ds=this->Internal->CurrentDataSet[foep];
2271 // ds == NULL could arrive is some cases when working in parallel
2272 vtkWarningMacro( "--- vtkMedReader::SetVTKFieldOnSupport: ds is NULL !!");
2276 VTKField& vtkfield=this->Internal->FieldCache[foep][fop];
2278 std::string name=vtkMedUtilities::SimplifyName(field->GetName());
2279 std::string vectorname = name+"_Vector";
2280 if(this->AnimationMode==Modes)
2282 double freq = cs.TimeOrFrequency;
2283 int index = this->GetFrequencyIndex(freq);
2284 name += string(" ") + vtkMedUtilities::GetModeKey(index, freq,
2285 this->AvailableTimes->GetNumberOfTuples()-1);
2286 vectorname += string(" ") + vtkMedUtilities::GetModeKey(index, freq,
2287 this->AvailableTimes->GetNumberOfTuples()-1);
2290 vtkfield.DataArray->SetName(name.c_str());
2291 vtkfield.DataArray->Squeeze();
2292 if(vtkfield.Vectors != NULL)
2294 vtkfield.Vectors->SetName(vectorname.c_str());
2295 vtkfield.Vectors->Squeeze();
2297 if(vtkfield.QuadratureIndexArray!=NULL)
2299 vtkfield.QuadratureIndexArray->Squeeze();
2302 if(foe->GetPointOrCell()==vtkMedUtilities::OnCell)
2304 if(field->GetFieldType()==vtkMedField::PointField)
2306 if(vtkfield.DataArray->GetNumberOfTuples()!=ds->GetNumberOfPoints())
2308 vtkDebugMacro("the data array " << vtkfield.DataArray->GetName()
2309 << " do not have the good number of tuples");
2312 ds->GetPointData()->AddArray(vtkfield.DataArray);
2313 if(vtkfield.Vectors != NULL && this->GenerateVectors)
2315 ds->GetPointData()->AddArray(vtkfield.Vectors);
2316 ds->GetPointData()->SetActiveVectors(vtkfield.Vectors->GetName());
2318 switch (vtkfield.DataArray->GetNumberOfComponents())
2321 ds->GetPointData()->SetActiveScalars(vtkfield.DataArray->GetName());
2324 ds->GetPointData()->SetActiveVectors(vtkfield.DataArray->GetName());
2327 // if the data set is only composed of VTK_VERTEX cells,
2328 // and no field called with the same name exist on cells,
2329 // map this field to cells too
2330 if(foe->GetVertexOnly()==1 && ds->GetCellData()->GetArray(
2331 vtkfield.DataArray->GetName())==NULL)
2333 ds->GetCellData()->AddArray(vtkfield.DataArray);
2334 if(vtkfield.Vectors != NULL && this->GenerateVectors)
2336 ds->GetCellData()->AddArray(vtkfield.Vectors);
2337 ds->GetCellData()->SetActiveVectors(vtkfield.Vectors->GetName());
2339 switch (vtkfield.DataArray->GetNumberOfComponents())
2342 ds->GetCellData()->SetActiveScalars(vtkfield.DataArray->GetName());
2345 ds->GetCellData()->SetActiveVectors(vtkfield.DataArray->GetName());
2350 if(field->GetFieldType()==vtkMedField::CellField)
2352 if((this->Internal->NumberOfPieces == 1) && vtkfield.DataArray->GetNumberOfTuples()!=ds->GetNumberOfCells() )
2354 vtkDebugMacro("the data array " << vtkfield.DataArray->GetName()
2355 << " do not have the good number of tuples"
2356 << " ncell=" << ds->GetNumberOfCells()
2357 << " ntuple=" << vtkfield.DataArray->GetNumberOfTuples());
2360 // In case we are in parallel and our process does not contain the data
2361 if(ds->GetNumberOfCells()==0)
2365 ds->GetCellData()->AddArray(vtkfield.DataArray);
2367 if(vtkfield.Vectors != NULL && this->GenerateVectors)
2369 ds->GetCellData()->AddArray(vtkfield.Vectors);
2370 ds->GetCellData()->SetActiveVectors(vtkfield.Vectors->GetName());
2372 switch (vtkfield.DataArray->GetNumberOfComponents())
2375 ds->GetCellData()->SetActiveScalars(vtkfield.DataArray->GetName());
2378 ds->GetCellData()->SetActiveVectors(vtkfield.DataArray->GetName());
2381 // if the data set is only composed of VTK_VERTEX cells,
2382 // and no field called with the same name exist on points,
2383 // map this field to points too
2384 if(foe->GetVertexOnly()==1 && ds->GetPointData()->GetArray(
2385 vtkfield.DataArray->GetName())==NULL)
2387 ds->GetPointData()->AddArray(vtkfield.DataArray);
2388 if(vtkfield.Vectors != NULL && this->GenerateVectors)
2390 ds->GetPointData()->AddArray(vtkfield.Vectors);
2391 ds->GetPointData()->SetActiveVectors(vtkfield.Vectors->GetName());
2393 switch (vtkfield.DataArray->GetNumberOfComponents())
2396 ds->GetPointData()->SetActiveScalars(vtkfield.DataArray->GetName());
2399 ds->GetPointData()->SetActiveVectors(vtkfield.DataArray->GetName());
2404 if(field->GetFieldType()==vtkMedField::QuadratureField ||
2405 field->GetFieldType()==vtkMedField::ElnoField )
2407 vtkIdType ncells=ds->GetNumberOfCells();
2408 vtkIdType nid=vtkfield.QuadratureIndexArray->GetNumberOfTuples();
2409 vtkIdType nda=vtkfield.DataArray->GetNumberOfTuples();
2410 if((nid!=ncells) && (this->Internal->NumberOfPieces == 1))
2413 "There should be as many quadrature index values as there are cells");
2420 vtkfield.DataArray->SetNumberOfTuples( 0 );
2421 vtkfield.DataArray->Squeeze();
2423 if (nid>ncells) // PROBABLY NOT NECESSARY
2425 vtkfield.QuadratureIndexArray->SetNumberOfTuples(ncells);
2426 int nquad=fop->GetNumberOfIntegrationPoint();
2427 vtkfield.DataArray->SetNumberOfTuples( nquad * ds->GetNumberOfCells() );
2428 vtkfield.DataArray->Squeeze();
2430 ds->GetFieldData()->AddArray(vtkfield.DataArray);
2431 ds->GetCellData()->AddArray(vtkfield.QuadratureIndexArray);
2433 nid=vtkfield.QuadratureIndexArray->GetNumberOfTuples();
2434 nda=vtkfield.DataArray->GetNumberOfTuples();
2436 if(vtkfield.Vectors != NULL && this->GenerateVectors)
2438 ds->GetFieldData()->AddArray(vtkfield.Vectors);
2445 if(vtkfield.DataArray->GetNumberOfTuples()!=ds->GetNumberOfPoints())
2447 vtkDebugMacro("the data array " << vtkfield.DataArray->GetName() << " do not have the good number of tuples");
2450 ds->GetPointData()->AddArray(vtkfield.DataArray);
2451 if(vtkfield.Vectors != NULL && this->GenerateVectors)
2453 ds->GetPointData()->AddArray(vtkfield.Vectors);
2454 ds->GetPointData()->SetActiveVectors(vtkfield.Vectors->GetName());
2456 switch (vtkfield.DataArray->GetNumberOfComponents())
2459 ds->GetPointData()->SetActiveScalars(vtkfield.DataArray->GetName());
2462 ds->GetPointData()->SetActiveVectors(vtkfield.DataArray->GetName());
2465 // all the VTK_VERTEX created cells have the same order than the points,
2466 // I can safely map the point array to the cells in this particular case.
2467 ds->GetCellData()->AddArray(vtkfield.DataArray);
2468 if(vtkfield.Vectors != NULL && this->GenerateVectors)
2470 ds->GetCellData()->AddArray(vtkfield.Vectors);
2471 ds->GetCellData()->SetActiveVectors(vtkfield.Vectors->GetName());
2473 switch (vtkfield.DataArray->GetNumberOfComponents())
2476 ds->GetCellData()->SetActiveScalars(vtkfield.DataArray->GetName());
2479 ds->GetCellData()->SetActiveVectors(vtkfield.DataArray->GetName());
2485 void vtkMedReader::InitializeQuadratureOffsets(vtkMedFieldOnProfile* fop,
2486 vtkMedFamilyOnEntityOnProfile* foep)
2488 vtkMedFamilyOnEntity* foe = foep->GetFamilyOnEntity();
2490 vtkMedFieldOverEntity* fieldOverEntity = fop->GetParentFieldOverEntity();
2491 vtkMedFieldStep* step = fieldOverEntity->GetParentStep();
2492 vtkMedField *field= step->GetParentField();
2494 // then I compute the quadrature key if needed, and try to find the quadrature offsets.
2495 if(this->Internal->QuadratureOffsetCache.find(foep)
2496 ==this->Internal->QuadratureOffsetCache.end())
2497 this->Internal->QuadratureOffsetCache[foep]=map<LocalizationKey,
2498 vtkSmartPointer<vtkIdTypeArray> > ();
2500 map<LocalizationKey, vtkSmartPointer<vtkIdTypeArray> >& quadOffsets=
2501 this->Internal->QuadratureOffsetCache[foep];
2503 LocalizationKey quadkey=this->GetLocalizationKey(fop);
2505 if(quadOffsets.find(quadkey)!=quadOffsets.end())
2506 {// the quadrature offset array has already been created
2510 vtkIdTypeArray* qoffsets=vtkIdTypeArray::New();
2511 quadOffsets[quadkey]=qoffsets;
2515 if(field->GetFieldType() == vtkMedField::ElnoField)
2517 qoffsets->GetInformation()->Set(vtkMedUtilities::ELNO(), 1);
2520 else if(field->GetFieldType() == vtkMedField::QuadratureField)
2522 qoffsets->GetInformation()->Set(vtkMedUtilities::ELGA(), 1);
2527 sstr<<"QuadraturePointOffset";
2530 qoffsets->SetName(sstr.str().c_str());
2532 vtkSmartPointer<vtkMedLocalization> loc=
2533 this->GetLocalization(fop->GetLocalizationName());
2537 if(fop->GetNumberOfIntegrationPoint()==1)
2538 {// cell-centered fields can be seen as quadrature fields with 1
2539 // quadrature point at the center
2540 vtkMedLocalization* center=vtkMedLocalization::New();
2543 center->BuildCenter(fieldOverEntity->GetEntity().GeometryType);
2545 else if(loc == NULL && field->GetFieldType() == vtkMedField::ElnoField)
2546 {// ELNO fields have no vtkMedLocalization,
2547 // I need to create a dummy one
2548 vtkMedLocalization* elnodef=vtkMedLocalization::New();
2551 elnodef->BuildELNO(fieldOverEntity->GetEntity().GeometryType);
2555 vtkErrorMacro("Cannot find localization of quadrature field "
2556 << field->GetName());
2559 this->AddQuadratureSchemeDefinition(qoffsets->GetInformation(), loc);
2562 void vtkMedReader::CreateVTKFieldOnSupport(vtkMedFieldOnProfile* fop,
2563 vtkMedFamilyOnEntityOnProfile* foep, int doCreateField)
2565 vtkMedFamilyOnEntity* foe = foep->GetFamilyOnEntity();
2566 vtkMedFieldOverEntity* fieldOverEntity = fop->GetParentFieldOverEntity();
2567 vtkMedFieldStep* step = fieldOverEntity->GetParentStep();
2568 vtkMedField* field = step->GetParentField();
2570 if(vtkMedUnstructuredGrid::SafeDownCast(
2571 foep->GetFamilyOnEntity()->GetParentGrid()) != NULL)
2573 fop->Load(MED_COMPACT_STMODE);
2577 fop->Load(MED_GLOBAL_STMODE);
2580 vtkMedIntArray* profids=NULL;
2582 vtkMedProfile* profile=fop->GetProfile();
2587 profids=profile->GetIds();
2590 VTKField& vtkfield=this->Internal->FieldCache[foep][fop];
2592 bool shallowok = true;
2594 // for structured grid, the values are loaded globally, and cells which are
2595 // not on the profile or not on the family are blanked out.
2596 // shallow copy is therefore always possible
2597 if(vtkMedUnstructuredGrid::SafeDownCast(
2598 foep->GetFamilyOnEntity()->GetParentGrid()) != NULL)
2600 shallowok = foep->CanShallowCopy(fop);
2605 vtkfield.DataArray = fop->GetData();
2609 vtkDataArray* data=vtkMedUtilities::NewArray(field->GetDataType());
2610 vtkfield.DataArray=data;
2612 vtkfield.DataArray->SetNumberOfComponents(field->GetNumberOfComponent());
2615 for(int comp=0; comp<field->GetNumberOfComponent(); comp++)
2617 vtkfield.DataArray->SetComponentName(comp, vtkMedUtilities::SimplifyName(
2618 field->GetComponentName()->GetValue(comp)).c_str());
2621 bool createOffsets=false;
2622 if(field->GetFieldType()==vtkMedField::QuadratureField ||
2623 field->GetFieldType()==vtkMedField::ElnoField )
2625 this->InitializeQuadratureOffsets(fop, foep);
2627 LocalizationKey quadKey = this->GetLocalizationKey(fop);
2628 vtkfield.QuadratureIndexArray
2629 =this->Internal->QuadratureOffsetCache[foep][quadKey];
2630 vtkDataSet* ds = this->Internal->CurrentDataSet[foep];
2632 vtkfield.DataArray->GetInformation()->Set(
2633 vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME(),
2634 vtkfield.QuadratureIndexArray->GetName());
2635 vtkfield.QuadratureIndexArray->GetInformation()->Set(
2636 vtkAbstractArray::GUI_HIDE(), 1);
2638 if(vtkfield.QuadratureIndexArray->GetNumberOfTuples()
2639 !=ds->GetNumberOfCells())
2641 vtkfield.QuadratureIndexArray->SetNumberOfTuples(0);
2651 // build the quadrature offset array if needed
2655 int nquad=fop->GetNumberOfIntegrationPoint();
2656 noffsets=fop->GetData()->GetNumberOfTuples()/nquad;
2658 vtkDataSet* ds=this->Internal->CurrentDataSet[foep];
2659 if(noffsets!=ds->GetNumberOfCells())
2662 "number of quadrature offsets (" << noffsets << ") must "
2663 << "match the number of cells (" << ds->GetNumberOfCells() << ")!");
2666 vtkfield.QuadratureIndexArray->SetNumberOfTuples(noffsets);
2667 for(vtkIdType id=0; id<noffsets; id++)
2669 vtkfield.QuadratureIndexArray->SetValue(id, nquad*id);
2674 else if(foe->GetPointOrCell() == vtkMedUtilities::OnCell
2675 && field->GetFieldType() != vtkMedField::PointField)
2677 // Cell-centered field on cell support
2679 if (field->GetFieldType()==vtkMedField::QuadratureField ||
2680 field->GetFieldType()==vtkMedField::ElnoField)
2682 nquad=fop->GetNumberOfIntegrationPoint();
2684 vtkMedIntArray* profids=NULL;
2687 profids=profile->GetIds();
2689 vtkIdType maxIndex=fop->GetData()->GetNumberOfTuples();
2691 vtkIdType quadIndex = 0;
2693 for (vtkIdType id = 0; id<maxIndex; id++)
2695 vtkIdType realIndex = (profids!=NULL ? profids->GetValue(id)-1 : id);
2696 if (!foep->KeepCell(realIndex))
2699 if (field->GetFieldType()==vtkMedField::QuadratureField ||
2700 field->GetFieldType()==vtkMedField::ElnoField)
2702 for (int q = 0; q<nquad; q++)
2704 vtkfield.DataArray->InsertNextTuple(nquad*realIndex+q,
2709 vtkfield.QuadratureIndexArray->InsertNextValue(quadIndex);
2715 vtkfield.DataArray->InsertNextTuple(id, fop->GetData());
2718 vtkfield.DataArray->Squeeze();
2719 vtkDataSet* ds = this->Internal->CurrentDataSet[foep];
2721 else if(foe->GetPointOrCell() == vtkMedUtilities::OnCell
2722 && field->GetFieldType() == vtkMedField::PointField)
2723 {// point field on cell support
2724 vtkMedIntArray* profids=NULL;
2726 vtkMedProfile* profile=fop->GetProfile();
2732 profids=profile->GetIds();
2735 vtkIdType maxId=fop->GetData()->GetNumberOfTuples();
2737 for(vtkIdType id=0; id<maxId; id++)
2739 // if I have a profile, then I should insert the value at the position given in the profile.
2740 vtkIdType destIndex;
2743 destIndex=profids->GetValue(id)-1; // -1 because med indices start at 1
2750 if(!foep->KeepPoint(destIndex))
2752 // if I use the med2VTKIndex, then the index to insert
2753 // this value is given by the map.
2754 destIndex = foep->GetVTKPointIndex(destIndex);
2755 vtkfield.DataArray->InsertTuple(destIndex, id, fop->GetData());
2757 vtkfield.DataArray->Squeeze();
2758 }// point field on cell support
2759 else if(foe->GetPointOrCell() == vtkMedUtilities::OnPoint &&
2760 field->GetFieldType() == vtkMedField::PointField)
2763 vtkIdType maxId = fop->GetData()->GetNumberOfTuples();
2765 for(vtkIdType id=0; id<maxId; id++)
2767 vtkIdType realIndex=id;
2770 realIndex=profids->GetValue(id)-1; // -1 because med indices start at 1
2773 if(!foep->KeepPoint(realIndex))
2776 vtkfield.DataArray->InsertNextTuple(fop->GetData()->GetTuple(id));
2778 vtkfield.DataArray->Squeeze();
2779 }// support on point
2781 // now generate the vector field if asked for
2782 if(this->GenerateVectors)
2784 int ncomp = vtkfield.DataArray->GetNumberOfComponents();
2785 if(ncomp > 1 && ncomp != 3)
2787 vtkDataArray* vectors = vtkfield.DataArray->NewInstance();
2788 vectors->SetNumberOfComponents(3);
2789 vectors->SetNumberOfTuples(vtkfield.DataArray->GetNumberOfTuples());
2790 vtkfield.Vectors = vectors;
2793 vectors->CopyInformation(vtkfield.DataArray->GetInformation());
2796 vectors->SetComponentName(2, "Z");
2798 vectors->SetComponentName(2, vtkfield.DataArray->GetComponentName(2));
2800 vectors->SetComponentName(1, vtkfield.DataArray->GetComponentName(1));
2801 vectors->SetComponentName(0, vtkfield.DataArray->GetComponentName(0));
2803 int tuplesize = (ncomp > 3? ncomp: 3);
2804 double* tuple = new double[tuplesize];
2805 for(int tid=0; tid<tuplesize; tid++)
2810 for(vtkIdType id=0; id < vtkfield.DataArray->GetNumberOfTuples(); id++)
2812 vtkfield.DataArray->GetTuple(id, tuple);
2813 vectors->SetTuple(id, tuple);
2819 int vtkMedReader::HasMeshAnyCellSelectedFamily(vtkMedMesh* mesh)
2821 int nfam = mesh->GetNumberOfCellFamily();
2822 for (int famid = 0; famid<nfam; famid++)
2824 vtkMedFamily* fam = mesh->GetFamily(famid);
2825 if (fam->GetPointOrCell()!=vtkMedUtilities::OnCell||!this->GetFamilyStatus(
2833 int vtkMedReader::HasMeshAnyPointSelectedFamily(vtkMedMesh* mesh)
2835 int nfam = mesh->GetNumberOfCellFamily();
2836 for (int famid = 0; famid<nfam; famid++)
2838 vtkMedFamily* fam = mesh->GetFamily(famid);
2839 if (fam->GetPointOrCell()!=vtkMedUtilities::OnPoint
2840 ||!this->GetFamilyStatus(mesh, fam))
2847 void vtkMedReader::BuildSIL(vtkMutableDirectedGraph* sil)
2854 vtkSmartPointer<vtkVariantArray> childEdge=
2855 vtkSmartPointer<vtkVariantArray>::New();
2856 childEdge->InsertNextValue(0);
2858 vtkSmartPointer<vtkVariantArray> crossEdge=
2859 vtkSmartPointer<vtkVariantArray>::New();
2860 crossEdge->InsertNextValue(1);
2862 // CrossEdge is an edge linking hierarchies.
2863 vtkUnsignedCharArray* crossEdgesArray=vtkUnsignedCharArray::New();
2864 crossEdgesArray->SetName("CrossEdges");
2865 sil->GetEdgeData()->AddArray(crossEdgesArray);
2866 crossEdgesArray->Delete();
2867 vtkstd::deque<vtkstd::string> names;
2869 // Now build the hierarchy.
2870 vtkIdType rootId=sil->AddVertex();
2871 names.push_back("SIL");
2873 // Add a global entry to encode global names for the families
2874 vtkIdType globalFamilyRoot=sil->AddChild(rootId, childEdge);
2875 names.push_back("Families");
2877 // Add a global entry to encode global names for the families
2878 vtkIdType globalGroupRoot=sil->AddChild(rootId, childEdge);
2879 names.push_back("Groups");
2881 // Add the groups subtree
2882 vtkIdType groupsRoot=sil->AddChild(rootId, childEdge);
2883 names.push_back("GroupTree");
2885 // Add a global entry to encode names for the cell types
2886 vtkIdType globalEntityRoot=sil->AddChild(rootId, childEdge);
2887 names.push_back("Entity");
2889 // Add the cell types subtree
2890 vtkIdType entityTypesRoot=sil->AddChild(rootId, childEdge);
2891 names.push_back("EntityTree");
2893 // this is the map that keep added cell types
2894 map<vtkMedEntity, vtkIdType> entityMap;
2895 map<med_entity_type, vtkIdType> typeMap;
2897 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
2898 meshfileit = this->Internal->MedFiles.begin();
2899 while(meshfileit != this->Internal->MedFiles.end())
2901 vtkMedFile* file = meshfileit->second;
2904 int numMeshes=file->GetNumberOfMesh();
2905 for(int meshIndex=0; meshIndex<numMeshes; meshIndex++)
2907 vtkMedMesh* mesh = file->GetMesh(meshIndex);
2908 vtkMedGrid* grid = this->FindGridStep(mesh);
2914 for(int entityIndex=0; entityIndex<grid->GetNumberOfEntityArray(); entityIndex++)
2916 vtkMedEntityArray* array=grid->GetEntityArray(entityIndex);
2917 vtkMedEntity entity = array->GetEntity();
2919 if(entityMap.find(entity)==entityMap.end())
2921 vtkIdType entityGlobalId=sil->AddChild(globalEntityRoot, childEdge);
2922 names.push_back(vtkMedUtilities::EntityKey(entity));
2925 if(typeMap.find(entity.EntityType)==typeMap.end())
2927 typeId=sil->AddChild(entityTypesRoot, childEdge);
2928 names.push_back(vtkMedUtilities::EntityName(entity.EntityType));
2929 typeMap[entity.EntityType]=typeId;
2933 typeId=typeMap[entity.EntityType];
2935 vtkIdType entityId=sil->AddChild(typeId, childEdge);
2936 names.push_back(entity.GeometryName);
2938 sil->AddEdge(entityId, entityGlobalId, crossEdge);
2940 entityMap[entity]=entityId;
2944 vtkIdType meshGroup = sil->AddChild(groupsRoot, childEdge);
2945 names.push_back(vtkMedUtilities::SimplifyName(mesh->GetName()));
2947 // add the two OnPoint and OnCell entries, for groups and for families
2948 vtkIdType meshCellGroups = sil->AddChild(meshGroup, childEdge);
2949 names.push_back(vtkMedUtilities::OnCellName);
2951 vtkIdType meshPointGroups = sil->AddChild(meshGroup, childEdge);
2952 names.push_back(vtkMedUtilities::OnPointName);
2954 // this maps will keep all added groups on nodes/cells of this mesh
2955 map<string, vtkIdType> nodeGroupMap;
2956 map<string, vtkIdType> cellGroupMap;
2959 for(int famIndex=0; famIndex<mesh->GetNumberOfFamily(); famIndex++)
2961 vtkMedFamily* family=mesh->GetFamily(famIndex);
2963 vtkIdType globalFamilyId = sil->AddChild(globalFamilyRoot, childEdge);
2964 names.push_back(vtkMedUtilities::FamilyKey(mesh->GetName(),
2965 family->GetPointOrCell(),
2966 family->GetName()));
2968 // family with Id 0 is both on cell and on points, so add it to both
2969 map<string, vtkIdType> & groupMap=(family->GetPointOrCell()
2970 ==vtkMedUtilities::OnPoint? nodeGroupMap: cellGroupMap);
2972 vtkIdType groupRootId =
2973 (family->GetPointOrCell()==vtkMedUtilities::OnPoint?
2974 meshPointGroups : meshCellGroups);
2976 // add all the groups of this family
2977 for(vtkIdType groupIndex=0; groupIndex<family->GetNumberOfGroup();
2980 vtkMedGroup* group=family->GetGroup(groupIndex);
2982 vtkIdType familyGroupId = sil->AddChild(globalFamilyId, childEdge);
2983 names.push_back(vtkMedUtilities::FamilyKey(
2984 mesh->GetName(), family->GetPointOrCell(),
2985 family->GetName()));
2987 vtkIdType groupGlobalId;
2988 if(groupMap.find(group->GetName())==groupMap.end())
2990 vtkIdType groupLocalId;
2991 groupLocalId=sil->AddChild(groupRootId, childEdge);
2992 names.push_back(vtkMedUtilities::SimplifyName(group->GetName()));
2994 groupGlobalId=sil->AddChild(globalGroupRoot, childEdge);
2995 names.push_back(vtkMedUtilities::GroupKey(
2996 mesh->GetName(), family->GetPointOrCell(),
2998 groupMap[group->GetName()]=groupGlobalId;
3000 sil->AddEdge(groupLocalId, groupGlobalId, crossEdge);
3002 vtkIdType groupId = groupMap[group->GetName()];
3003 sil->AddEdge(familyGroupId, groupId, childEdge);
3010 // This array is used to assign names to nodes.
3011 vtkStringArray* namesArray=vtkStringArray::New();
3012 namesArray->SetName("Names");
3013 namesArray->SetNumberOfTuples(sil->GetNumberOfVertices());
3014 sil->GetVertexData()->AddArray(namesArray);
3015 namesArray->Delete();
3016 vtkstd::deque<vtkstd::string>::iterator iter;
3018 for(cc=0, iter=names.begin(); iter!=names.end(); ++iter, ++cc)
3020 namesArray->SetValue(cc, (*iter).c_str());
3024 void vtkMedReader::ClearMedSupports()
3026 this->Internal->DataSetCache.clear();
3027 //this->Internal->Med2VTKPointIndex.clear();
3029 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
3030 meshfileit = this->Internal->MedFiles.begin();
3031 while(meshfileit != this->Internal->MedFiles.end())
3033 vtkMedFile* file = meshfileit->second;
3035 int numMeshes=file->GetNumberOfMesh();
3036 for(int meshIndex=0; meshIndex<numMeshes; meshIndex++)
3038 vtkMedMesh* mesh=file->GetMesh(meshIndex);
3039 mesh->ClearMedSupports();
3042 int numProf = file->GetNumberOfProfile();
3043 for (int prof = 0; prof<numProf; prof++)
3045 vtkMedProfile* profile = file->GetProfile(prof);
3046 if (profile->GetIds()!=NULL)
3047 profile->GetIds()->Initialize();
3052 void vtkMedReader::ClearMedFields()
3054 this->Internal->FieldCache.clear();
3055 this->Internal->QuadOffsetKey.clear();
3056 this->Internal->QuadratureOffsetCache.clear();
3058 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
3059 fieldfileit = this->Internal->MedFiles.begin();
3060 while(fieldfileit != this->Internal->MedFiles.end())
3062 vtkMedFile* file = fieldfileit->second;
3065 int numFields=file->GetNumberOfField();
3066 for(int ff=0; ff<numFields; ff++)
3068 vtkMedField* field=file->GetField(ff);
3069 int nstep=field->GetNumberOfFieldStep();
3070 for(int sid=0; sid<nstep; sid++)
3072 vtkMedFieldStep* step = field->GetFieldStep(sid);
3073 for(int id=0; id<step->GetNumberOfFieldOverEntity(); id++)
3075 vtkMedFieldOverEntity * fieldOverEntity=step->GetFieldOverEntity(id);
3076 for(int pid = 0; pid < fieldOverEntity->GetNumberOfFieldOnProfile(); pid++)
3078 vtkMedFieldOnProfile* fop = fieldOverEntity->GetFieldOnProfile(pid);
3079 if(fop->GetData() != NULL)
3088 void vtkMedReader::ClearCaches(int when)
3093 this->Internal->CurrentDataSet.clear();
3094 this->Internal->DataSetCache.clear();
3095 this->Internal->FieldCache.clear();
3096 this->Internal->UsedSupports.clear();
3097 this->Internal->QuadratureOffsetCache.clear();
3098 this->Internal->QuadOffsetKey.clear();
3099 //this->Internal->Med2VTKPointIndex.clear();
3102 this->Internal->CurrentDataSet.clear();
3103 this->Internal->UsedSupports.clear();
3104 if(this->CacheStrategy==CacheNothing)
3106 this->ClearMedSupports();
3107 this->ClearMedFields();
3109 else if(this->CacheStrategy==CacheGeometry)
3111 this->ClearMedFields();
3114 case AfterCreateMedSupports:
3115 // TODO : clear the unused supports and associated cached datasets and fields
3117 case EndBuildVTKSupports:
3120 if(this->CacheStrategy==CacheNothing)
3122 this->ClearMedSupports();
3123 this->ClearMedFields();
3125 else if(this->CacheStrategy==CacheGeometry && this->AnimationMode != Modes)
3127 this->ClearMedFields();
3133 void vtkMedReader::PrintSelf(ostream& os, vtkIndent indent)
3135 PRINT_STRING(os, indent, FileName);
3136 PRINT_IVAR(os, indent, AnimationMode);
3137 PRINT_IVAR(os, indent, TimeIndexForIterations);
3138 PRINT_OBJECT(os, indent, PointFields);
3139 PRINT_OBJECT(os, indent, CellFields);
3140 PRINT_OBJECT(os, indent, QuadratureFields);
3141 PRINT_OBJECT(os, indent, ElnoFields);
3142 PRINT_OBJECT(os, indent, Groups);
3143 PRINT_OBJECT(os, indent, Entities);
3144 PRINT_IVAR(os, indent, CacheStrategy);
3145 this->Superclass::PrintSelf(os, indent);