1 // Copyright (C) 2010-2011 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 "vtkMultiTimeStepAlgorithm.h"
60 #include "vtkUnstructuredGrid.h"
62 #include "vtkPointData.h"
63 #include "vtkCellData.h"
64 #include "vtkFieldData.h"
65 #include "vtkInformationQuadratureSchemeDefinitionVectorKey.h"
66 #include "vtkQuadratureSchemeDefinition.h"
67 #include "vtkCellType.h"
68 #include "vtkCellArray.h"
69 #include "vtkDoubleArray.h"
70 #include "vtkConfigure.h"
71 #include "vtkMultiProcessController.h"
72 #include "vtkCommunicator.h"
74 #include "vtkSMDoubleVectorProperty.h"
75 #include "vtkInformationDataObjectKey.h"
89 vtkSmartPointer<vtkDataArray> DataArray;
90 vtkSmartPointer<vtkDataArray> Vectors;
91 vtkSmartPointer<vtkIdTypeArray> QuadratureIndexArray;
94 class vtkMedListOfFieldSteps : public std::list<vtkMedFieldStep*>
97 typedef int LocalizationKey;
99 class vtkMedReader::vtkMedReaderInternal
103 int CurrentPieceNumber;
105 double UpdateTimeStep;
106 vtkTimeStamp FileNameMTime;
107 vtkTimeStamp MetaDataMTime;
108 vtkTimeStamp GroupSelectionMTime;
109 vtkTimeStamp FamilySelectionMTime;
111 int RealAnimationMode;
112 vtkMedSelection* Families;
114 // this stores the aggregation of all compute steps from
115 // both meshes and fields.
116 std::map<med_float, std::set<med_int> > GlobalComputeStep;
118 // Store the vtkMutableDirectedGraph that represents links between family, groups and cell types
119 vtkMutableDirectedGraph* SIL;
121 // this map is used to keep clean data sets in the cache, without any field.
122 // for each support, store the vtkDataSet
123 map<vtkMedFamilyOnEntityOnProfile*, vtkSmartPointer<vtkDataSet> > DataSetCache;
125 // this is the current dataset for the given support.
126 map<vtkMedFamilyOnEntityOnProfile*, vtkDataSet*> CurrentDataSet;
128 // for each support, cache the VTK arrays that correspond to a given field at the given step.
129 map<vtkMedFamilyOnEntityOnProfile*, map<vtkMedFieldOnProfile*, VTKField> > FieldCache;
130 //map<vtkMedFamilyOnEntity*, map<vtkMedFieldOnProfile*, bool> > FieldMatchCache;
132 // This list keep tracks of all the currently selected supports
133 set<vtkMedFamilyOnEntityOnProfile*> UsedSupports;
135 // this map keeps for each support, the quadrature offset array so that
136 // different fields on the same support can use
137 // the same offset array, provided they use the same gauss points definitions
138 map<vtkMedFamilyOnEntityOnProfile*,
139 map<LocalizationKey, vtkSmartPointer<vtkIdTypeArray> > >
140 QuadratureOffsetCache;
142 map<vtkMedFamilyOnEntityOnProfile*,
143 map<vtkMedFieldOnProfile*, LocalizationKey> > QuadOffsetKey;
145 std::map<std::string, vtkSmartPointer<vtkMedFile> > MedFiles;
147 vtkMedReaderInternal()
149 this->SIL=vtkMutableDirectedGraph::New();
150 this->SILUpdateStamp=-1;
151 this->RealAnimationMode=vtkMedReader::PhysicalTime;
152 this->Families=vtkMedSelection::New();
153 this->FamilySelectionMTime.Modified();
154 this->GroupSelectionMTime.Modified();
156 ~vtkMedReaderInternal()
159 this->Families->Delete();
162 void ClearSupportCache(vtkMedFamilyOnEntityOnProfile* foep)
164 //this->Med2VTKPointIndex.erase(foep);
165 this->QuadratureOffsetCache.erase(foep);
166 //this->FieldMatchCache.erase(foe);
167 this->FieldCache.erase(foep);
168 this->CurrentDataSet.erase(foep);
169 this->DataSetCache.erase(foep);
172 vtkIdType GetChild(vtkIdType parent, const vtkStdString childName)
174 vtkStringArray* names=vtkStringArray::SafeDownCast(
175 this->SIL->GetVertexData()->GetArray("Names"));
178 vtkIdType nedges=this->SIL->GetOutDegree(parent);
179 for(vtkIdType id=0; id<nedges; id++)
181 vtkOutEdgeType edge=this->SIL->GetOutEdge(parent, id);
182 if(names->GetValue(edge.Target)==childName)
188 vtkIdType GetGroupId(const char* key)
190 std::string meshname, celltypename, groupname;
191 vtkMedUtilities::SplitGroupKey(key, meshname, celltypename, groupname);
192 vtkIdType root=GetChild(0, "SIL");
195 vtkIdType mesh=GetChild(root, meshname);
198 vtkIdType type=GetChild(mesh, celltypename);
201 return GetChild(type, groupname);
207 //vtkCxxRevisionMacro(vtkMedReader, "$Revision$");
208 vtkStandardNewMacro(vtkMedReader);
210 vtkMedReader::vtkMedReader()
213 this->SetNumberOfInputPorts(0);
214 this->PointFields=vtkMedSelection::New();
215 this->CellFields=vtkMedSelection::New();
216 this->QuadratureFields=vtkMedSelection::New();
217 this->ElnoFields=vtkMedSelection::New();
218 this->Entities=vtkMedSelection::New();
219 this->Groups=vtkMedSelection::New();
220 this->Frequencies=vtkMedSelection::New();
221 this->AnimationMode=Default;
222 this->TimeIndexForIterations=0;
223 this->CacheStrategy=CacheGeometry;
224 this->Internal=new vtkMedReaderInternal;
225 this->TimePrecision=0.00001;
226 this->AvailableTimes=vtkDoubleArray::New();
227 this->GenerateVectors = 0;
230 vtkMedReader::~vtkMedReader()
232 this->SetFileName(NULL);
233 this->PointFields->Delete();
234 this->CellFields->Delete();
235 this->QuadratureFields->Delete();
236 this->ElnoFields->Delete();
237 this->Entities->Delete();
238 this->Groups->Delete();
239 this->Frequencies->Delete();
240 delete this->Internal;
241 this->AvailableTimes->Delete();
244 int vtkMedReader::GetSILUpdateStamp()
246 return this->Internal->SILUpdateStamp;
249 void vtkMedReader::SetFileName(const char* fname)
252 if(fname==this->FileName)
254 if(fname&&this->FileName&&!strcmp(fname, this->FileName))
258 delete[] this->FileName;
261 size_t fnl=strlen(fname)+1;
262 char* dst=new char[fnl];
263 const char* src=fname;
278 this->Internal->MedFiles.clear();
279 this->Internal->FileNameMTime.Modified();
283 int vtkMedReader::CanReadFile(const char* fname)
285 // the factory give a driver only when it can read the file version,
286 // or it returns a NULL pointer.
287 vtkSmartPointer<vtkMedFile> file=vtkSmartPointer<vtkMedFile>::New();
288 file->SetFileName(fname);
290 if(!file->CreateDriver())
296 int vtkMedReader::RequestInformation(vtkInformation *request,
297 vtkInformationVector **inputVector, vtkInformationVector *outputVector)
299 vtkInformation* outInfo = outputVector->GetInformationObject(0);
300 outInfo->Set(vtkStreamingDemandDrivenPipeline::MAXIMUM_NUMBER_OF_PIECES(),-1);
302 if(this->Internal->MetaDataMTime <= this->Internal->FileNameMTime)
304 this->ClearCaches(Initialize);
306 vtkMedFile* file=vtkMedFile::New();
307 file->SetFileName(this->FileName);
308 this->Internal->MedFiles[this->FileName] = file;
311 std::list<vtkMedFile*> file_stack;
312 file_stack.push_back(file);
314 // if the file cannot create a driver, this means that the filename is not
315 // valid, or that I do not know how to read this file.
316 while (file_stack.size() > 0)
318 vtkMedFile* file = file_stack.front();
319 file_stack.pop_front();
320 // This reads information from disk
321 file->ReadInformation();
323 // add all files linked to in the current file to the files to read.
324 for (int linkid=0; linkid<file->GetNumberOfLink(); linkid++)
326 vtkMedLink* link = file->GetLink(linkid);
327 const char* filename = link->GetFullLink(file->GetFileName());
328 if(this->Internal->MedFiles.find(filename) == this->Internal->MedFiles.end())
330 vtkMedFile* newfile = vtkMedFile::New();
331 newfile->SetFileName(filename);
332 this->Internal->MedFiles[filename] = newfile;
333 file_stack.push_back(newfile);
339 // This computes some meta information, like which field use which
340 // support, but do not read large data from disk.
343 // This computes the initial global id of each cell type.
344 this->InitializeCellGlobalIds();
346 this->ClearSelections();
348 this->BuildSIL(this->Internal->SIL);
349 this->Internal->SILUpdateStamp++;
351 this->GatherComputeSteps();
353 this->Internal->MetaDataMTime.Modified();
356 outInfo->Set(vtkDataObject::SIL(), this->Internal->SIL);
357 request->AppendUnique(vtkExecutive::KEYS_TO_COPY(),
358 vtkDataObject::SIL());
359 request->AppendUnique(vtkExecutive::KEYS_TO_COPY(),
360 vtkMedUtilities::BLOCK_NAME());
361 this->AdvertiseTime(outInfo);
365 int vtkMedReader::RequestData(vtkInformation *request,
366 vtkInformationVector **inputVector, vtkInformationVector *outputVector)
368 if(this->FileName==NULL)
370 vtkWarningMacro( << "FileName must be set and meta data loaded");
374 vtkInformation *info=outputVector->GetInformationObject(0);
376 vtkMultiBlockDataSet *output=vtkMultiBlockDataSet::SafeDownCast(info->Get(
377 vtkDataObject::DATA_OBJECT()));
379 if (info->Has(vtkStreamingDemandDrivenPipeline::UPDATE_NUMBER_OF_PIECES()))
381 this->Internal->NumberOfPieces=info->Get(
382 vtkStreamingDemandDrivenPipeline::UPDATE_NUMBER_OF_PIECES());
386 vtkMultiProcessController* controller =
387 vtkMultiProcessController::GetGlobalController();
390 this->Internal->NumberOfPieces=controller->GetNumberOfProcesses();
394 this->Internal->NumberOfPieces = 1;
397 if (info->Has(vtkStreamingDemandDrivenPipeline::UPDATE_PIECE_NUMBER()))
399 this->Internal->CurrentPieceNumber=info->Get(
400 vtkStreamingDemandDrivenPipeline::UPDATE_PIECE_NUMBER());
404 this->Internal->CurrentPieceNumber=0;
405 vtkMultiProcessController* controller =
406 vtkMultiProcessController::GetGlobalController();
409 this->Internal->CurrentPieceNumber= controller->GetLocalProcessId();
413 this->Internal->CurrentPieceNumber=0;
418 vtkStreamingDemandDrivenPipeline::UPDATE_NUMBER_OF_GHOST_LEVELS()))
420 this->Internal->GhostLevel=info->Get(
421 vtkStreamingDemandDrivenPipeline::UPDATE_NUMBER_OF_GHOST_LEVELS());
425 this->Internal->GhostLevel=0;
428 //vtkStreamingDemandDrivenPipeline::UPDATE_TIME_STEPS()
429 /* if (info->Has(vtkStreamingDemandDrivenPipeline::UPDATE_TIME_STEPS()))
431 this->Internal->UpdateTimeStep=info->Get(
432 vtkStreamingDemandDrivenPipeline::UPDATE_TIME_STEPS())[0];
436 this->Internal->UpdateTimeStep=0;
439 if (info->Has(vtkStreamingDemandDrivenPipeline::UPDATE_TIME_STEP()))
441 this->Internal->UpdateTimeStep=info->Get(
442 vtkStreamingDemandDrivenPipeline::UPDATE_TIME_STEP());
446 this->Internal->UpdateTimeStep=0;
449 this->InitializeParallelRead();
450 output->Initialize();
452 this->ChooseRealAnimationMode();
454 std::list<vtkMedDriver::FileOpen> openlist;
455 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
456 fileit = this->Internal->MedFiles.begin();
457 while(fileit != this->Internal->MedFiles.end())
459 vtkMedFile* file = fileit->second;
460 openlist.push_back(vtkMedDriver::FileOpen(file->GetMedDriver()));
464 // clear the dataset cache of unneeded geometry
465 this->ClearCaches(StartRequest);
467 // This call create the vtkMedSupports, but do not create the corresponding vtkDataSet;
468 this->CreateMedSupports();
469 this->ClearCaches(AfterCreateMedSupports);
470 // This call creates the actual vtkDataSet that corresponds to each support
473 int maxprogress=2*this->Internal->UsedSupports.size();
476 for(set<vtkMedFamilyOnEntityOnProfile*>::iterator it=
477 this->Internal->UsedSupports.begin(); it
478 !=this->Internal->UsedSupports.end(); it++)
481 vtkMedFamilyOnEntityOnProfile* foep = *it;
482 sstr<<"Support : "<<vtkMedUtilities::SimplifyName(
483 foep->GetFamilyOnEntity()->GetFamily()->GetName());
484 this->SetProgressText(sstr.str().c_str());
485 int doBuildSupportField = 1;
487 this->BuildVTKSupport(foep, doBuildSupportField);
488 this->UpdateProgress((float) progress/((float) maxprogress-1));
493 this->ClearCaches(EndBuildVTKSupports);
494 // This call maps the fields to the supports
495 for(set<vtkMedFamilyOnEntityOnProfile*>::iterator it=
496 this->Internal->UsedSupports.begin(); it
497 !=this->Internal->UsedSupports.end(); it++)
499 vtkMedFamilyOnEntityOnProfile* foep = *it;
500 if((foep->GetValid() == 0) && (this->Internal->NumberOfPieces == 1))
503 sstr<<"Loading fields on "<<vtkMedUtilities::SimplifyName(
504 foep->GetFamilyOnEntity()->GetFamily()->GetName());
505 this->SetProgressText(sstr.str().c_str());
507 this->MapFieldsOnSupport(*it, doMapField);
508 this->UpdateProgress((float) progress/((float) maxprogress-1));
513 // This call clean up caches (what is actually done depends of the CacheStrategy)
514 this->ClearCaches(EndRequest);
518 void vtkMedReader::InitializeCellGlobalIds()
520 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
521 fileit = this->Internal->MedFiles.begin();
522 while(fileit != this->Internal->MedFiles.end())
524 vtkMedFile* file = fileit->second;
526 for(int m=0; m<file->GetNumberOfMesh(); m++)
528 vtkMedMesh* mesh = file->GetMesh(m);
529 med_int nstep = mesh->GetNumberOfGridStep();
530 for(med_int stepid=0; stepid<nstep; stepid++)
532 vtkMedGrid* grid = mesh->GetGridStep(stepid);
533 grid->InitializeCellGlobalIds();
539 // Method to create the filters for the MED parallel read functions
540 // It is defined here as we have all information for initialization
541 void vtkMedReader::InitializeParallelRead()
543 // If there is only one process for reading no need to enter here
544 if (this->Internal->NumberOfPieces <= 1)
549 // FIRST: Generate filters for the cells
551 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
552 meshfit = this->Internal->MedFiles.begin();
553 while(meshfit != this->Internal->MedFiles.end())
555 vtkMedFile* meshfile = meshfit->second;
557 med_idt pFileID = meshfile->GetMedDriver()->GetParallelFileId();
559 for(int mid=0; mid<meshfile->GetNumberOfMesh(); mid++)
561 vtkMedMesh* mesh = meshfile->GetMesh(mid);
562 for(int gid=0; gid<mesh->GetNumberOfGridStep(); gid++)
564 vtkMedGrid* grid = mesh->GetGridStep(gid);
565 // read point family data and create EntityArrays
567 for(int eid=0; eid < grid->GetNumberOfEntityArray(); eid++)
569 vtkMedEntityArray* array = grid->GetEntityArray(eid);
571 // Next continue is to avoid to create filters for the
572 // points, at the moment we charge the points in all nodes
573 if (array->GetEntity().GeometryType == MED_POINT1) // !MED_NODE
576 med_int nbofconstituentpervalue = vtkMedUtilities::GetNumberOfNodes(
577 array->GetEntity().GeometryType);
578 if (nbofconstituentpervalue == -1)
579 vtkErrorMacro("Still not implemented for MED_POLYGON and MED_POLYHEDRON"); // à gerer
581 // Calculating block sizes
582 int nEntity = array->GetNumberOfEntity();
583 int block_size = ( nEntity / this->Internal->NumberOfPieces );
584 med_size start = block_size * this->Internal->CurrentPieceNumber + 1;
585 med_size stride = block_size;
587 med_size blocksize = block_size;
588 med_size lastblocksize = (nEntity % this->Internal->NumberOfPieces);
589 if ((this->Internal->NumberOfPieces ==
590 this->Internal->CurrentPieceNumber+1) && (lastblocksize != 0))
592 blocksize += lastblocksize;
593 stride += lastblocksize;
597 vtkMedFilter *filter = vtkMedFilter::New();
598 filter->SetFilterSizes( start, stride, count, blocksize, lastblocksize );
599 array->SetFilter(filter);
605 // SECOND: Filters for the Fields
607 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
609 // link the FieldOnProfile with the profiles
610 fieldfileit = this->Internal->MedFiles.begin();
611 while(fieldfileit != this->Internal->MedFiles.end())
613 vtkMedFile* fieldfile = fieldfileit->second;
615 med_idt pFileID = fieldfile->GetMedDriver()->GetParallelFileId();
617 for(int fid=0; fid<fieldfile->GetNumberOfField(); fid++)
619 vtkMedField* field = fieldfile->GetField(fid);
621 if (field->GetFieldType() == vtkMedField::CellField)
623 for(int sid = 0; sid< field->GetNumberOfFieldStep(); sid++)
625 vtkMedFieldStep* step = field->GetFieldStep(sid);
627 for(int foeid = 0; foeid < step->GetNumberOfFieldOverEntity(); foeid++)
628 // TODO : seul le premier pas de temps est dispo au debut
630 vtkMedFieldOverEntity* fieldOverEntity = step->GetFieldOverEntity(foeid);
632 for(int fopid = 0; fopid < fieldOverEntity->GetNumberOfFieldOnProfile(); fopid++)
634 vtkMedFieldOnProfile* fop = fieldOverEntity->GetFieldOnProfile(fopid);
635 // Here implement the filters as before:
636 // 1- Modify vtkMedFieldOnProfile to contain a filter
637 // 2- Create the filters here only if they are on CELLs (use GetFieldType)
638 med_int nbofconstituentpervalue = field->GetNumberOfComponent();
640 int nVectors = fop->GetNumberOfValues();
642 int block_size = ( nVectors / this->Internal->NumberOfPieces );
643 int start = block_size * this->Internal->CurrentPieceNumber + 1;
644 int stride = block_size;
646 int blocksize = block_size;
647 int lastblocksize = (nVectors % this->Internal->NumberOfPieces);
648 if ((this->Internal->NumberOfPieces ==
649 this->Internal->CurrentPieceNumber+1) && (lastblocksize != 0))
651 blocksize += lastblocksize;
652 stride += lastblocksize;
656 vtkMedFilter *filter = vtkMedFilter::New();
657 filter->SetFilterSizes( start, stride, count, blocksize, lastblocksize );
658 fop->SetFilter(filter);
668 void vtkMedReader::LinkMedInfo()
670 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
672 // link the FieldOnProfile with the profiles
673 fieldfileit = this->Internal->MedFiles.begin();
674 while(fieldfileit != this->Internal->MedFiles.end())
676 vtkMedFile* fieldfile = fieldfileit->second;
679 for(int fid=0; fid<fieldfile->GetNumberOfField(); fid++)
681 vtkMedField* field = fieldfile->GetField(fid);
683 for(int sid = 0; sid< field->GetNumberOfFieldStep(); sid++)
685 vtkMedFieldStep* step = field->GetFieldStep(sid);
687 for(int foeid = 0; foeid < step->GetNumberOfFieldOverEntity(); foeid++)
689 vtkMedFieldOverEntity* fieldOverEntity = step->GetFieldOverEntity(foeid);
691 for(int fopid = 0; fopid < fieldOverEntity->GetNumberOfFieldOnProfile(); fopid++)
693 vtkMedFieldOnProfile* fop = fieldOverEntity->GetFieldOnProfile(fopid);
695 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
696 profilefileit = this->Internal->MedFiles.begin();
697 while(profilefileit != this->Internal->MedFiles.end() && fop->GetProfile() == NULL)
699 vtkMedFile* profilefile = profilefileit->second;
702 for(int pid = 0; pid < profilefile->GetNumberOfProfile(); pid++)
704 vtkMedProfile *profile = profilefile->GetProfile(pid);
705 if(strcmp(profile->GetName(), fop->GetProfileName()) == 0)
707 fop->SetProfile(profile);
718 // first, add a familyOnEntityOnProfile to all FamilyOnEntity with a NULL
719 // profile. This is used if no field is mapped to this support.
720 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
721 meshfit = this->Internal->MedFiles.begin();
722 while(meshfit != this->Internal->MedFiles.end())
724 vtkMedFile* meshfile = meshfit->second;
727 for(int mid=0; mid<meshfile->GetNumberOfMesh(); mid++)
729 vtkMedMesh* mesh = meshfile->GetMesh(mid);
731 for(int gid=0; gid<mesh->GetNumberOfGridStep(); gid++)
733 vtkMedGrid* grid = mesh->GetGridStep(gid);
734 // read point family data and create EntityArrays
736 for(int eid=0; eid < grid->GetNumberOfEntityArray(); eid++)
738 vtkMedEntityArray* array = grid->GetEntityArray(eid);
740 for(int fid=0; fid < array->GetNumberOfFamilyOnEntity(); fid++)
742 vtkMedFamilyOnEntity* foe = array->GetFamilyOnEntity(fid);
743 if(foe->GetFamilyOnEntityOnProfile((vtkMedProfile*)NULL) == NULL)
745 vtkMedFamilyOnEntityOnProfile* foep =
746 vtkMedFamilyOnEntityOnProfile::New();
747 foep->SetFamilyOnEntity(foe);
748 foep->SetProfile(NULL);
749 foe->AddFamilyOnEntityOnProfile(foep);
758 fieldfileit = this->Internal->MedFiles.begin();
759 while(fieldfileit != this->Internal->MedFiles.end())
761 vtkMedFile* fieldfile = fieldfileit->second;
764 for(int fieldid=0; fieldid < fieldfile->GetNumberOfField(); fieldid++)
766 vtkMedField* field = fieldfile->GetField(fieldid);
768 for(int fstepid=0; fstepid < field->GetNumberOfFieldStep(); fstepid++)
770 vtkMedFieldStep* step = field->GetFieldStep(fstepid);
772 vtkMedComputeStep meshcs = step->GetMeshComputeStep();
774 for(int foeid=0; foeid<step->GetNumberOfFieldOverEntity() ;foeid++)
776 vtkMedFieldOverEntity* fieldOverEntity = step->GetFieldOverEntity(foeid);
777 const vtkMedEntity& fieldentity = fieldOverEntity->GetEntity();
780 fopid < fieldOverEntity->GetNumberOfFieldOnProfile(); fopid++)
782 vtkMedFieldOnProfile* fop =
783 fieldOverEntity->GetFieldOnProfile(fopid);
785 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
786 meshfileit = this->Internal->MedFiles.begin();
787 while(meshfileit != this->Internal->MedFiles.end())
789 vtkMedFile* meshfile = meshfileit->second;
792 if(field->GetLocal() == 1 && (meshfile != fieldfile))
795 for(int mid=0; mid<meshfile->GetNumberOfMesh(); mid++)
797 vtkMedMesh* mesh = meshfile->GetMesh(mid);
799 // the field must be on this mesh.
800 if(strcmp(mesh->GetName(),
801 field->GetMeshName()) != 0)
804 vtkMedGrid* grid = mesh->GetGridStep(meshcs);
807 vtkErrorMacro("the field " << field->GetName()
808 << " at step iteration:"
809 << step->GetComputeStep().IterationIt
811 << step->GetComputeStep().TimeIt
812 << " uses mesh at step "
813 << meshcs.TimeIt << " " << meshcs.IterationIt
814 << "which does not exists!");
818 for(int eid=0; eid < grid->GetNumberOfEntityArray(); eid++)
820 vtkMedEntityArray* array = grid->GetEntityArray(eid);
822 // if the support is on points,
823 // the field must also be on points
824 if(array->GetEntity().EntityType == MED_NODE &&
825 fieldentity.EntityType != MED_NODE)
828 if(array->GetEntity().EntityType != MED_NODE &&
829 fieldentity.EntityType == MED_NODE)
832 // for fields not on points, the geometry type
833 // of the support must match
834 if(array->GetEntity().EntityType != MED_NODE &&
835 array->GetEntity().GeometryType != fieldentity.GeometryType)
838 for(int fid = 0; fid < array->GetNumberOfFamilyOnEntity(); fid++)
840 vtkMedFamilyOnEntity* foe = array->GetFamilyOnEntity(fid);
841 if(foe->GetFamilyOnEntityOnProfile(fop->GetProfile()) == NULL)
843 vtkMedFamilyOnEntityOnProfile* foep =
844 vtkMedFamilyOnEntityOnProfile::New();
845 foep->SetProfile(fop->GetProfile());
846 foep->SetFamilyOnEntity(foe);
847 foe->AddFamilyOnEntityOnProfile(foep);
850 // also add the family on entity with no profile.
851 if(foe->GetFamilyOnEntityOnProfile((vtkMedProfile*)NULL) == NULL)
853 vtkMedFamilyOnEntityOnProfile* foep =
854 vtkMedFamilyOnEntityOnProfile::New();
855 foep->SetProfile(NULL);
856 foep->SetFamilyOnEntity(foe);
857 foe->AddFamilyOnEntityOnProfile(foep);
870 // Now, link localizations and interpolations
871 fieldfileit = this->Internal->MedFiles.begin();
872 while(fieldfileit != this->Internal->MedFiles.end())
874 vtkMedFile* fieldfile = fieldfileit->second;
877 for(int locid = 0; locid < fieldfile->GetNumberOfLocalization(); locid++)
879 vtkMedLocalization* loc = fieldfile->GetLocalization(locid);
881 for(int fid = 0; fid < fieldfile->GetNumberOfField() &&
882 loc->GetInterpolation() == NULL; fid++)
884 vtkMedField* field = fieldfile->GetField(fid);
885 for(int interpid = 0; interpid < field->GetNumberOfInterpolation();
888 vtkMedInterpolation* interp = field->GetInterpolation(interpid);
889 if(strcmp(loc->GetInterpolationName(),
890 interp->GetName()) == 0)
892 loc->SetInterpolation(interp);
900 // now that the interpolation is set, build the shape functions.
901 fieldfileit = this->Internal->MedFiles.begin();
902 while(fieldfileit != this->Internal->MedFiles.end())
904 vtkMedFile* fieldfile = fieldfileit->second;
907 for(int locid = 0; locid < fieldfile->GetNumberOfLocalization(); locid++)
909 vtkMedLocalization* loc = fieldfile->GetLocalization(locid);
910 loc->BuildShapeFunction();
914 // set the supportmesh pointer in the structural element
915 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
916 fileit = this->Internal->MedFiles.begin();
917 while(fileit != this->Internal->MedFiles.end())
919 vtkMedFile* file = fileit->second;
922 for(int structelemit = 0;
923 structelemit < file->GetNumberOfStructElement();
926 vtkMedStructElement* structElem =
927 file->GetStructElement(structelemit);
929 for(int supportmeshit = 0;
930 supportmeshit < file->GetNumberOfSupportMesh();
933 vtkMedMesh* supportMesh =
934 file->GetSupportMesh(supportmeshit);
936 if(strcmp(supportMesh->GetName(), structElem->GetName()) == 0 )
938 structElem->SetSupportMesh(supportMesh);
945 // set the pointer to the profile used by the constant attributes
946 fileit = this->Internal->MedFiles.begin();
947 while(fileit != this->Internal->MedFiles.end())
949 vtkMedFile* file = fileit->second;
952 for(int structelemit = 0;
953 structelemit < file->GetNumberOfStructElement();
956 vtkMedStructElement* structElem =
957 file->GetStructElement(structelemit);
959 for(int cstattit = 0; cstattit < structElem->GetNumberOfConstantAttribute(); cstattit++)
961 vtkMedConstantAttribute* cstatt = structElem->GetConstantAttribute(cstattit);
964 profit < file->GetNumberOfProfile();
967 vtkMedProfile* profile =
968 file->GetProfile(profit);
970 if(strcmp(profile->GetName(), cstatt->GetProfileName()) == 0 )
972 cstatt->SetProfile(profile);
980 meshfit = this->Internal->MedFiles.begin();
981 while(meshfit != this->Internal->MedFiles.end())
983 vtkMedFile* meshfile = meshfit->second;
986 for(int mid=0; mid<meshfile->GetNumberOfMesh(); mid++)
988 vtkMedMesh* mesh = meshfile->GetMesh(mid);
990 for(int gid=0; gid<mesh->GetNumberOfGridStep(); gid++)
992 vtkMedGrid* grid = mesh->GetGridStep(gid);
993 // read point family data and create EntityArrays
995 for(int eid=0; eid < grid->GetNumberOfEntityArray(); eid++)
997 vtkMedEntityArray* array = grid->GetEntityArray(eid);
998 if(array->GetEntity().EntityType != MED_STRUCT_ELEMENT)
1001 for(int structelemit = 0; structelemit < meshfile->GetNumberOfStructElement(); structelemit++)
1003 vtkMedStructElement* structelem = meshfile->GetStructElement(structelemit);
1004 if(structelem->GetGeometryType() == array->GetEntity().GeometryType)
1006 array->SetStructElement(structelem);
1016 int vtkMedReader::GetFrequencyArrayStatus(const char* name)
1018 return this->Frequencies->GetKeyStatus(name);
1021 void vtkMedReader::SetFrequencyArrayStatus(const char* name, int status)
1023 if(this->Frequencies->GetKeyStatus(name) == status)
1028 this->Frequencies->SetKeyStatus(name, status);
1033 int vtkMedReader::GetNumberOfFrequencyArrays()
1035 return this->Frequencies->GetNumberOfKey();
1038 const char* vtkMedReader::GetFrequencyArrayName(int index)
1040 return this->Frequencies->GetKey(index);
1045 bool operator()(pair<double, int> i, pair<double, int> j)
1047 if(i.first!=j.first)
1048 return (i.first<j.first);
1049 return i.second<j.second;
1053 vtkDoubleArray* vtkMedReader::GetAvailableTimes()
1055 this->AvailableTimes->Initialize();
1056 this->AvailableTimes->SetNumberOfComponents(1);
1058 std::set<std::string> newFrequencies;
1061 std::map<med_float, std::set<med_int> >::iterator it =
1062 this->Internal->GlobalComputeStep.begin();
1063 while(it != this->Internal->GlobalComputeStep.end())
1065 double time = it->first;
1066 this->AvailableTimes->InsertNextValue(time);
1067 string name = vtkMedUtilities::GetModeKey(tid, time, this->Internal->GlobalComputeStep.size()-1);
1068 this->Frequencies->AddKey(name.c_str());
1069 newFrequencies.insert(name);
1074 // now check if old frequencies have been removed
1075 for(int f = 0; f < this->Frequencies->GetNumberOfKey(); f++)
1077 const char* name = this->Frequencies->GetKey(f);
1078 if(newFrequencies.find(name) == newFrequencies.end())
1080 this->Frequencies->RemoveKeyByIndex(f);
1085 return this->AvailableTimes;
1088 void vtkMedReader::ChooseRealAnimationMode()
1090 if(this->AnimationMode!=Default)
1092 this->Internal->RealAnimationMode=this->AnimationMode;
1096 // if there is exactly one physical time and more than one iteration
1097 // set the animation mode to iteration, else default to physical time.
1098 if (this->Internal->GlobalComputeStep.size() == 1 &&
1099 this->Internal->GlobalComputeStep[0].size() > 1)
1101 this->Internal->RealAnimationMode=Iteration;
1105 this->Internal->RealAnimationMode=PhysicalTime;
1108 int vtkMedReader::GetEntityStatus(const vtkMedEntity& entity)
1110 if (entity.EntityType==MED_NODE)
1112 if(entity.EntityType == MED_DESCENDING_FACE
1113 || entity.EntityType == MED_DESCENDING_EDGE)
1116 return this->Entities->GetKeyStatus(vtkMedUtilities::EntityKey(entity).c_str());
1119 int vtkMedReader::GetFamilyStatus(vtkMedMesh* mesh, vtkMedFamily* family)
1124 if(this->Internal->GroupSelectionMTime > this->Internal->FamilySelectionMTime)
1126 this->SelectFamiliesFromGroups();
1129 int status = this->Internal->Families->GetKeyStatus(vtkMedUtilities::FamilyKey(
1130 mesh->GetName(), family->GetPointOrCell(),
1131 family->GetName()).c_str());
1136 int vtkMedReader::IsMeshSelected(vtkMedMesh* mesh)
1138 for(int fam=0; fam<mesh->GetNumberOfPointFamily(); fam++)
1140 if(this->GetFamilyStatus(mesh, mesh->GetPointFamily(fam))!=0)
1144 for(int fam=0; fam<mesh->GetNumberOfCellFamily(); fam++)
1146 if(this->GetFamilyStatus(mesh, mesh->GetCellFamily(fam))!=0)
1152 void vtkMedReader::GatherComputeSteps()
1154 this->Internal->GlobalComputeStep.clear();
1156 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
1157 fieldfileit = this->Internal->MedFiles.begin();
1158 while(fieldfileit != this->Internal->MedFiles.end())
1160 vtkMedFile* file = fieldfileit->second;
1163 // first loop over all fields to gather their compute steps
1164 for(int fieldId=0; fieldId<file->GetNumberOfField(); fieldId++)
1166 vtkMedField* field=file->GetField(fieldId);
1168 for(int stepId=0; stepId<field->GetNumberOfFieldStep(); stepId++)
1170 vtkMedFieldStep* step=field->GetFieldStep(stepId);
1171 const vtkMedComputeStep& cs = step->GetComputeStep();
1172 this->Internal->GlobalComputeStep[cs.TimeOrFrequency].insert(cs.IterationIt);
1176 // then loop over all meshes to gather their grid steps too.
1177 // for meshes, do not add the MED_UNDEF_DT time
1178 for(int meshId=0; meshId<file->GetNumberOfMesh(); meshId++)
1180 vtkMedMesh* mesh=file->GetMesh(meshId);
1182 for(int stepId=0; stepId<mesh->GetNumberOfGridStep(); stepId++)
1184 vtkMedGrid* grid=mesh->GetGridStep(stepId);
1185 const vtkMedComputeStep& cs = grid->GetComputeStep();
1186 if(cs.TimeOrFrequency != MED_UNDEF_DT || cs.TimeIt != MED_NO_DT)
1188 this->Internal->GlobalComputeStep[cs.TimeOrFrequency].insert(cs.IterationIt);
1193 if(this->Internal->GlobalComputeStep.size() == 0)
1195 this->Internal->GlobalComputeStep[MED_UNDEF_DT].insert(MED_NO_IT);
1199 int vtkMedReader::IsFieldSelected(vtkMedField* field)
1201 return this->IsPointFieldSelected(field)||this->IsCellFieldSelected(field)
1202 ||this->IsQuadratureFieldSelected(field) || this->IsElnoFieldSelected(field);
1205 int vtkMedReader::IsPointFieldSelected(vtkMedField* field)
1207 return field->GetFieldType()==vtkMedField::PointField
1208 &&this->GetPointFieldArrayStatus(vtkMedUtilities::SimplifyName(
1209 field->GetName()).c_str());
1212 int vtkMedReader::IsCellFieldSelected(vtkMedField* field)
1214 return field->GetFieldType()==vtkMedField::CellField
1215 &&this->GetCellFieldArrayStatus(vtkMedUtilities::SimplifyName(
1216 field->GetName()).c_str());
1219 vtkMedProfile* vtkMedReader::GetProfile(const char* pname)
1221 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
1222 fileit = this->Internal->MedFiles.begin();
1223 while(fileit != this->Internal->MedFiles.end())
1225 vtkMedFile* file = fileit->second;
1227 vtkMedProfile* profile = file->GetProfile(pname);
1234 vtkMedLocalization* vtkMedReader::GetLocalization(const char* lname)
1236 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
1237 fileit = this->Internal->MedFiles.begin();
1238 while(fileit != this->Internal->MedFiles.end())
1240 vtkMedFile* file = fileit->second;
1242 vtkMedLocalization* loc = file->GetLocalization(lname);
1250 int vtkMedReader::GetLocalizationKey(vtkMedFieldOnProfile* fop)
1252 vtkMedLocalization* def=this->GetLocalization(fop->GetLocalizationName());
1254 // This is not a quadrature field with explicit definition.
1255 // There are two possible cases : either the intergration point is
1256 // at the center of the cell
1257 //1 quadrature point at the cell center
1259 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
1260 fileit = this->Internal->MedFiles.begin();
1261 while(fileit != this->Internal->MedFiles.end())
1263 vtkMedFile* file = fileit->second;
1266 if(def && def->GetParentFile() == file)
1267 return nloc + def->GetMedIterator() - 1;
1269 nloc += file->GetNumberOfLocalization();
1273 if(fop->GetNumberOfIntegrationPoint()==1)
1274 return nloc + 1 + fop->GetParentFieldOverEntity()->GetEntity().GeometryType;
1276 // or it is an elno field (field stored on nodes of the cells,
1277 // but with discontinuities at the vertices)
1278 return -fop->GetParentFieldOverEntity()->GetEntity().GeometryType;//ELNO
1281 int vtkMedReader::IsQuadratureFieldSelected(vtkMedField* field)
1283 return field->GetFieldType()==vtkMedField::QuadratureField
1284 &&this->GetQuadratureFieldArrayStatus(vtkMedUtilities::SimplifyName(
1285 field->GetName()).c_str());
1288 int vtkMedReader::IsElnoFieldSelected(vtkMedField* field)
1290 return field->GetFieldType()==vtkMedField::ElnoField
1291 &&this->GetElnoFieldArrayStatus(vtkMedUtilities::SimplifyName(
1292 field->GetName()).c_str());
1296 // Give the animation steps to the pipeline
1297 void vtkMedReader::AdvertiseTime(vtkInformation* info)
1299 this->ChooseRealAnimationMode();
1301 if(this->Internal->RealAnimationMode==PhysicalTime)
1303 // I advertise the union of all times available
1304 // in all selected fields and meshes
1305 set<double> timeset;
1307 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
1308 fieldfileit = this->Internal->MedFiles.begin();
1309 while(fieldfileit != this->Internal->MedFiles.end())
1311 vtkMedFile* file = fieldfileit->second;
1314 // first loop over all fields to gather their compute steps
1315 for(int fieldId=0; fieldId<file->GetNumberOfField(); fieldId++)
1317 vtkMedField* field=file->GetField(fieldId);
1319 if(!this->IsFieldSelected(field))
1322 field->GatherFieldTimes(timeset);
1325 // then loop over all meshes to gather their grid steps too.
1326 for(int meshId=0; meshId<file->GetNumberOfMesh(); meshId++)
1328 vtkMedMesh* mesh=file->GetMesh(meshId);
1330 if(!this->IsMeshSelected(mesh))
1333 mesh->GatherGridTimes(timeset);
1337 if(timeset.size() > 0)
1339 // remove MED_UNDEF_DT if there are other time step
1340 if(timeset.size() > 1)
1341 timeset.erase(MED_UNDEF_DT);
1343 vector<double> times;
1344 set<double>::iterator it = timeset.begin();
1345 while(it != timeset.end())
1347 times.push_back(*it);
1350 sort(times.begin(), times.end());
1352 info->Set(vtkStreamingDemandDrivenPipeline::TIME_STEPS(), ×[0],
1354 double timeRange[2];
1355 timeRange[0]=times[0];
1356 timeRange[1]=times[times.size()-1];
1357 info->Set(vtkStreamingDemandDrivenPipeline::TIME_RANGE(), &timeRange[0],
1362 info->Remove(vtkStreamingDemandDrivenPipeline::TIME_STEPS());
1363 info->Remove(vtkStreamingDemandDrivenPipeline::TIME_RANGE());
1366 else if(this->Internal->RealAnimationMode==Iteration)
1368 // I advertise the union of all iterations available at the given
1369 // Time for all selected fields.
1370 set<med_int> iterationsets;
1371 med_float time = MED_UNDEF_DT;
1372 if(this->TimeIndexForIterations >= 0 &&
1373 this->TimeIndexForIterations <
1374 this->AvailableTimes->GetNumberOfTuples())
1376 time = this->AvailableTimes->
1377 GetValue((vtkIdType)this->TimeIndexForIterations);
1380 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
1381 fieldfileit = this->Internal->MedFiles.begin();
1382 while(fieldfileit != this->Internal->MedFiles.end())
1384 vtkMedFile* file = fieldfileit->second;
1387 for(int fieldId=0; fieldId<file->GetNumberOfField(); fieldId++)
1389 vtkMedField* field=file->GetField(fieldId);
1390 if(!this->IsFieldSelected(field))
1393 field->GatherFieldIterations(time, iterationsets);
1395 // then loop over all meshes to gather their grid steps too.
1396 for(int meshId=0; meshId<file->GetNumberOfMesh(); meshId++)
1398 vtkMedMesh* mesh=file->GetMesh(meshId);
1400 if(!this->IsMeshSelected(mesh))
1403 mesh->GatherGridIterations(time, iterationsets);
1407 if(iterationsets.size()>0)
1409 // remove MED_NO_IT if there are other available iterations.
1410 if(iterationsets.size()>1)
1411 iterationsets.erase(MED_NO_IT);
1413 vector<double> iterations;
1414 set<med_int>::iterator it=iterationsets.begin();
1415 while(it!=iterationsets.end())
1417 iterations.push_back((double)*it);
1420 sort(iterations.begin(), iterations.end());
1421 info->Set(vtkStreamingDemandDrivenPipeline::TIME_STEPS(), &iterations[0],
1423 double timeRange[2];
1424 timeRange[0]=iterations[0];
1425 timeRange[1]=iterations[iterations.size()-1];
1426 info->Set(vtkStreamingDemandDrivenPipeline::TIME_RANGE(), &timeRange[0],
1431 info->Remove(vtkStreamingDemandDrivenPipeline::TIME_STEPS());
1432 info->Remove(vtkStreamingDemandDrivenPipeline::TIME_RANGE());
1437 info->Remove(vtkStreamingDemandDrivenPipeline::TIME_STEPS());
1438 info->Remove(vtkStreamingDemandDrivenPipeline::TIME_RANGE());
1442 vtkIdType vtkMedReader::GetFrequencyIndex(double freq)
1444 return this->AvailableTimes->LookupValue(freq);
1447 int vtkMedReader::RequestDataObject(vtkInformation* request,
1448 vtkInformationVector** vtkNotUsed(inputVector), vtkInformationVector* outputVector)
1450 vtkInformation *info = outputVector->GetInformationObject(0);
1451 if (vtkMultiBlockDataSet::SafeDownCast(
1452 info->Get(vtkDataObject::DATA_OBJECT())))
1454 // The output is already created
1459 vtkMultiBlockDataSet* output=vtkMultiBlockDataSet::New();
1460 this->GetExecutive()->SetOutputData(0, output);
1462 this->GetOutputPortInformation(0)->Set(vtkDataObject::DATA_EXTENT_TYPE(),
1463 output->GetExtentType());
1468 void vtkMedReader::ClearSelections()
1470 this->PointFields->Initialize();
1471 this->CellFields->Initialize();
1472 this->QuadratureFields->Initialize();
1473 this->ElnoFields->Initialize();
1475 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
1476 fieldfileit = this->Internal->MedFiles.begin();
1477 while(fieldfileit != this->Internal->MedFiles.end())
1479 vtkMedFile* file = fieldfileit->second;
1482 for(int index=0; index < file->GetNumberOfField(); index++)
1484 vtkMedField* field = file->GetField(index);
1485 switch(field->GetFieldType())
1487 case vtkMedField::PointField :
1488 this->PointFields->AddKey(vtkMedUtilities::SimplifyName(
1489 field->GetName()).c_str());
1491 case vtkMedField::CellField :
1492 this->CellFields->AddKey(vtkMedUtilities::SimplifyName(
1493 field->GetName()).c_str());
1495 case vtkMedField::QuadratureField :
1496 this->QuadratureFields->AddKey(vtkMedUtilities::SimplifyName(
1497 field->GetName()).c_str());
1499 case vtkMedField::ElnoField :
1500 this->ElnoFields->AddKey(vtkMedUtilities::SimplifyName(
1501 field->GetName()).c_str());
1506 this->Internal->Families->Initialize();
1507 this->Groups->Initialize();
1508 for(int meshIndex=0; meshIndex < file->GetNumberOfMesh(); meshIndex++)
1510 vtkMedMesh* mesh = file->GetMesh(meshIndex);
1511 for(int famIndex=0; famIndex<mesh->GetNumberOfPointFamily(); famIndex++)
1513 vtkMedFamily* fam=mesh->GetPointFamily(famIndex);
1515 int ng=fam->GetNumberOfGroup();
1516 for(int gindex=0; gindex<ng; gindex++)
1518 vtkMedGroup* group=fam->GetGroup(gindex);
1519 string gname=vtkMedUtilities::GroupKey(mesh->GetName(),
1520 fam->GetPointOrCell(), group->GetName());
1522 this->Groups->AddKey(gname.c_str());
1523 this->Groups->SetKeyStatus(gname.c_str(), 0);
1526 for(int famIndex=0; famIndex<mesh->GetNumberOfCellFamily(); famIndex++)
1528 vtkMedFamily* fam=mesh->GetCellFamily(famIndex);
1530 int ng=fam->GetNumberOfGroup();
1531 for(int gindex=0; gindex<ng; gindex++)
1533 vtkMedGroup* group=fam->GetGroup(gindex);
1534 string gname=vtkMedUtilities::GroupKey(mesh->GetName(),
1535 fam->GetPointOrCell(), group->GetName());
1537 this->Groups->AddKey(gname.c_str());
1538 this->Groups->SetKeyStatus(gname.c_str(), 1);
1542 this->Internal->GroupSelectionMTime.Modified();
1544 for(int meshIndex=0; meshIndex< file->GetNumberOfMesh(); meshIndex++)
1546 if(file->GetMesh(meshIndex)->GetNumberOfGridStep() == 0)
1549 vtkMedGrid* grid=file->GetMesh(meshIndex)->GetGridStep(0);
1551 for(int entityIndex=0; entityIndex<grid->GetNumberOfEntityArray();
1554 vtkMedEntityArray* array=grid->GetEntityArray(entityIndex);
1555 string name=vtkMedUtilities::EntityKey(array->GetEntity());
1556 this->Entities->AddKey(name.c_str());
1563 void vtkMedReader::SelectFamiliesFromGroups()
1565 this->Internal->Families->Initialize();
1567 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
1568 meshfileit = this->Internal->MedFiles.begin();
1569 while(meshfileit != this->Internal->MedFiles.end())
1571 vtkMedFile* file = meshfileit->second;
1574 for(int meshIndex=0; meshIndex < file->GetNumberOfMesh(); meshIndex++)
1576 vtkMedMesh* mesh = file->GetMesh(meshIndex);
1577 for(int famIndex=0; famIndex<mesh->GetNumberOfFamily(); famIndex++)
1579 vtkMedFamily* fam=mesh->GetFamily(famIndex);
1580 string name=vtkMedUtilities::FamilyKey(mesh->GetName(),
1581 fam->GetPointOrCell(), fam->GetName());
1583 this->Internal->Families->SetKeyStatus(name.c_str(), 0);
1585 for(int gindex=0; gindex<fam->GetNumberOfGroup(); gindex++)
1587 vtkMedGroup* group=fam->GetGroup(gindex);
1588 string gname=vtkMedUtilities::GroupKey(mesh->GetName(),
1589 fam->GetPointOrCell(), group->GetName());
1590 int state=this->Groups->GetKeyStatus(gname.c_str());
1594 this->SetFamilyStatus(name.c_str(), 1);
1601 this->Internal->FamilySelectionMTime.Modified();
1604 int vtkMedReader::GetNumberOfPointFieldArrays()
1606 return this->PointFields->GetNumberOfKey();
1610 vtkMedReader::GetPointFieldArrayName(int index)
1612 return this->PointFields->GetKey(index);
1615 int vtkMedReader::GetPointFieldArrayStatus(const char* name)
1617 return this->PointFields->GetKeyStatus(name);
1620 void vtkMedReader::SetPointFieldArrayStatus(const char* name, int status)
1622 if(this->PointFields->KeyExists(name)&&this->PointFields->GetKeyStatus(
1626 this->PointFields->SetKeyStatus(name, status);
1631 int vtkMedReader::GetNumberOfCellFieldArrays()
1633 return this->CellFields->GetNumberOfKey();
1637 vtkMedReader::GetCellFieldArrayName(int index)
1639 return this->CellFields->GetKey(index);
1642 int vtkMedReader::GetCellFieldArrayStatus(const char* name)
1644 return this->CellFields->GetKeyStatus(name);
1647 void vtkMedReader::SetCellFieldArrayStatus(const char* name, int status)
1649 if(this->CellFields->KeyExists(name)&&this->CellFields->GetKeyStatus(
1653 this->CellFields->SetKeyStatus(name, status);
1658 int vtkMedReader::GetNumberOfQuadratureFieldArrays()
1660 return this->QuadratureFields->GetNumberOfKey();
1663 const char* vtkMedReader::GetQuadratureFieldArrayName(int index)
1665 return this->QuadratureFields->GetKey(index);
1668 int vtkMedReader::GetQuadratureFieldArrayStatus(const char* name)
1670 return this->QuadratureFields->GetKeyStatus(name);
1673 void vtkMedReader::SetQuadratureFieldArrayStatus(const char* name, int status)
1675 if(this->QuadratureFields->KeyExists(name)
1676 &&this->QuadratureFields->GetKeyStatus(name)==status)
1679 this->QuadratureFields->SetKeyStatus(name, status);
1684 int vtkMedReader::GetNumberOfElnoFieldArrays()
1686 return this->ElnoFields->GetNumberOfKey();
1689 const char* vtkMedReader::GetElnoFieldArrayName(int index)
1691 return this->ElnoFields->GetKey(index);
1694 int vtkMedReader::GetElnoFieldArrayStatus(const char* name)
1696 return this->ElnoFields->GetKeyStatus(name);
1699 void vtkMedReader::SetElnoFieldArrayStatus(const char* name, int status)
1701 if(this->ElnoFields->KeyExists(name)
1702 &&this->ElnoFields->GetKeyStatus(name)==status)
1705 this->ElnoFields->SetKeyStatus(name, status);
1710 void vtkMedReader::SetEntityStatus(const char* name, int status)
1712 if(this->Entities->KeyExists(name)&&this->Entities->GetKeyStatus(name)
1716 this->Entities->SetKeyStatus(name, status);
1721 void vtkMedReader::SetFamilyStatus(const char* name, int status)
1723 if(this->Internal->Families->KeyExists(name)
1724 &&this->Internal->Families->GetKeyStatus(name)==status)
1727 this->Internal->Families->SetKeyStatus(name, status);
1730 void vtkMedReader::SetGroupStatus(const char* name, int status)
1733 if(this->Groups->KeyExists(name)&&this->Groups->GetKeyStatus(name)
1737 this->Groups->SetKeyStatus(name, status);
1741 this->Internal->GroupSelectionMTime.Modified();
1744 int vtkMedReader::GetGroupStatus(const char* key)
1746 return this->Groups->GetKeyStatus(key);
1749 void vtkMedReader::AddQuadratureSchemeDefinition(vtkInformation* info,
1750 vtkMedLocalization* loc)
1752 if(info==NULL||loc==NULL)
1755 vtkInformationQuadratureSchemeDefinitionVectorKey *key=
1756 vtkQuadratureSchemeDefinition::DICTIONARY();
1758 vtkQuadratureSchemeDefinition* def=vtkQuadratureSchemeDefinition::New();
1759 int cellType=vtkMedUtilities::GetVTKCellType(loc->GetGeometryType());
1760 def->Initialize(cellType, vtkMedUtilities::GetNumberOfPoint(
1761 loc->GetGeometryType()), loc->GetNumberOfQuadraturePoint(),
1762 (double*)loc->GetShapeFunction()->GetVoidPointer(0),
1763 (double*)loc->GetWeights()->GetVoidPointer(0));
1764 key->Set(info, def, cellType);
1769 void vtkMedReader::LoadConnectivity(vtkMedEntityArray* array)
1771 vtkMedGrid* grid = array->GetParentGrid();
1772 array->LoadConnectivity();
1773 if (array->GetConnectivity()==MED_NODAL||vtkMedUtilities::GetDimension(
1774 array->GetEntity().GeometryType)<2
1775 || grid->GetParentMesh()->GetMeshType() == MED_STRUCTURED_MESH)
1778 vtkMedEntity subentity;
1779 subentity.EntityType = vtkMedUtilities::GetSubType(array->GetEntity().EntityType);
1781 vtkMedUnstructuredGrid* ugrid = vtkMedUnstructuredGrid::SafeDownCast(grid);
1785 for(int index=0; index<vtkMedUtilities::GetNumberOfSubEntity(
1786 array->GetEntity().GeometryType); index++)
1788 subentity.GeometryType = vtkMedUtilities::GetSubGeometry(
1789 array->GetEntity().GeometryType, index);
1790 vtkMedEntityArray* subarray=ugrid->GetEntityArray(subentity);
1794 vtkErrorMacro("DESC connectivity used, but sub types do not exist in file.");
1797 subarray->LoadConnectivity();
1801 vtkDataSet* vtkMedReader::CreateUnstructuredGridForPointSupport(
1802 vtkMedFamilyOnEntityOnProfile* foep)
1804 vtkUnstructuredGrid* vtkgrid = vtkUnstructuredGrid::New();
1805 foep->ComputeIntersection(NULL);
1806 vtkMedFamilyOnEntity* foe = foep->GetFamilyOnEntity();
1807 vtkMedGrid* medgrid=foe->GetParentGrid();
1808 vtkMedUnstructuredGrid* medugrid=vtkMedUnstructuredGrid::SafeDownCast(
1810 vtkMedCurvilinearGrid* medcgrid=vtkMedCurvilinearGrid::SafeDownCast(
1813 medgrid->LoadCoordinates();
1815 vtkIdType npts=medgrid->GetNumberOfPoints();
1817 bool shallowCopy= (medugrid != NULL || medcgrid!=NULL);
1818 if(medgrid->GetParentMesh()->GetSpaceDimension()!=3)
1824 shallowCopy = foep->CanShallowCopyPointField(NULL);
1827 vtkDataArray* coords = NULL;
1829 if(medugrid != NULL)
1830 coords = medugrid->GetCoordinates();
1831 if(medcgrid != NULL)
1832 coords = medcgrid->GetCoordinates();
1835 vtkIdType numberOfPoints;
1836 vtkPoints* points=vtkPoints::New(coords->GetDataType());
1837 vtkgrid->SetPoints(points);
1840 vtkIdTypeArray* pointGlobalIds=vtkIdTypeArray::New();
1841 pointGlobalIds->SetName("MED_POINT_ID");
1842 pointGlobalIds->SetNumberOfComponents(1);
1843 vtkgrid->GetPointData()->SetGlobalIds(pointGlobalIds);
1844 pointGlobalIds->Delete();
1848 vtkgrid->GetPoints()->SetData(coords);
1849 numberOfPoints=npts;
1851 pointGlobalIds->SetNumberOfTuples(numberOfPoints);
1852 vtkIdType* ptr=pointGlobalIds->GetPointer(0);
1853 for(int pid=0; pid<numberOfPoints; pid++)
1858 vtkIdType currentIndex=0;
1860 for(vtkIdType index=0; index<medgrid->GetNumberOfPoints(); index++)
1862 if (!foep->KeepPoint(index))
1867 double coord[3]={0.0, 0.0, 0.0};
1868 double * tuple=medgrid->GetCoordTuple(index);
1869 for(int dim=0; dim<medgrid->GetParentMesh()->GetSpaceDimension()&&dim<3; dim++)
1871 coord[dim]=tuple[dim];
1873 vtkgrid->GetPoints()->InsertPoint(currentIndex, coord);
1874 pointGlobalIds->InsertNextValue(index+1);
1877 vtkgrid->GetPoints()->Squeeze();
1878 pointGlobalIds->Squeeze();
1879 numberOfPoints=currentIndex;
1882 // now create the VTK_VERTEX cells
1883 for(vtkIdType id=0; id<numberOfPoints; id++)
1885 vtkgrid->InsertNextCell(VTK_VERTEX, 1, &id);
1889 // in this particular case, the global ids of the cells is the same as the global ids of the points.
1890 vtkgrid->GetCellData()->SetGlobalIds(vtkgrid->GetPointData()->GetGlobalIds());
1895 vtkMedGrid* vtkMedReader::FindGridStep(vtkMedMesh* mesh)
1897 if(this->Internal->RealAnimationMode == vtkMedReader::PhysicalTime)
1899 vtkMedComputeStep cs;
1900 cs.TimeOrFrequency = this->Internal->UpdateTimeStep;
1901 return mesh->FindGridStep(cs, vtkMedReader::PhysicalTime);
1903 else if(this->Internal->RealAnimationMode == vtkMedReader::Modes)
1905 vtkMedComputeStep cs;
1906 cs.IterationIt = MED_NO_IT;
1907 cs.TimeIt = MED_NO_DT;
1908 cs.TimeOrFrequency = MED_NO_DT;
1909 return mesh->FindGridStep(cs, vtkMedReader::Modes);
1913 vtkMedComputeStep cs;
1914 // the time is set by choosing its index in the global
1915 // array giving the available times : this->TimeIndexForIterations
1916 cs.TimeOrFrequency = (med_int)this->AvailableTimes->GetValue(
1917 (vtkIdType)this->TimeIndexForIterations);
1918 // the iteration is asked by the pipeline
1919 cs.IterationIt = (med_int)this->Internal->UpdateTimeStep;
1920 return mesh->FindGridStep(cs, vtkMedReader::Iteration);
1925 void vtkMedReader::CreateMedSupports()
1927 this->Internal->UsedSupports.clear();
1929 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
1930 meshfileit = this->Internal->MedFiles.begin();
1931 while(meshfileit != this->Internal->MedFiles.end())
1933 vtkMedFile* file = meshfileit->second;
1936 for(int meshIndex=0; meshIndex<file->GetNumberOfMesh(); meshIndex++)
1938 vtkMedMesh* mesh = file->GetMesh(meshIndex);
1939 vtkMedGrid* grid = this->FindGridStep(mesh);
1943 for(int entityIndex=0; entityIndex<grid->GetNumberOfEntityArray();
1946 vtkMedEntityArray* array=grid->GetEntityArray(entityIndex);
1947 if(this->GetEntityStatus(array->GetEntity())==0)
1952 file->GetMedDriver()->LoadFamilyIds(array);
1953 for(int foeIndex=0; foeIndex<array->GetNumberOfFamilyOnEntity();
1956 vtkMedFamilyOnEntity* foe=array->GetFamilyOnEntity(foeIndex);
1957 vtkMedFamily* family=foe->GetFamily();
1958 if(this->GetFamilyStatus(mesh, family)==0)
1961 // now, I look over all non-point fields to see which profiles
1962 // have to be used on points.
1963 bool selectedSupport = false;
1965 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
1966 fieldfileit = this->Internal->MedFiles.begin();
1967 while(fieldfileit != this->Internal->MedFiles.end())
1969 vtkMedFile* fieldfile = fieldfileit->second;
1972 for(int fieldId=0; fieldId<fieldfile->GetNumberOfField(); fieldId++)
1974 vtkMedField* field=fieldfile->GetField(fieldId);
1976 if (!this->IsFieldSelected(field))
1979 vtkMedListOfFieldSteps steps;
1981 this->GatherFieldSteps(field, steps);
1983 vtkMedListOfFieldSteps::iterator it=steps.begin();
1984 while(it!=steps.end())
1986 vtkMedFieldStep *step = *it;
1987 step->LoadInformation();
1990 for(int eid=0; eid<step->GetNumberOfFieldOverEntity(); eid++)
1992 vtkMedFieldOverEntity* fieldOverEntity =
1993 step->GetFieldOverEntity(eid);
1995 for(int pid = 0; pid < fieldOverEntity->GetNumberOfFieldOnProfile(); pid++)
1997 vtkMedFieldOnProfile* fop =
1998 fieldOverEntity->GetFieldOnProfile(pid);
1999 vtkMedProfile* prof = fop->GetProfile();
2000 vtkMedFamilyOnEntityOnProfile* foep =
2001 foe->GetFamilyOnEntityOnProfile(prof);
2004 this->Internal->UsedSupports.insert(foep);
2005 selectedSupport = true;
2012 // If no field use this family on entity, I nevertheless create the
2013 // support, with an empty profile.
2014 if(!selectedSupport)
2016 vtkMedFamilyOnEntityOnProfile* foep =
2017 foe->GetFamilyOnEntityOnProfile((vtkMedProfile*)NULL);
2020 foep = vtkMedFamilyOnEntityOnProfile::New();
2021 foep->SetFamilyOnEntity(foe);
2022 foep->SetProfile(NULL);
2023 foe->AddFamilyOnEntityOnProfile(foep);
2026 this->Internal->UsedSupports.insert(foep);
2034 bool vtkMedReader::BuildVTKSupport(
2035 vtkMedFamilyOnEntityOnProfile* foep,
2039 vtkMedFamilyOnEntity* foe = foep->GetFamilyOnEntity();
2042 vtkMultiProcessController* controller = vtkMultiProcessController::GetGlobalController();
2043 if (controller != NULL)
2045 numProc = controller->GetNumberOfProcesses();
2048 if ((foep->GetValid() == 0) && numProc == 1)
2053 vtkMedGrid* grid=foe->GetParentGrid();
2055 vtkMedEntityArray* array=foe->GetEntityArray();
2056 vtkMedMesh* mesh=grid->GetParentMesh();
2057 vtkSmartPointer<vtkStringArray> path = vtkSmartPointer<vtkStringArray>::New();
2058 string meshName=vtkMedUtilities::SimplifyName(mesh->GetName());
2059 path->InsertNextValue(meshName);
2060 std::string finalName;
2062 if(foe->GetPointOrCell()==vtkMedUtilities::OnPoint)
2064 path->InsertNextValue(vtkMedUtilities::OnPointName);
2065 finalName=vtkMedUtilities::SimplifyName(foe->GetFamily()->GetName());
2069 path->InsertNextValue(vtkMedUtilities::OnCellName);
2070 path->InsertNextValue(vtkMedUtilities::SimplifyName(foe->GetFamily()->GetName()));
2071 finalName=vtkMedUtilities::EntityKey(array->GetEntity());
2074 if(foep->GetProfile() != NULL)
2076 path->InsertNextValue(finalName);
2077 finalName = foep->GetProfile()->GetName();
2080 ostringstream progressBarTxt;
2081 for(int depth=0; depth<path->GetNumberOfValues(); depth++)
2083 progressBarTxt<<path->GetValue(depth)<<" ";
2085 progressBarTxt<<finalName;
2086 SetProgressText(progressBarTxt.str().c_str());
2088 vtkDataSet* cachedDataSet = NULL;
2089 if(this->Internal->DataSetCache.find(foep)!=this->Internal->DataSetCache.end())
2091 cachedDataSet = this->Internal->DataSetCache[foep];
2095 vtkDataSet* dataset = NULL;
2098 if(foe->GetPointOrCell()==vtkMedUtilities::OnPoint)
2100 dataset = this->CreateUnstructuredGridForPointSupport(foep);
2104 dataset = foep->GetFamilyOnEntity()->GetParentGrid()->
2105 CreateVTKDataSet(foep);
2114 this->Internal->DataSetCache[foep]=dataset;
2115 cachedDataSet = dataset;
2120 vtkMultiBlockDataSet* root=vtkMedUtilities::GetParent(this->GetOutput(), path);
2121 int nb=root->GetNumberOfBlocks();
2123 if(cachedDataSet != NULL)
2125 vtkDataSet* realDataSet=cachedDataSet->NewInstance();
2126 root->SetBlock(nb, realDataSet);
2127 realDataSet->Delete();
2129 root->GetMetaData(nb)->Set(vtkCompositeDataSet::NAME(), finalName.c_str());
2130 realDataSet->ShallowCopy(cachedDataSet);
2132 this->Internal->DataSetCache[foep]=cachedDataSet;
2133 this->Internal->CurrentDataSet[foep]=realDataSet;
2135 path->InsertNextValue(finalName);
2136 path->SetName("BLOCK_NAME");
2137 realDataSet->GetFieldData()->AddArray(path);
2138 realDataSet->GetInformation()->Remove(vtkMedUtilities::BLOCK_NAME());
2139 for(int depth=0; depth<path->GetNumberOfValues(); depth++)
2141 realDataSet->GetInformation()->Set(vtkMedUtilities::BLOCK_NAME(),
2142 path->GetValue(depth), depth);
2147 void vtkMedReader::MapFieldOnSupport(vtkMedFieldOnProfile* fop,
2148 vtkMedFamilyOnEntityOnProfile* foep,
2151 bool cached = false;
2153 if(this->Internal->FieldCache.find(foep)
2154 !=this->Internal->FieldCache.end())
2156 map<vtkMedFieldOnProfile*, VTKField>& fieldCache =
2157 this->Internal->FieldCache[foep];
2158 if(fieldCache.find(fop)!=fieldCache.end())
2166 this->CreateVTKFieldOnSupport(fop, foep, doCreateField);
2168 this->SetVTKFieldOnSupport(fop, foep);
2171 void vtkMedReader::MapFieldsOnSupport(vtkMedFamilyOnEntityOnProfile* foep,
2174 // now loop over all fields to map it to the created grids
2175 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
2176 fieldfileit = this->Internal->MedFiles.begin();
2177 while(fieldfileit != this->Internal->MedFiles.end())
2179 vtkMedFile* file = fieldfileit->second;
2182 for(int fieldId=0; fieldId<file->GetNumberOfField(); fieldId++)
2184 vtkMedField* field=file->GetField(fieldId);
2186 if(strcmp(foep->GetFamilyOnEntity()->
2187 GetParentGrid()->GetParentMesh()->GetName(),
2188 field->GetMeshName()) != 0)
2191 if(strcmp(foep->GetFamilyOnEntity()->
2192 GetParentGrid()->GetParentMesh()->GetName(),
2193 field->GetMeshName()) != 0)
2196 if (!this->IsFieldSelected(field))
2199 vtkMedListOfFieldSteps steps;
2201 this->GatherFieldSteps(field, steps);
2203 vtkMedListOfFieldSteps::iterator it=steps.begin();
2204 while(it!=steps.end())
2206 vtkMedFieldStep *step = *it;
2207 step->LoadInformation();
2210 for(int eid=0; eid<step->GetNumberOfFieldOverEntity(); eid++)
2212 vtkMedFieldOverEntity* fieldOverEntity = step->GetFieldOverEntity(eid);
2213 for(int pid = 0; pid < fieldOverEntity->GetNumberOfFieldOnProfile(); pid++)
2215 vtkMedFieldOnProfile* fop = fieldOverEntity->GetFieldOnProfile(pid);
2216 if(foep->CanMapField(fop))
2218 this->MapFieldOnSupport(fop, foep, doCreateField);
2227 void vtkMedReader::GatherFieldSteps(vtkMedField* field,
2228 vtkMedListOfFieldSteps& steps)
2230 if(this->Internal->RealAnimationMode == vtkMedReader::PhysicalTime)
2232 vtkMedComputeStep cs;
2233 cs.TimeOrFrequency = this->Internal->UpdateTimeStep;
2234 vtkMedFieldStep* fs =
2235 field->FindFieldStep(cs, vtkMedReader::PhysicalTime);
2237 steps.push_back(fs);
2239 else if(this->Internal->RealAnimationMode == vtkMedReader::Iteration)
2241 vtkMedComputeStep cs;
2242 cs.IterationIt = (med_int)this->Internal->UpdateTimeStep;
2243 cs.TimeOrFrequency = (med_int)this->AvailableTimes->GetValue(
2244 (vtkIdType)this->TimeIndexForIterations);
2245 vtkMedFieldStep* fs =
2246 field->FindFieldStep(cs, vtkMedReader::Iteration);
2248 steps.push_back(fs);
2252 for(int modeid = 0; modeid < this->Frequencies->GetNumberOfKey(); modeid++)
2254 if(this->Frequencies->GetKeyStatus(
2255 this->Frequencies->GetKey(modeid)) == 0)
2260 vtkMedComputeStep cs;
2261 cs.TimeOrFrequency = this->AvailableTimes->GetValue(modeid);
2262 vtkMedFieldStep* fs =
2263 field->FindFieldStep(cs, vtkMedReader::PhysicalTime);
2265 steps.push_back(fs);
2270 void vtkMedReader::SetVTKFieldOnSupport(vtkMedFieldOnProfile* fop,
2271 vtkMedFamilyOnEntityOnProfile* foep)
2273 vtkMedFamilyOnEntity* foe = foep->GetFamilyOnEntity();
2274 vtkMedFieldOverEntity* fieldOverEntity = fop->GetParentFieldOverEntity();
2275 vtkMedFieldStep* step = fieldOverEntity->GetParentStep();
2276 vtkMedField* field = step->GetParentField();
2278 const vtkMedComputeStep& cs = step->GetComputeStep();
2280 vtkDataSet* ds=this->Internal->CurrentDataSet[foep];
2283 // ds == NULL could arrive is some cases when working in parallel
2284 vtkWarningMacro( "--- vtkMedReader::SetVTKFieldOnSupport: ds is NULL !!");
2288 VTKField& vtkfield=this->Internal->FieldCache[foep][fop];
2290 std::string name=vtkMedUtilities::SimplifyName(field->GetName());
2291 std::string vectorname = name+"_Vector";
2292 if(this->AnimationMode==Modes)
2294 double freq = cs.TimeOrFrequency;
2295 int index = this->GetFrequencyIndex(freq);
2296 name += string(" ") + vtkMedUtilities::GetModeKey(index, freq,
2297 this->AvailableTimes->GetNumberOfTuples()-1);
2298 vectorname += string(" ") + vtkMedUtilities::GetModeKey(index, freq,
2299 this->AvailableTimes->GetNumberOfTuples()-1);
2302 vtkfield.DataArray->SetName(name.c_str());
2303 vtkfield.DataArray->Squeeze();
2304 if(vtkfield.Vectors != NULL)
2306 vtkfield.Vectors->SetName(vectorname.c_str());
2307 vtkfield.Vectors->Squeeze();
2309 if(vtkfield.QuadratureIndexArray!=NULL)
2311 vtkfield.QuadratureIndexArray->Squeeze();
2314 if(foe->GetPointOrCell()==vtkMedUtilities::OnCell)
2316 if(field->GetFieldType()==vtkMedField::PointField)
2318 if(vtkfield.DataArray->GetNumberOfTuples()!=ds->GetNumberOfPoints())
2320 vtkDebugMacro("the data array " << vtkfield.DataArray->GetName()
2321 << " do not have the good number of tuples");
2324 ds->GetPointData()->AddArray(vtkfield.DataArray);
2325 if(vtkfield.Vectors != NULL && this->GenerateVectors)
2327 ds->GetPointData()->AddArray(vtkfield.Vectors);
2328 ds->GetPointData()->SetActiveVectors(vtkfield.Vectors->GetName());
2330 switch (vtkfield.DataArray->GetNumberOfComponents())
2333 ds->GetPointData()->SetActiveScalars(vtkfield.DataArray->GetName());
2336 ds->GetPointData()->SetActiveVectors(vtkfield.DataArray->GetName());
2339 // if the data set is only composed of VTK_VERTEX cells,
2340 // and no field called with the same name exist on cells,
2341 // map this field to cells too
2342 if(foe->GetVertexOnly()==1&&ds->GetCellData()->GetArray(
2343 vtkfield.DataArray->GetName())==NULL)
2345 ds->GetCellData()->AddArray(vtkfield.DataArray);
2346 if(vtkfield.Vectors != NULL && this->GenerateVectors)
2348 ds->GetCellData()->AddArray(vtkfield.Vectors);
2349 ds->GetCellData()->SetActiveVectors(vtkfield.Vectors->GetName());
2351 switch (vtkfield.DataArray->GetNumberOfComponents())
2354 ds->GetCellData()->SetActiveScalars(vtkfield.DataArray->GetName());
2357 ds->GetCellData()->SetActiveVectors(vtkfield.DataArray->GetName());
2362 if(field->GetFieldType()==vtkMedField::CellField)
2364 if((this->Internal->NumberOfPieces == 1) && vtkfield.DataArray->GetNumberOfTuples()!=ds->GetNumberOfCells() )
2366 vtkDebugMacro("the data array " << vtkfield.DataArray->GetName()
2367 << " do not have the good number of tuples"
2368 << " ncell=" << ds->GetNumberOfCells()
2369 << " ntuple=" << vtkfield.DataArray->GetNumberOfTuples());
2372 // In case we are in parallel and our process does not contain the data
2373 if(ds->GetNumberOfCells()==0)
2377 ds->GetCellData()->AddArray(vtkfield.DataArray);
2379 if(vtkfield.Vectors != NULL && this->GenerateVectors)
2381 ds->GetCellData()->AddArray(vtkfield.Vectors);
2382 ds->GetCellData()->SetActiveVectors(vtkfield.Vectors->GetName());
2384 switch (vtkfield.DataArray->GetNumberOfComponents())
2387 ds->GetCellData()->SetActiveScalars(vtkfield.DataArray->GetName());
2390 ds->GetCellData()->SetActiveVectors(vtkfield.DataArray->GetName());
2393 // if the data set is only composed of VTK_VERTEX cells,
2394 // and no field called with the same name exist on points,
2395 // map this field to points too
2396 if(foe->GetVertexOnly()==1 && ds->GetPointData()->GetArray(
2397 vtkfield.DataArray->GetName())==NULL)
2399 ds->GetPointData()->AddArray(vtkfield.DataArray);
2400 if(vtkfield.Vectors != NULL && this->GenerateVectors)
2402 ds->GetPointData()->AddArray(vtkfield.Vectors);
2403 ds->GetPointData()->SetActiveVectors(vtkfield.Vectors->GetName());
2405 switch (vtkfield.DataArray->GetNumberOfComponents())
2408 ds->GetPointData()->SetActiveScalars(vtkfield.DataArray->GetName());
2411 ds->GetPointData()->SetActiveVectors(vtkfield.DataArray->GetName());
2416 if(field->GetFieldType()==vtkMedField::QuadratureField ||
2417 field->GetFieldType()==vtkMedField::ElnoField )
2419 vtkIdType ncells=ds->GetNumberOfCells();
2420 vtkIdType nid=vtkfield.QuadratureIndexArray->GetNumberOfTuples();
2421 vtkIdType nda=vtkfield.DataArray->GetNumberOfTuples();
2422 if((nid!=ncells) && (this->Internal->NumberOfPieces == 1))
2425 "There should be as many quadrature index values as there are cells");
2432 vtkfield.DataArray->SetNumberOfTuples( 0 );
2433 vtkfield.DataArray->Squeeze();
2435 if (nid>ncells) // PROBABLY NOT NECESSARY
2437 vtkfield.QuadratureIndexArray->SetNumberOfTuples(ncells);
2438 int nquad=fop->GetNumberOfIntegrationPoint();
2439 vtkfield.DataArray->SetNumberOfTuples( nquad * ds->GetNumberOfCells() );
2440 vtkfield.DataArray->Squeeze();
2442 ds->GetFieldData()->AddArray(vtkfield.DataArray);
2443 ds->GetCellData()->AddArray(vtkfield.QuadratureIndexArray);
2445 nid=vtkfield.QuadratureIndexArray->GetNumberOfTuples();
2446 nda=vtkfield.DataArray->GetNumberOfTuples();
2448 if(vtkfield.Vectors != NULL && this->GenerateVectors)
2450 ds->GetFieldData()->AddArray(vtkfield.Vectors);
2457 if(vtkfield.DataArray->GetNumberOfTuples()!=ds->GetNumberOfPoints())
2459 vtkDebugMacro("the data array " << vtkfield.DataArray->GetName() << " do not have the good number of tuples");
2462 ds->GetPointData()->AddArray(vtkfield.DataArray);
2463 if(vtkfield.Vectors != NULL && this->GenerateVectors)
2465 ds->GetPointData()->AddArray(vtkfield.Vectors);
2466 ds->GetPointData()->SetActiveVectors(vtkfield.Vectors->GetName());
2468 switch (vtkfield.DataArray->GetNumberOfComponents())
2471 ds->GetPointData()->SetActiveScalars(vtkfield.DataArray->GetName());
2474 ds->GetPointData()->SetActiveVectors(vtkfield.DataArray->GetName());
2477 // all the VTK_VERTEX created cells have the same order than the points,
2478 // I can safely map the point array to the cells in this particular case.
2479 ds->GetCellData()->AddArray(vtkfield.DataArray);
2480 if(vtkfield.Vectors != NULL && this->GenerateVectors)
2482 ds->GetCellData()->AddArray(vtkfield.Vectors);
2483 ds->GetCellData()->SetActiveVectors(vtkfield.Vectors->GetName());
2485 switch (vtkfield.DataArray->GetNumberOfComponents())
2488 ds->GetCellData()->SetActiveScalars(vtkfield.DataArray->GetName());
2491 ds->GetCellData()->SetActiveVectors(vtkfield.DataArray->GetName());
2497 void vtkMedReader::InitializeQuadratureOffsets(vtkMedFieldOnProfile* fop,
2498 vtkMedFamilyOnEntityOnProfile* foep)
2500 vtkMedFamilyOnEntity* foe = foep->GetFamilyOnEntity();
2502 vtkMedFieldOverEntity* fieldOverEntity = fop->GetParentFieldOverEntity();
2503 vtkMedFieldStep* step = fieldOverEntity->GetParentStep();
2504 vtkMedField *field= step->GetParentField();
2506 // then I compute the quadrature key if needed, and try to find the quadrature offsets.
2507 if(this->Internal->QuadratureOffsetCache.find(foep)
2508 ==this->Internal->QuadratureOffsetCache.end())
2509 this->Internal->QuadratureOffsetCache[foep]=map<LocalizationKey,
2510 vtkSmartPointer<vtkIdTypeArray> > ();
2512 map<LocalizationKey, vtkSmartPointer<vtkIdTypeArray> >& quadOffsets=
2513 this->Internal->QuadratureOffsetCache[foep];
2515 LocalizationKey quadkey=this->GetLocalizationKey(fop);
2517 if(quadOffsets.find(quadkey)!=quadOffsets.end())
2518 {// the quadrature offset array has already been created
2522 vtkIdTypeArray* qoffsets=vtkIdTypeArray::New();
2523 quadOffsets[quadkey]=qoffsets;
2527 if(field->GetFieldType() == vtkMedField::ElnoField)
2529 qoffsets->GetInformation()->Set(vtkMedUtilities::ELNO(), 1);
2532 else if(field->GetFieldType() == vtkMedField::QuadratureField)
2534 qoffsets->GetInformation()->Set(vtkMedUtilities::ELGA(), 1);
2539 sstr<<"QuadraturePointOffset";
2542 qoffsets->SetName(sstr.str().c_str());
2544 vtkSmartPointer<vtkMedLocalization> loc=
2545 this->GetLocalization(fop->GetLocalizationName());
2549 if(fop->GetNumberOfIntegrationPoint()==1)
2550 {// cell-centered fields can be seen as quadrature fields with 1
2551 // quadrature point at the center
2552 vtkMedLocalization* center=vtkMedLocalization::New();
2555 center->BuildCenter(fieldOverEntity->GetEntity().GeometryType);
2557 else if(loc == NULL && field->GetFieldType() == vtkMedField::ElnoField)
2558 {// ELNO fields have no vtkMedLocalization,
2559 // I need to create a dummy one
2560 vtkMedLocalization* elnodef=vtkMedLocalization::New();
2563 elnodef->BuildELNO(fieldOverEntity->GetEntity().GeometryType);
2567 vtkErrorMacro("Cannot find localization of quadrature field "
2568 << field->GetName());
2571 this->AddQuadratureSchemeDefinition(qoffsets->GetInformation(), loc);
2574 void vtkMedReader::CreateVTKFieldOnSupport(vtkMedFieldOnProfile* fop,
2575 vtkMedFamilyOnEntityOnProfile* foep, int doCreateField)
2577 vtkMedFamilyOnEntity* foe = foep->GetFamilyOnEntity();
2578 vtkMedFieldOverEntity* fieldOverEntity = fop->GetParentFieldOverEntity();
2579 vtkMedFieldStep* step = fieldOverEntity->GetParentStep();
2580 vtkMedField* field = step->GetParentField();
2582 if(vtkMedUnstructuredGrid::SafeDownCast(
2583 foep->GetFamilyOnEntity()->GetParentGrid()) != NULL)
2585 fop->Load(MED_COMPACT_STMODE);
2589 fop->Load(MED_GLOBAL_STMODE);
2592 vtkMedIntArray* profids=NULL;
2594 vtkMedProfile* profile=fop->GetProfile();
2599 profids=profile->GetIds();
2602 VTKField& vtkfield=this->Internal->FieldCache[foep][fop];
2604 bool shallowok = true;
2606 // for structured grid, the values are loaded globally, and cells which are
2607 // not on the profile or not on the family are blanked out.
2608 // shallow copy is therefore always possible
2609 if(vtkMedUnstructuredGrid::SafeDownCast(
2610 foep->GetFamilyOnEntity()->GetParentGrid()) != NULL)
2612 shallowok = foep->CanShallowCopy(fop);
2617 vtkfield.DataArray = fop->GetData();
2621 vtkDataArray* data=vtkMedUtilities::NewArray(field->GetDataType());
2622 vtkfield.DataArray=data;
2624 vtkfield.DataArray->SetNumberOfComponents(field->GetNumberOfComponent());
2627 for(int comp=0; comp<field->GetNumberOfComponent(); comp++)
2629 vtkfield.DataArray->SetComponentName(comp, vtkMedUtilities::SimplifyName(
2630 field->GetComponentName()->GetValue(comp)).c_str());
2633 bool createOffsets=false;
2634 if(field->GetFieldType()==vtkMedField::QuadratureField ||
2635 field->GetFieldType()==vtkMedField::ElnoField )
2637 this->InitializeQuadratureOffsets(fop, foep);
2639 LocalizationKey quadKey = this->GetLocalizationKey(fop);
2640 vtkfield.QuadratureIndexArray
2641 =this->Internal->QuadratureOffsetCache[foep][quadKey];
2642 vtkDataSet* ds = this->Internal->CurrentDataSet[foep];
2644 vtkfield.DataArray->GetInformation()->Set(
2645 vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME(),
2646 vtkfield.QuadratureIndexArray->GetName());
2647 vtkfield.QuadratureIndexArray->GetInformation()->Set(
2648 vtkAbstractArray::GUI_HIDE(), 1);
2650 if(vtkfield.QuadratureIndexArray->GetNumberOfTuples()
2651 !=ds->GetNumberOfCells())
2653 vtkfield.QuadratureIndexArray->SetNumberOfTuples(0);
2663 // build the quadrature offset array if needed
2667 int nquad=fop->GetNumberOfIntegrationPoint();
2668 noffsets=fop->GetData()->GetNumberOfTuples()/nquad;
2670 vtkDataSet* ds=this->Internal->CurrentDataSet[foep];
2671 if(noffsets!=ds->GetNumberOfCells())
2674 "number of quadrature offsets (" << noffsets << ") must "
2675 << "match the number of cells (" << ds->GetNumberOfCells() << ")!");
2678 vtkfield.QuadratureIndexArray->SetNumberOfTuples(noffsets);
2679 for(vtkIdType id=0; id<noffsets; id++)
2681 vtkfield.QuadratureIndexArray->SetValue(id, nquad*id);
2686 else if(foe->GetPointOrCell() == vtkMedUtilities::OnCell
2687 && field->GetFieldType() != vtkMedField::PointField)
2689 // Cell-centered field on cell support
2691 if (field->GetFieldType()==vtkMedField::QuadratureField ||
2692 field->GetFieldType()==vtkMedField::ElnoField)
2694 nquad=fop->GetNumberOfIntegrationPoint();
2696 vtkMedIntArray* profids=NULL;
2699 profids=profile->GetIds();
2701 vtkIdType maxIndex=fop->GetData()->GetNumberOfTuples();
2703 vtkIdType quadIndex = 0;
2705 for (vtkIdType id = 0; id<maxIndex; id++)
2707 vtkIdType realIndex = (profids!=NULL ? profids->GetValue(id)-1 : id);
2708 if (!foep->KeepCell(realIndex))
2711 if (field->GetFieldType()==vtkMedField::QuadratureField ||
2712 field->GetFieldType()==vtkMedField::ElnoField)
2714 for (int q = 0; q<nquad; q++)
2716 vtkfield.DataArray->InsertNextTuple(nquad*realIndex+q,
2721 vtkfield.QuadratureIndexArray->InsertNextValue(quadIndex);
2727 vtkfield.DataArray->InsertNextTuple(id, fop->GetData());
2730 vtkfield.DataArray->Squeeze();
2731 vtkDataSet* ds = this->Internal->CurrentDataSet[foep];
2733 else if(foe->GetPointOrCell() == vtkMedUtilities::OnCell
2734 && field->GetFieldType() == vtkMedField::PointField)
2735 {// point field on cell support
2736 vtkMedIntArray* profids=NULL;
2738 vtkMedProfile* profile=fop->GetProfile();
2744 profids=profile->GetIds();
2747 vtkIdType maxId=fop->GetData()->GetNumberOfTuples();
2749 for(vtkIdType id=0; id<maxId; id++)
2751 // if I have a profile, then I should insert the value at the position given in the profile.
2752 vtkIdType destIndex;
2755 destIndex=profids->GetValue(id)-1; // -1 because med indices start at 1
2762 if(!foep->KeepPoint(destIndex))
2764 // if I use the med2VTKIndex, then the index to insert
2765 // this value is given by the map.
2766 destIndex = foep->GetVTKPointIndex(destIndex);
2767 vtkfield.DataArray->InsertTuple(destIndex, id, fop->GetData());
2769 vtkfield.DataArray->Squeeze();
2770 }// point field on cell support
2771 else if(foe->GetPointOrCell() == vtkMedUtilities::OnPoint &&
2772 field->GetFieldType() == vtkMedField::PointField)
2775 vtkIdType maxId = fop->GetData()->GetNumberOfTuples();
2777 for(vtkIdType id=0; id<maxId; id++)
2779 vtkIdType realIndex=id;
2782 realIndex=profids->GetValue(id)-1; // -1 because med indices start at 1
2785 if(!foep->KeepPoint(realIndex))
2788 vtkfield.DataArray->InsertNextTuple(fop->GetData()->GetTuple(realIndex));
2790 vtkfield.DataArray->Squeeze();
2791 }// support on point
2793 // now generate the vector field if asked for
2794 if(this->GenerateVectors)
2796 int ncomp = vtkfield.DataArray->GetNumberOfComponents();
2797 if(ncomp > 1 && ncomp != 3)
2799 vtkDataArray* vectors = vtkfield.DataArray->NewInstance();
2800 vectors->SetNumberOfComponents(3);
2801 vectors->SetNumberOfTuples(vtkfield.DataArray->GetNumberOfTuples());
2802 vtkfield.Vectors = vectors;
2805 vectors->CopyInformation(vtkfield.DataArray->GetInformation());
2808 vectors->SetComponentName(2, "Z");
2810 vectors->SetComponentName(2, vtkfield.DataArray->GetComponentName(2));
2812 vectors->SetComponentName(1, vtkfield.DataArray->GetComponentName(1));
2813 vectors->SetComponentName(0, vtkfield.DataArray->GetComponentName(0));
2815 int tuplesize = (ncomp > 3? ncomp: 3);
2816 double* tuple = new double[tuplesize];
2817 for(int tid=0; tid<tuplesize; tid++)
2822 for(vtkIdType id=0; id < vtkfield.DataArray->GetNumberOfTuples(); id++)
2824 vtkfield.DataArray->GetTuple(id, tuple);
2825 vectors->SetTuple(id, tuple);
2831 int vtkMedReader::HasMeshAnyCellSelectedFamily(vtkMedMesh* mesh)
2833 int nfam = mesh->GetNumberOfCellFamily();
2834 for (int famid = 0; famid<nfam; famid++)
2836 vtkMedFamily* fam = mesh->GetFamily(famid);
2837 if (fam->GetPointOrCell()!=vtkMedUtilities::OnCell||!this->GetFamilyStatus(
2845 int vtkMedReader::HasMeshAnyPointSelectedFamily(vtkMedMesh* mesh)
2847 int nfam = mesh->GetNumberOfCellFamily();
2848 for (int famid = 0; famid<nfam; famid++)
2850 vtkMedFamily* fam = mesh->GetFamily(famid);
2851 if (fam->GetPointOrCell()!=vtkMedUtilities::OnPoint
2852 ||!this->GetFamilyStatus(mesh, fam))
2859 void vtkMedReader::BuildSIL(vtkMutableDirectedGraph* sil)
2866 vtkSmartPointer<vtkVariantArray> childEdge=
2867 vtkSmartPointer<vtkVariantArray>::New();
2868 childEdge->InsertNextValue(0);
2870 vtkSmartPointer<vtkVariantArray> crossEdge=
2871 vtkSmartPointer<vtkVariantArray>::New();
2872 crossEdge->InsertNextValue(1);
2874 // CrossEdge is an edge linking hierarchies.
2875 vtkUnsignedCharArray* crossEdgesArray=vtkUnsignedCharArray::New();
2876 crossEdgesArray->SetName("CrossEdges");
2877 sil->GetEdgeData()->AddArray(crossEdgesArray);
2878 crossEdgesArray->Delete();
2879 /*vtk*/std::deque</*vtk*/std::string> names;
2881 // Now build the hierarchy.
2882 vtkIdType rootId=sil->AddVertex();
2883 names.push_back("SIL");
2885 // Add a global entry to encode global names for the families
2886 vtkIdType globalFamilyRoot=sil->AddChild(rootId, childEdge);
2887 names.push_back("Families");
2889 // Add a global entry to encode global names for the families
2890 vtkIdType globalGroupRoot=sil->AddChild(rootId, childEdge);
2891 names.push_back("Groups");
2893 // Add the groups subtree
2894 vtkIdType groupsRoot=sil->AddChild(rootId, childEdge);
2895 names.push_back("GroupTree");
2897 // Add a global entry to encode names for the cell types
2898 vtkIdType globalEntityRoot=sil->AddChild(rootId, childEdge);
2899 names.push_back("Entity");
2901 // Add the cell types subtree
2902 vtkIdType entityTypesRoot=sil->AddChild(rootId, childEdge);
2903 names.push_back("EntityTree");
2905 // this is the map that keep added cell types
2906 map<vtkMedEntity, vtkIdType> entityMap;
2907 map<med_entity_type, vtkIdType> typeMap;
2909 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
2910 meshfileit = this->Internal->MedFiles.begin();
2911 while(meshfileit != this->Internal->MedFiles.end())
2913 vtkMedFile* file = meshfileit->second;
2916 int numMeshes=file->GetNumberOfMesh();
2917 for(int meshIndex=0; meshIndex<numMeshes; meshIndex++)
2919 vtkMedMesh* mesh = file->GetMesh(meshIndex);
2920 vtkMedGrid* grid = this->FindGridStep(mesh);
2926 for(int entityIndex=0; entityIndex<grid->GetNumberOfEntityArray(); entityIndex++)
2928 vtkMedEntityArray* array=grid->GetEntityArray(entityIndex);
2929 vtkMedEntity entity = array->GetEntity();
2931 if(entityMap.find(entity)==entityMap.end())
2933 vtkIdType entityGlobalId=sil->AddChild(globalEntityRoot, childEdge);
2934 names.push_back(vtkMedUtilities::EntityKey(entity));
2937 if(typeMap.find(entity.EntityType)==typeMap.end())
2939 typeId=sil->AddChild(entityTypesRoot, childEdge);
2940 names.push_back(vtkMedUtilities::EntityName(entity.EntityType));
2941 typeMap[entity.EntityType]=typeId;
2945 typeId=typeMap[entity.EntityType];
2947 vtkIdType entityId=sil->AddChild(typeId, childEdge);
2948 names.push_back(entity.GeometryName);
2950 sil->AddEdge(entityId, entityGlobalId, crossEdge);
2952 entityMap[entity]=entityId;
2956 vtkIdType meshGroup = sil->AddChild(groupsRoot, childEdge);
2957 names.push_back(vtkMedUtilities::SimplifyName(mesh->GetName()));
2959 // add the two OnPoint and OnCell entries, for groups and for families
2960 vtkIdType meshCellGroups = sil->AddChild(meshGroup, childEdge);
2961 names.push_back(vtkMedUtilities::OnCellName);
2963 vtkIdType meshPointGroups = sil->AddChild(meshGroup, childEdge);
2964 names.push_back(vtkMedUtilities::OnPointName);
2966 // this maps will keep all added groups on nodes/cells of this mesh
2967 map<string, vtkIdType> nodeGroupMap;
2968 map<string, vtkIdType> cellGroupMap;
2971 for(int famIndex=0; famIndex<mesh->GetNumberOfFamily(); famIndex++)
2973 vtkMedFamily* family=mesh->GetFamily(famIndex);
2975 vtkIdType globalFamilyId = sil->AddChild(globalFamilyRoot, childEdge);
2976 names.push_back(vtkMedUtilities::FamilyKey(mesh->GetName(),
2977 family->GetPointOrCell(),
2978 family->GetName()));
2980 // family with Id 0 is both on cell and on points, so add it to both
2981 map<string, vtkIdType> & groupMap=(family->GetPointOrCell()
2982 ==vtkMedUtilities::OnPoint? nodeGroupMap: cellGroupMap);
2984 vtkIdType groupRootId =
2985 (family->GetPointOrCell()==vtkMedUtilities::OnPoint?
2986 meshPointGroups : meshCellGroups);
2988 // add all the groups of this family
2989 for(vtkIdType groupIndex=0; groupIndex<family->GetNumberOfGroup();
2992 vtkMedGroup* group=family->GetGroup(groupIndex);
2994 vtkIdType familyGroupId = sil->AddChild(globalFamilyId, childEdge);
2995 names.push_back(vtkMedUtilities::FamilyKey(
2996 mesh->GetName(), family->GetPointOrCell(),
2997 family->GetName()));
2999 vtkIdType groupGlobalId;
3000 if(groupMap.find(group->GetName())==groupMap.end())
3002 vtkIdType groupLocalId;
3003 groupLocalId=sil->AddChild(groupRootId, childEdge);
3004 names.push_back(vtkMedUtilities::SimplifyName(group->GetName()));
3006 groupGlobalId=sil->AddChild(globalGroupRoot, childEdge);
3007 names.push_back(vtkMedUtilities::GroupKey(
3008 mesh->GetName(), family->GetPointOrCell(),
3010 groupMap[group->GetName()]=groupGlobalId;
3012 sil->AddEdge(groupLocalId, groupGlobalId, crossEdge);
3014 vtkIdType groupId = groupMap[group->GetName()];
3015 sil->AddEdge(familyGroupId, groupId, childEdge);
3022 // This array is used to assign names to nodes.
3023 vtkStringArray* namesArray=vtkStringArray::New();
3024 namesArray->SetName("Names");
3025 namesArray->SetNumberOfTuples(sil->GetNumberOfVertices());
3026 sil->GetVertexData()->AddArray(namesArray);
3027 namesArray->Delete();
3028 /*vtk*/std::deque</*vtk*/std::string>::iterator iter;
3030 for(cc=0, iter=names.begin(); iter!=names.end(); ++iter, ++cc)
3032 namesArray->SetValue(cc, (*iter).c_str());
3036 void vtkMedReader::ClearMedSupports()
3038 this->Internal->DataSetCache.clear();
3039 //this->Internal->Med2VTKPointIndex.clear();
3041 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
3042 meshfileit = this->Internal->MedFiles.begin();
3043 while(meshfileit != this->Internal->MedFiles.end())
3045 vtkMedFile* file = meshfileit->second;
3047 int numMeshes=file->GetNumberOfMesh();
3048 for(int meshIndex=0; meshIndex<numMeshes; meshIndex++)
3050 vtkMedMesh* mesh=file->GetMesh(meshIndex);
3051 mesh->ClearMedSupports();
3054 int numProf = file->GetNumberOfProfile();
3055 for (int prof = 0; prof<numProf; prof++)
3057 vtkMedProfile* profile = file->GetProfile(prof);
3058 if (profile->GetIds()!=NULL)
3059 profile->GetIds()->Initialize();
3064 void vtkMedReader::ClearMedFields()
3066 this->Internal->FieldCache.clear();
3067 this->Internal->QuadOffsetKey.clear();
3068 this->Internal->QuadratureOffsetCache.clear();
3070 std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
3071 fieldfileit = this->Internal->MedFiles.begin();
3072 while(fieldfileit != this->Internal->MedFiles.end())
3074 vtkMedFile* file = fieldfileit->second;
3077 int numFields=file->GetNumberOfField();
3078 for(int ff=0; ff<numFields; ff++)
3080 vtkMedField* field=file->GetField(ff);
3081 int nstep=field->GetNumberOfFieldStep();
3082 for(int sid=0; sid<nstep; sid++)
3084 vtkMedFieldStep* step = field->GetFieldStep(sid);
3085 for(int id=0; id<step->GetNumberOfFieldOverEntity(); id++)
3087 vtkMedFieldOverEntity * fieldOverEntity=step->GetFieldOverEntity(id);
3088 for(int pid = 0; pid < fieldOverEntity->GetNumberOfFieldOnProfile(); pid++)
3090 vtkMedFieldOnProfile* fop = fieldOverEntity->GetFieldOnProfile(pid);
3091 if(fop->GetData() != NULL)
3100 void vtkMedReader::ClearCaches(int when)
3105 this->Internal->CurrentDataSet.clear();
3106 this->Internal->DataSetCache.clear();
3107 this->Internal->FieldCache.clear();
3108 this->Internal->UsedSupports.clear();
3109 this->Internal->QuadratureOffsetCache.clear();
3110 this->Internal->QuadOffsetKey.clear();
3111 //this->Internal->Med2VTKPointIndex.clear();
3114 this->Internal->CurrentDataSet.clear();
3115 this->Internal->UsedSupports.clear();
3116 if(this->CacheStrategy==CacheNothing)
3118 this->ClearMedSupports();
3119 this->ClearMedFields();
3121 else if(this->CacheStrategy==CacheGeometry)
3123 this->ClearMedFields();
3126 case AfterCreateMedSupports:
3127 // TODO : clear the unused supports and associated cached datasets and fields
3129 case EndBuildVTKSupports:
3132 if(this->CacheStrategy==CacheNothing)
3134 this->ClearMedSupports();
3135 this->ClearMedFields();
3137 else if(this->CacheStrategy==CacheGeometry && this->AnimationMode != Modes)
3139 this->ClearMedFields();
3145 void vtkMedReader::PrintSelf(ostream& os, vtkIndent indent)
3147 PRINT_STRING(os, indent, FileName);
3148 PRINT_IVAR(os, indent, AnimationMode);
3149 PRINT_IVAR(os, indent, TimeIndexForIterations);
3150 PRINT_OBJECT(os, indent, PointFields);
3151 PRINT_OBJECT(os, indent, CellFields);
3152 PRINT_OBJECT(os, indent, QuadratureFields);
3153 PRINT_OBJECT(os, indent, ElnoFields);
3154 PRINT_OBJECT(os, indent, Groups);
3155 PRINT_OBJECT(os, indent, Entities);
3156 PRINT_IVAR(os, indent, CacheStrategy);
3157 this->Superclass::PrintSelf(os, indent);