Back to index

salome-paravis  6.5.0
vtkMedReader.cxx
Go to the documentation of this file.
00001 // Copyright (C) 2010-2012  CEA/DEN, EDF R&D
00002 //
00003 // This library is free software; you can redistribute it and/or
00004 // modify it under the terms of the GNU Lesser General Public
00005 // License as published by the Free Software Foundation; either
00006 // version 2.1 of the License.
00007 //
00008 // This library is distributed in the hope that it will be useful,
00009 // but WITHOUT ANY WARRANTY; without even the implied warranty of
00010 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00011 // Lesser General Public License for more details.
00012 //
00013 // You should have received a copy of the GNU Lesser General Public
00014 // License along with this library; if not, write to the Free Software
00015 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
00016 //
00017 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
00018 //
00019 
00020 #include "vtkMedReader.h"
00021 
00022 #include "vtkMedFile.h"
00023 #include "vtkMedDriver.h"
00024 #include "vtkMedUtilities.h"
00025 #include "vtkMedFamily.h"
00026 #include "vtkMedGroup.h"
00027 #include "vtkMedMesh.h"
00028 #include "vtkMedUnstructuredGrid.h"
00029 #include "vtkMedCurvilinearGrid.h"
00030 #include "vtkMedRegularGrid.h"
00031 #include "vtkMedEntityArray.h"
00032 #include "vtkMedField.h"
00033 #include "vtkMedFieldStep.h"
00034 #include "vtkMedFieldOverEntity.h"
00035 #include "vtkMedProfile.h"
00036 #include "vtkMedIntArray.h"
00037 #include "vtkMedLocalization.h"
00038 #include "vtkMedFamilyOnEntity.h"
00039 #include "vtkMedFieldOnProfile.h"
00040 #include "vtkMedSelection.h"
00041 #include "vtkMedFamilyOnEntityOnProfile.h"
00042 #include "vtkMedLink.h"
00043 #include "vtkMedInterpolation.h"
00044 #include "vtkMedStructElement.h"
00045 #include "vtkMedConstantAttribute.h"
00046 
00047 #include "vtkObjectFactory.h"
00048 #include "vtkMutableDirectedGraph.h"
00049 #include "vtkStringArray.h"
00050 #include "vtkDataSetAttributes.h"
00051 #include "vtkUnsignedCharArray.h"
00052 #include "vtkInformation.h"
00053 #include "vtkInformationVector.h"
00054 #include "vtkSmartPointer.h"
00055 #include "vtkVariantArray.h"
00056 #include "vtkMultiBlockDataSet.h"
00057 #include "vtkExecutive.h"
00058 #include "vtkStreamingDemandDrivenPipeline.h"
00059 #include "vtkUnstructuredGrid.h"
00060 #include "vtkMath.h"
00061 #include "vtkPointData.h"
00062 #include "vtkCellData.h"
00063 #include "vtkFieldData.h"
00064 #include "vtkInformationQuadratureSchemeDefinitionVectorKey.h"
00065 #include "vtkQuadratureSchemeDefinition.h"
00066 #include "vtkCellType.h"
00067 #include "vtkCellArray.h"
00068 #include "vtkDoubleArray.h"
00069 #include "vtkConfigure.h"
00070 #include "vtkMultiProcessController.h"
00071 #include "vtkCommunicator.h"
00072 
00073 #include "vtkSMDoubleVectorProperty.h"
00074 #include "vtkInformationDataObjectKey.h"
00075 
00076 #include <deque>
00077 #include <map>
00078 #include <string>
00079 #include <sstream>
00080 #include <vector>
00081 #include <list>
00082 #include <set>
00083 #include <algorithm>
00084 using namespace std;
00085 
00086 struct VTKField
00087 {
00088   vtkSmartPointer<vtkDataArray> DataArray;
00089   vtkSmartPointer<vtkDataArray> Vectors;
00090   vtkSmartPointer<vtkIdTypeArray> QuadratureIndexArray;
00091 };
00092 
00093 class vtkMedListOfFieldSteps : public std::list<vtkMedFieldStep*>
00094 {};
00095 
00096 typedef int LocalizationKey;
00097 
00098 class vtkMedReader::vtkMedReaderInternal
00099 {
00100 public:
00101   int NumberOfPieces;
00102   int CurrentPieceNumber;
00103   int GhostLevel;
00104   double UpdateTimeStep;
00105   vtkTimeStamp FileNameMTime;
00106   vtkTimeStamp MetaDataMTime;
00107   vtkTimeStamp GroupSelectionMTime;
00108   vtkTimeStamp FamilySelectionMTime;
00109   int SILUpdateStamp;
00110   int RealAnimationMode;
00111   vtkMedSelection* Families;
00112 
00113   // this stores the aggregation of all compute steps from
00114   // both meshes and fields.
00115   std::map<med_float, std::set<med_int> > GlobalComputeStep;
00116 
00117   // Store the vtkMutableDirectedGraph that represents links between family, groups and cell types
00118   vtkMutableDirectedGraph* SIL;
00119 
00120   // this map is used to keep clean data sets in the cache, without any field.
00121   // for each support, store the vtkDataSet
00122   map<vtkMedFamilyOnEntityOnProfile*, vtkSmartPointer<vtkDataSet> > DataSetCache;
00123 
00124   // this is the current dataset for the given support.
00125   map<vtkMedFamilyOnEntityOnProfile*, vtkDataSet*> CurrentDataSet;
00126 
00127   // for each support, cache the VTK arrays that correspond to a given field at the given step.
00128   map<vtkMedFamilyOnEntityOnProfile*, map<vtkMedFieldOnProfile*, VTKField> > FieldCache;
00129   //map<vtkMedFamilyOnEntity*, map<vtkMedFieldOnProfile*, bool> > FieldMatchCache;
00130 
00131   // This list keep tracks of all the currently selected supports
00132   set<vtkMedFamilyOnEntityOnProfile*> UsedSupports;
00133 
00134   // this map keeps for each support, the quadrature offset array so that
00135   // different fields on the same support can use
00136   // the same offset array, provided they use the same gauss points definitions
00137   map<vtkMedFamilyOnEntityOnProfile*,
00138       map<LocalizationKey, vtkSmartPointer<vtkIdTypeArray> > >
00139       QuadratureOffsetCache;
00140 
00141   map<vtkMedFamilyOnEntityOnProfile*,
00142       map<vtkMedFieldOnProfile*, LocalizationKey> > QuadOffsetKey;
00143 
00144   std::map<std::string, vtkSmartPointer<vtkMedFile> > MedFiles;
00145 
00146   vtkMedReaderInternal()
00147   {
00148     this->SIL=vtkMutableDirectedGraph::New();
00149     this->SILUpdateStamp=-1;
00150     this->RealAnimationMode=vtkMedReader::PhysicalTime;
00151     this->Families=vtkMedSelection::New();
00152     this->FamilySelectionMTime.Modified();
00153     this->GroupSelectionMTime.Modified();
00154   }
00155   ~vtkMedReaderInternal()
00156   {
00157     this->SIL->Delete();
00158     this->Families->Delete();
00159   }
00160 
00161   void ClearSupportCache(vtkMedFamilyOnEntityOnProfile* foep)
00162   {
00163     //this->Med2VTKPointIndex.erase(foep);
00164     this->QuadratureOffsetCache.erase(foep);
00165     //this->FieldMatchCache.erase(foe);
00166     this->FieldCache.erase(foep);
00167     this->CurrentDataSet.erase(foep);
00168     this->DataSetCache.erase(foep);
00169   }
00170 
00171   vtkIdType GetChild(vtkIdType parent, const vtkStdString childName)
00172   {
00173     vtkStringArray* names=vtkStringArray::SafeDownCast(
00174         this->SIL->GetVertexData()->GetArray("Names"));
00175     if(names==NULL)
00176       return -1;
00177     vtkIdType nedges=this->SIL->GetOutDegree(parent);
00178     for(vtkIdType id=0; id<nedges; id++)
00179       {
00180       vtkOutEdgeType edge=this->SIL->GetOutEdge(parent, id);
00181       if(names->GetValue(edge.Target)==childName)
00182         return edge.Target;
00183       }
00184     return -1;
00185   }
00186 
00187   vtkIdType GetGroupId(const char* key)
00188   {
00189     vtkstd::string meshname, celltypename, groupname;
00190     vtkMedUtilities::SplitGroupKey(key, meshname, celltypename, groupname);
00191     vtkIdType root=GetChild(0, "SIL");
00192     if(root==-1)
00193       return -1;
00194     vtkIdType mesh=GetChild(root, meshname);
00195     if(mesh==-1)
00196       return -1;
00197     vtkIdType type=GetChild(mesh, celltypename);
00198     if(type==-1)
00199       return -1;
00200     return GetChild(type, groupname);
00201 
00202   }
00203 
00204 };
00205 
00206 vtkCxxRevisionMacro(vtkMedReader, "$Revision: 1.1.4.15.2.2 $");
00207 vtkStandardNewMacro(vtkMedReader);
00208 
00209 vtkMedReader::vtkMedReader()
00210 {
00211   this->FileName=NULL;
00212   this->SetNumberOfInputPorts(0);
00213   this->PointFields=vtkMedSelection::New();
00214   this->CellFields=vtkMedSelection::New();
00215   this->QuadratureFields=vtkMedSelection::New();
00216   this->ElnoFields=vtkMedSelection::New();
00217   this->Entities=vtkMedSelection::New();
00218   this->Groups=vtkMedSelection::New();
00219   this->Frequencies=vtkMedSelection::New();
00220   this->AnimationMode=Default;
00221   this->TimeIndexForIterations=0;
00222   this->CacheStrategy=CacheGeometry;
00223   this->Internal=new vtkMedReaderInternal;
00224   this->TimePrecision=0.00001;
00225   this->AvailableTimes=vtkDoubleArray::New();
00226   this->GenerateVectors = 0;
00227 }
00228 
00229 vtkMedReader::~vtkMedReader()
00230 {
00231   this->SetFileName(NULL);
00232   this->PointFields->Delete();
00233   this->CellFields->Delete();
00234   this->QuadratureFields->Delete();
00235   this->ElnoFields->Delete();
00236   this->Entities->Delete();
00237   this->Groups->Delete();
00238   this->Frequencies->Delete();
00239   delete this->Internal;
00240   this->AvailableTimes->Delete();
00241 }
00242 
00243 int vtkMedReader::GetSILUpdateStamp()
00244 {
00245   return this->Internal->SILUpdateStamp;
00246 }
00247 
00248 void vtkMedReader::SetFileName(const char* fname)
00249 {
00250   int modified=0;
00251   if(fname==this->FileName)
00252     return;
00253   if(fname&&this->FileName&&!strcmp(fname, this->FileName))
00254     return;
00255   modified=1;
00256   if(this->FileName)
00257     delete[] this->FileName;
00258   if (fname)
00259     {
00260     size_t fnl=strlen(fname)+1;
00261     char* dst=new char[fnl];
00262     const char* src=fname;
00263     this->FileName=dst;
00264     do
00265       {
00266       *dst++=*src++;
00267       }
00268     while (--fnl);
00269     }
00270   else
00271     {
00272     this->FileName=0;
00273     }
00274   if (modified)
00275     {
00276     this->Modified();
00277     this->Internal->MedFiles.clear();
00278     this->Internal->FileNameMTime.Modified();
00279     }
00280 }
00281 
00282 int vtkMedReader::CanReadFile(const char* fname)
00283 {
00284   // the factory give a driver only when it can read the file version,
00285   // or it returns a NULL pointer.
00286   vtkSmartPointer<vtkMedFile> file=vtkSmartPointer<vtkMedFile>::New();
00287   file->SetFileName(fname);
00288 
00289   if(!file->CreateDriver())
00290     return 0;
00291 
00292   return 1;
00293 }
00294 
00295 int vtkMedReader::RequestInformation(vtkInformation *request,
00296     vtkInformationVector **inputVector, vtkInformationVector *outputVector)
00297 {
00298   vtkInformation* outInfo = outputVector->GetInformationObject(0);
00299   outInfo->Set(vtkStreamingDemandDrivenPipeline::MAXIMUM_NUMBER_OF_PIECES(),-1);
00300 
00301   if(this->Internal->MetaDataMTime <= this->Internal->FileNameMTime)
00302     {
00303     this->ClearCaches(Initialize);
00304 
00305     vtkMedFile* file=vtkMedFile::New();
00306     file->SetFileName(this->FileName);
00307     this->Internal->MedFiles[this->FileName] = file;
00308     file->Delete();
00309 
00310     std::list<vtkMedFile*> file_stack;
00311     file_stack.push_back(file);
00312 
00313     // if the file cannot create a driver, this means that the filename is not
00314     // valid, or that I do not know how to read this file.
00315     while (file_stack.size() > 0)
00316       {
00317       vtkMedFile* file = file_stack.front();
00318       file_stack.pop_front();
00319       // This reads information from disk
00320       file->ReadInformation();
00321 
00322       // add all files linked to in the current file to the files to read.
00323       for (int linkid=0; linkid<file->GetNumberOfLink(); linkid++)
00324         {
00325         vtkMedLink* link = file->GetLink(linkid);
00326         const char* filename = link->GetFullLink(file->GetFileName());
00327         if(this->Internal->MedFiles.find(filename) == this->Internal->MedFiles.end())
00328           {
00329           vtkMedFile* newfile = vtkMedFile::New();
00330           newfile->SetFileName(filename);
00331           this->Internal->MedFiles[filename] = newfile;
00332           file_stack.push_back(newfile);
00333           newfile->Delete();
00334           }
00335         }
00336       }
00337 
00338     // This computes some meta information, like which field use which
00339     // support, but do not read large data from disk.
00340     this->LinkMedInfo();
00341 
00342     // This computes the initial global id of each cell type.
00343     this->InitializeCellGlobalIds();
00344 
00345     this->ClearSelections();
00346 
00347     this->BuildSIL(this->Internal->SIL);
00348     this->Internal->SILUpdateStamp++;
00349 
00350     this->GatherComputeSteps();
00351 
00352     this->Internal->MetaDataMTime.Modified();
00353     }
00354 
00355   outInfo->Set(vtkDataObject::SIL(), this->Internal->SIL);
00356   request->AppendUnique(vtkExecutive::KEYS_TO_COPY(),
00357                         vtkDataObject::SIL());
00358   request->AppendUnique(vtkExecutive::KEYS_TO_COPY(),
00359                         vtkMedUtilities::BLOCK_NAME());
00360   this->AdvertiseTime(outInfo);
00361   return 1;
00362 }
00363 
00364 int vtkMedReader::RequestData(vtkInformation *request,
00365     vtkInformationVector **inputVector, vtkInformationVector *outputVector)
00366 {
00367   if(this->FileName==NULL)
00368     {
00369     vtkWarningMacro( << "FileName must be set and meta data loaded");
00370     return 0;
00371     }
00372 
00373   vtkInformation *info=outputVector->GetInformationObject(0);
00374 
00375   vtkMultiBlockDataSet *output=vtkMultiBlockDataSet::SafeDownCast(info->Get(
00376       vtkDataObject::DATA_OBJECT()));
00377 
00378   if (info->Has(vtkStreamingDemandDrivenPipeline::UPDATE_NUMBER_OF_PIECES()))
00379     {
00380     this->Internal->NumberOfPieces=info->Get(
00381         vtkStreamingDemandDrivenPipeline::UPDATE_NUMBER_OF_PIECES());
00382     }
00383   else
00384     {
00385     vtkMultiProcessController* controller =
00386               vtkMultiProcessController::GetGlobalController();
00387     if(controller)
00388       {
00389       this->Internal->NumberOfPieces=controller->GetNumberOfProcesses();
00390       }
00391     else
00392       {
00393       this->Internal->NumberOfPieces = 1;
00394       }
00395     }
00396   if (info->Has(vtkStreamingDemandDrivenPipeline::UPDATE_PIECE_NUMBER()))
00397     {
00398     this->Internal->CurrentPieceNumber=info->Get(
00399         vtkStreamingDemandDrivenPipeline::UPDATE_PIECE_NUMBER());
00400     }
00401   else
00402     {
00403     this->Internal->CurrentPieceNumber=0;
00404     vtkMultiProcessController* controller =
00405               vtkMultiProcessController::GetGlobalController();
00406     if(controller)
00407       {
00408       this->Internal->CurrentPieceNumber= controller->GetLocalProcessId();
00409       }
00410     else
00411       {
00412       this->Internal->CurrentPieceNumber=0;
00413       }
00414     }
00415 
00416   if (info->Has(
00417       vtkStreamingDemandDrivenPipeline::UPDATE_NUMBER_OF_GHOST_LEVELS()))
00418     {
00419     this->Internal->GhostLevel=info->Get(
00420         vtkStreamingDemandDrivenPipeline::UPDATE_NUMBER_OF_GHOST_LEVELS());
00421     }
00422   else
00423     {
00424     this->Internal->GhostLevel=0;
00425     }
00426 
00427   if (info->Has(vtkStreamingDemandDrivenPipeline::UPDATE_TIME_STEPS()))
00428     {
00429     this->Internal->UpdateTimeStep=info->Get(
00430         vtkStreamingDemandDrivenPipeline::UPDATE_TIME_STEPS())[0];
00431     }
00432   else
00433     {
00434     this->Internal->UpdateTimeStep=0;
00435     }
00436 
00437   this->InitializeParallelRead();
00438   output->Initialize();
00439 
00440   this->ChooseRealAnimationMode();
00441 
00442   std::list<vtkMedDriver::FileOpen> openlist;
00443   std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
00444       fileit = this->Internal->MedFiles.begin();
00445   while(fileit != this->Internal->MedFiles.end())
00446     {
00447     vtkMedFile* file = fileit->second;
00448     openlist.push_back(vtkMedDriver::FileOpen(file->GetMedDriver()));
00449     fileit++;
00450     }
00451 
00452   // clear the dataset cache of unneeded geometry
00453   this->ClearCaches(StartRequest);
00454 
00455   // This call create the vtkMedSupports, but do not create the corresponding vtkDataSet;
00456   this->CreateMedSupports();
00457   this->ClearCaches(AfterCreateMedSupports);
00458   // This call creates the actual vtkDataSet that corresponds to each support
00459   int supportId = 0;
00460   int progress=0;
00461   int maxprogress=2*this->Internal->UsedSupports.size();
00462   supportId = 0;
00463   int it_counter = 0;
00464   for(set<vtkMedFamilyOnEntityOnProfile*>::iterator it=
00465       this->Internal->UsedSupports.begin(); it
00466       !=this->Internal->UsedSupports.end(); it++)
00467     {
00468     ostringstream sstr;
00469     vtkMedFamilyOnEntityOnProfile* foep = *it;
00470     sstr<<"Support : "<<vtkMedUtilities::SimplifyName(
00471         foep->GetFamilyOnEntity()->GetFamily()->GetName());
00472     this->SetProgressText(sstr.str().c_str());
00473     int doBuildSupportField = 1;
00474     it_counter++;
00475     this->BuildVTKSupport(foep, doBuildSupportField);
00476     this->UpdateProgress((float) progress/((float) maxprogress-1));
00477     progress++;
00478     supportId++;
00479     }
00480 
00481   this->ClearCaches(EndBuildVTKSupports);
00482   // This call maps the fields to the supports
00483   for(set<vtkMedFamilyOnEntityOnProfile*>::iterator it=
00484       this->Internal->UsedSupports.begin(); it
00485       !=this->Internal->UsedSupports.end(); it++)
00486     {
00487     vtkMedFamilyOnEntityOnProfile* foep = *it;
00488     if((foep->GetValid() == 0) && (this->Internal->NumberOfPieces == 1))
00489       continue;
00490     ostringstream sstr;
00491     sstr<<"Loading fields on "<<vtkMedUtilities::SimplifyName(
00492         foep->GetFamilyOnEntity()->GetFamily()->GetName());
00493     this->SetProgressText(sstr.str().c_str());
00494     int doMapField = 1;
00495     this->MapFieldsOnSupport(*it, doMapField);
00496     this->UpdateProgress((float) progress/((float) maxprogress-1));
00497     progress++;
00498     supportId++;
00499     }
00500 
00501   // This call clean up caches (what is actually done depends of the CacheStrategy)
00502   this->ClearCaches(EndRequest);
00503   return 1;
00504 }
00505 
00506 void vtkMedReader::InitializeCellGlobalIds()
00507 {
00508   std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
00509       fileit = this->Internal->MedFiles.begin();
00510   while(fileit != this->Internal->MedFiles.end())
00511     {
00512     vtkMedFile* file = fileit->second;
00513     fileit++;
00514     for(int m=0; m<file->GetNumberOfMesh(); m++)
00515       {
00516       vtkMedMesh* mesh = file->GetMesh(m);
00517       med_int nstep = mesh->GetNumberOfGridStep();
00518       for(med_int stepid=0; stepid<nstep; stepid++)
00519         {
00520         vtkMedGrid* grid = mesh->GetGridStep(stepid);
00521         grid->InitializeCellGlobalIds();
00522         }
00523       }
00524     }
00525 }
00526 
00527 // Method to create the filters for the MED parallel read functions
00528 // It is defined here as we have all information for initialization
00529 void vtkMedReader::InitializeParallelRead()
00530 {
00531   // If there is only one process for reading no need to enter here
00532   if (this->Internal->NumberOfPieces <= 1)
00533     {
00534     return;
00535     }
00536 
00537   // FIRST: Generate filters for the cells
00538 
00539   std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
00540   meshfit = this->Internal->MedFiles.begin();
00541   while(meshfit != this->Internal->MedFiles.end())
00542     {
00543     vtkMedFile* meshfile = meshfit->second;
00544     meshfit++;
00545     med_idt pFileID = meshfile->GetMedDriver()->GetParallelFileId();
00546 
00547     for(int mid=0; mid<meshfile->GetNumberOfMesh(); mid++)
00548       {
00549       vtkMedMesh* mesh = meshfile->GetMesh(mid);
00550       for(int gid=0; gid<mesh->GetNumberOfGridStep(); gid++)
00551         {
00552         vtkMedGrid* grid = mesh->GetGridStep(gid);
00553         // read point family data and create EntityArrays
00554 
00555         for(int eid=0; eid < grid->GetNumberOfEntityArray(); eid++)
00556          {
00557           vtkMedEntityArray* array = grid->GetEntityArray(eid);
00558 
00559           // Next continue is to avoid to create filters for the
00560           // points, at the moment we charge the points in all nodes
00561           if (array->GetEntity().GeometryType == MED_POINT1) // !MED_NODE
00562             continue;
00563 
00564           med_int nbofconstituentpervalue = vtkMedUtilities::GetNumberOfNodes(
00565                                             array->GetEntity().GeometryType);
00566           if (nbofconstituentpervalue == -1)
00567             vtkErrorMacro("Still not implemented for MED_POLYGON and MED_POLYHEDRON"); // à gerer
00568 
00569           // Calculating block sizes
00570           int nEntity = array->GetNumberOfEntity();
00571           int block_size = ( nEntity / this->Internal->NumberOfPieces );
00572           med_size    start  = block_size * this->Internal->CurrentPieceNumber + 1;
00573           med_size    stride = block_size;
00574           med_size    count  = 1;
00575           med_size    blocksize = block_size;
00576           med_size    lastblocksize = (nEntity % this->Internal->NumberOfPieces);
00577           if ((this->Internal->NumberOfPieces ==
00578               this->Internal->CurrentPieceNumber+1) && (lastblocksize != 0))
00579             {
00580             blocksize += lastblocksize;
00581             stride    += lastblocksize;
00582             }
00583           lastblocksize = 0;
00584 
00585           vtkMedFilter *filter = vtkMedFilter::New();
00586           filter->SetFilterSizes( start, stride, count, blocksize, lastblocksize );
00587           array->SetFilter(filter);
00588          }//entity array
00589         }// grid step
00590       }//mesh
00591     }//mesh file
00592 
00593   // SECOND: Filters for the Fields
00594 
00595   std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
00596       fieldfileit;
00597   // link the FieldOnProfile with the profiles
00598   fieldfileit = this->Internal->MedFiles.begin();
00599   while(fieldfileit != this->Internal->MedFiles.end())
00600     {
00601     vtkMedFile* fieldfile = fieldfileit->second;
00602     fieldfileit++;
00603     med_idt pFileID = fieldfile->GetMedDriver()->GetParallelFileId();
00604 
00605     for(int fid=0; fid<fieldfile->GetNumberOfField(); fid++)
00606       {
00607       vtkMedField* field = fieldfile->GetField(fid);
00608 
00609       if (field->GetFieldType() == vtkMedField::CellField)
00610       {
00611       for(int sid = 0; sid< field->GetNumberOfFieldStep(); sid++)
00612         {
00613         vtkMedFieldStep* step = field->GetFieldStep(sid);
00614 
00615         for(int foeid = 0; foeid < step->GetNumberOfFieldOverEntity(); foeid++)
00616         // TODO : seul le premier pas de temps est dispo au debut
00617           {
00618           vtkMedFieldOverEntity* fieldOverEntity = step->GetFieldOverEntity(foeid);
00619 
00620           for(int fopid = 0; fopid < fieldOverEntity->GetNumberOfFieldOnProfile(); fopid++)
00621             {
00622             vtkMedFieldOnProfile* fop = fieldOverEntity->GetFieldOnProfile(fopid);
00623             // Here implement the filters as before:
00624             // 1- Modify vtkMedFieldOnProfile to contain a filter
00625             // 2- Create the filters here only if they are on CELLs (use GetFieldType)
00626             med_int nbofconstituentpervalue = field->GetNumberOfComponent();
00627 
00628             int nVectors = fop->GetNumberOfValues();
00629 
00630             int block_size = ( nVectors / this->Internal->NumberOfPieces );
00631             int    start  = block_size * this->Internal->CurrentPieceNumber + 1;
00632             int    stride = block_size;
00633             int    count  = 1;
00634             int    blocksize = block_size;
00635             int    lastblocksize = (nVectors % this->Internal->NumberOfPieces);
00636             if ((this->Internal->NumberOfPieces ==
00637                  this->Internal->CurrentPieceNumber+1) && (lastblocksize != 0))
00638               {
00639               blocksize += lastblocksize;
00640               stride    += lastblocksize;
00641               }
00642             lastblocksize = 0;
00643 
00644             vtkMedFilter *filter = vtkMedFilter::New();
00645             filter->SetFilterSizes( start, stride, count, blocksize, lastblocksize );
00646             fop->SetFilter(filter);
00647             }
00648           }
00649         }
00650       } // end IF
00651       }
00652     }
00653 
00654 }
00655 
00656 void  vtkMedReader::LinkMedInfo()
00657 {
00658   std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
00659       fieldfileit;
00660   // link the FieldOnProfile with the profiles
00661   fieldfileit = this->Internal->MedFiles.begin();
00662   while(fieldfileit != this->Internal->MedFiles.end())
00663     {
00664     vtkMedFile* fieldfile = fieldfileit->second;
00665     fieldfileit++;
00666 
00667     for(int fid=0; fid<fieldfile->GetNumberOfField(); fid++)
00668       {
00669       vtkMedField* field = fieldfile->GetField(fid);
00670 
00671       for(int sid = 0; sid< field->GetNumberOfFieldStep(); sid++)
00672         {
00673         vtkMedFieldStep* step = field->GetFieldStep(sid);
00674 
00675         for(int foeid = 0; foeid < step->GetNumberOfFieldOverEntity(); foeid++)
00676           {
00677           vtkMedFieldOverEntity* fieldOverEntity = step->GetFieldOverEntity(foeid);
00678 
00679           for(int fopid = 0; fopid < fieldOverEntity->GetNumberOfFieldOnProfile(); fopid++)
00680             {
00681             vtkMedFieldOnProfile* fop = fieldOverEntity->GetFieldOnProfile(fopid);
00682 
00683             std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
00684                 profilefileit = this->Internal->MedFiles.begin();
00685             while(profilefileit != this->Internal->MedFiles.end() && fop->GetProfile() == NULL)
00686               {
00687               vtkMedFile* profilefile = profilefileit->second;
00688               profilefileit++;
00689 
00690               for(int pid = 0; pid < profilefile->GetNumberOfProfile(); pid++)
00691                 {
00692                 vtkMedProfile *profile = profilefile->GetProfile(pid);
00693                 if(strcmp(profile->GetName(), fop->GetProfileName()) == 0)
00694                   {
00695                   fop->SetProfile(profile);
00696                   break;
00697                   }
00698                 }
00699               }
00700             }
00701           }
00702         }
00703       }
00704     }
00705 
00706   // first, add a familyOnEntityOnProfile to all FamilyOnEntity with a NULL
00707   // profile. This is used if no field is mapped to this support.
00708   std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
00709       meshfit = this->Internal->MedFiles.begin();
00710   while(meshfit != this->Internal->MedFiles.end())
00711     {
00712     vtkMedFile* meshfile = meshfit->second;
00713     meshfit++;
00714 
00715     for(int mid=0; mid<meshfile->GetNumberOfMesh(); mid++)
00716       {
00717       vtkMedMesh* mesh = meshfile->GetMesh(mid);
00718 
00719       for(int gid=0; gid<mesh->GetNumberOfGridStep(); gid++)
00720         {
00721         vtkMedGrid* grid = mesh->GetGridStep(gid);
00722         // read point family data and create EntityArrays
00723 
00724         for(int eid=0; eid < grid->GetNumberOfEntityArray(); eid++)
00725           {
00726           vtkMedEntityArray* array = grid->GetEntityArray(eid);
00727 
00728           for(int fid=0; fid < array->GetNumberOfFamilyOnEntity(); fid++)
00729             {
00730             vtkMedFamilyOnEntity* foe = array->GetFamilyOnEntity(fid);
00731             if(foe->GetFamilyOnEntityOnProfile((vtkMedProfile*)NULL) == NULL)
00732               {
00733               vtkMedFamilyOnEntityOnProfile* foep =
00734                   vtkMedFamilyOnEntityOnProfile::New();
00735               foep->SetFamilyOnEntity(foe);
00736               foep->SetProfile(NULL);
00737               foe->AddFamilyOnEntityOnProfile(foep);
00738               foep->Delete();
00739               }
00740             }//family on entity
00741           }//entity array
00742         }// grid step
00743       }//mesh
00744     }//mesh file
00745 
00746   fieldfileit = this->Internal->MedFiles.begin();
00747   while(fieldfileit != this->Internal->MedFiles.end())
00748     {
00749     vtkMedFile* fieldfile = fieldfileit->second;
00750     fieldfileit++;
00751 
00752     for(int fieldid=0; fieldid < fieldfile->GetNumberOfField(); fieldid++)
00753       {
00754       vtkMedField* field = fieldfile->GetField(fieldid);
00755 
00756       for(int fstepid=0; fstepid < field->GetNumberOfFieldStep(); fstepid++)
00757         {
00758         vtkMedFieldStep* step = field->GetFieldStep(fstepid);
00759 
00760         vtkMedComputeStep meshcs = step->GetMeshComputeStep();
00761 
00762         for(int foeid=0; foeid<step->GetNumberOfFieldOverEntity() ;foeid++)
00763           {
00764           vtkMedFieldOverEntity* fieldOverEntity = step->GetFieldOverEntity(foeid);
00765           const vtkMedEntity& fieldentity = fieldOverEntity->GetEntity();
00766 
00767           for (int fopid = 0;
00768                fopid < fieldOverEntity->GetNumberOfFieldOnProfile(); fopid++)
00769             {
00770             vtkMedFieldOnProfile* fop =
00771                 fieldOverEntity->GetFieldOnProfile(fopid);
00772 
00773             std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
00774                 meshfileit = this->Internal->MedFiles.begin();
00775             while(meshfileit != this->Internal->MedFiles.end())
00776               {
00777               vtkMedFile* meshfile = meshfileit->second;
00778               meshfileit++;
00779 
00780               if(field->GetLocal() == 1 && (meshfile != fieldfile))
00781                 continue;
00782 
00783               for(int mid=0; mid<meshfile->GetNumberOfMesh(); mid++)
00784                 {
00785                 vtkMedMesh* mesh = meshfile->GetMesh(mid);
00786 
00787                 // the field must be on this mesh.
00788                 if(strcmp(mesh->GetName(),
00789                           field->GetMeshName()) != 0)
00790                   continue;
00791 
00792                 vtkMedGrid* grid = mesh->GetGridStep(meshcs);
00793                 if(grid == NULL)
00794                   {
00795                   vtkErrorMacro("the field " << field->GetName()
00796                                 << " at step iteration:"
00797                                 << step->GetComputeStep().IterationIt
00798                                 << " and time "
00799                                 << step->GetComputeStep().TimeIt
00800                                 << " uses mesh at step "
00801                                 << meshcs.TimeIt << " " << meshcs.IterationIt
00802                                 << "which does not exists!");
00803                   continue;
00804                   }
00805 
00806                 for(int eid=0; eid < grid->GetNumberOfEntityArray(); eid++)
00807                   {
00808                   vtkMedEntityArray* array = grid->GetEntityArray(eid);
00809 
00810                   // if the support is on points,
00811                   // the field must also be on points
00812                   if(array->GetEntity().EntityType == MED_NODE &&
00813                      fieldentity.EntityType != MED_NODE)
00814                     continue;
00815 
00816                   if(array->GetEntity().EntityType != MED_NODE &&
00817                      fieldentity.EntityType == MED_NODE)
00818                     continue;
00819 
00820                   // for fields not on points, the geometry type
00821                   // of the support must match
00822                   if(array->GetEntity().EntityType != MED_NODE &&
00823                      array->GetEntity().GeometryType != fieldentity.GeometryType)
00824                     continue;
00825 
00826                   for(int fid = 0; fid < array->GetNumberOfFamilyOnEntity(); fid++)
00827                     {
00828                     vtkMedFamilyOnEntity* foe = array->GetFamilyOnEntity(fid);
00829                     if(foe->GetFamilyOnEntityOnProfile(fop->GetProfile()) == NULL)
00830                       {
00831                       vtkMedFamilyOnEntityOnProfile* foep =
00832                           vtkMedFamilyOnEntityOnProfile::New();
00833                       foep->SetProfile(fop->GetProfile());
00834                       foep->SetFamilyOnEntity(foe);
00835                       foe->AddFamilyOnEntityOnProfile(foep);
00836                       foep->Delete();
00837                       }
00838                     // also add the family on entity with no profile.
00839                     if(foe->GetFamilyOnEntityOnProfile((vtkMedProfile*)NULL) == NULL)
00840                       {
00841                       vtkMedFamilyOnEntityOnProfile* foep =
00842                           vtkMedFamilyOnEntityOnProfile::New();
00843                       foep->SetProfile(NULL);
00844                       foep->SetFamilyOnEntity(foe);
00845                       foe->AddFamilyOnEntityOnProfile(foep);
00846                       foep->Delete();
00847                       }
00848                     }//familyOnEntity
00849                   }//entityArray
00850                 }//mesh
00851               }//mesh file
00852             }//field on profile
00853           }//fieldOverEntity
00854         }//fiedstep
00855       }// fieldid
00856     }//fieldfileit
00857 
00858   // Now, link localizations and interpolations
00859   fieldfileit = this->Internal->MedFiles.begin();
00860   while(fieldfileit != this->Internal->MedFiles.end())
00861     {
00862     vtkMedFile* fieldfile = fieldfileit->second;
00863     fieldfileit++;
00864 
00865     for(int locid = 0; locid < fieldfile->GetNumberOfLocalization(); locid++)
00866       {
00867       vtkMedLocalization* loc = fieldfile->GetLocalization(locid);
00868 
00869       for(int fid = 0; fid < fieldfile->GetNumberOfField() &&
00870                     loc->GetInterpolation() == NULL; fid++)
00871         {
00872         vtkMedField* field = fieldfile->GetField(fid);
00873         for(int interpid = 0; interpid < field->GetNumberOfInterpolation();
00874         interpid++)
00875           {
00876           vtkMedInterpolation* interp = field->GetInterpolation(interpid);
00877           if(strcmp(loc->GetInterpolationName(),
00878                     interp->GetName()) == 0)
00879             {
00880             loc->SetInterpolation(interp);
00881             break;
00882             }
00883           }
00884         }
00885       }
00886     }
00887 
00888   // now that the interpolation is set, build the  shape functions.
00889   fieldfileit = this->Internal->MedFiles.begin();
00890   while(fieldfileit != this->Internal->MedFiles.end())
00891     {
00892     vtkMedFile* fieldfile = fieldfileit->second;
00893     fieldfileit++;
00894 
00895     for(int locid = 0; locid < fieldfile->GetNumberOfLocalization(); locid++)
00896       {
00897       vtkMedLocalization* loc = fieldfile->GetLocalization(locid);
00898       loc->BuildShapeFunction();
00899       }
00900     }
00901 
00902   // set the supportmesh pointer in the structural element
00903   std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
00904       fileit = this->Internal->MedFiles.begin();
00905   while(fileit != this->Internal->MedFiles.end())
00906     {
00907     vtkMedFile* file = fileit->second;
00908     fileit++;
00909 
00910     for(int structelemit = 0;
00911         structelemit < file->GetNumberOfStructElement();
00912         structelemit++)
00913       {
00914       vtkMedStructElement* structElem =
00915           file->GetStructElement(structelemit);
00916 
00917       for(int supportmeshit = 0;
00918           supportmeshit < file->GetNumberOfSupportMesh();
00919           supportmeshit++)
00920         {
00921         vtkMedMesh* supportMesh =
00922             file->GetSupportMesh(supportmeshit);
00923 
00924         if(strcmp(supportMesh->GetName(), structElem->GetName()) == 0 )
00925           {
00926           structElem->SetSupportMesh(supportMesh);
00927           break;
00928           }
00929         }
00930       }
00931     }
00932 
00933   // set the pointer to the profile used by the constant attributes
00934   fileit = this->Internal->MedFiles.begin();
00935   while(fileit != this->Internal->MedFiles.end())
00936   {
00937   vtkMedFile* file = fileit->second;
00938   fileit++;
00939 
00940   for(int structelemit = 0;
00941       structelemit < file->GetNumberOfStructElement();
00942       structelemit++)
00943     {
00944     vtkMedStructElement* structElem =
00945         file->GetStructElement(structelemit);
00946 
00947     for(int cstattit = 0; cstattit < structElem->GetNumberOfConstantAttribute(); cstattit++)
00948       {
00949       vtkMedConstantAttribute* cstatt = structElem->GetConstantAttribute(cstattit);
00950 
00951       for(int profit = 0;
00952           profit < file->GetNumberOfProfile();
00953           profit++)
00954         {
00955         vtkMedProfile* profile =
00956             file->GetProfile(profit);
00957 
00958         if(strcmp(profile->GetName(), cstatt->GetProfileName()) == 0 )
00959           {
00960           cstatt->SetProfile(profile);
00961           break;
00962           }
00963         }
00964       }
00965     }
00966   }
00967 
00968   meshfit = this->Internal->MedFiles.begin();
00969   while(meshfit != this->Internal->MedFiles.end())
00970   {
00971   vtkMedFile* meshfile = meshfit->second;
00972   meshfit++;
00973 
00974   for(int mid=0; mid<meshfile->GetNumberOfMesh(); mid++)
00975     {
00976     vtkMedMesh* mesh = meshfile->GetMesh(mid);
00977 
00978     for(int gid=0; gid<mesh->GetNumberOfGridStep(); gid++)
00979       {
00980       vtkMedGrid* grid = mesh->GetGridStep(gid);
00981       // read point family data and create EntityArrays
00982 
00983       for(int eid=0; eid < grid->GetNumberOfEntityArray(); eid++)
00984         {
00985         vtkMedEntityArray* array = grid->GetEntityArray(eid);
00986         if(array->GetEntity().EntityType != MED_STRUCT_ELEMENT)
00987           continue;
00988 
00989         for(int structelemit = 0; structelemit < meshfile->GetNumberOfStructElement(); structelemit++)
00990           {
00991           vtkMedStructElement* structelem = meshfile->GetStructElement(structelemit);
00992           if(structelem->GetGeometryType() == array->GetEntity().GeometryType)
00993             {
00994             array->SetStructElement(structelem);
00995             break;
00996             }
00997           }
00998         }
00999       }
01000     }
01001   }
01002 }
01003 
01004 int vtkMedReader::GetFrequencyArrayStatus(const char* name)
01005 {
01006   return this->Frequencies->GetKeyStatus(name);
01007 }
01008 
01009 void vtkMedReader::SetFrequencyArrayStatus(const char* name, int status)
01010 {
01011   if(this->Frequencies->GetKeyStatus(name) == status)
01012     {
01013     return;
01014     }
01015 
01016   this->Frequencies->SetKeyStatus(name, status);
01017 
01018   this->Modified();
01019 }
01020 
01021 int vtkMedReader::GetNumberOfFrequencyArrays()
01022 {
01023   return this->Frequencies->GetNumberOfKey();
01024 }
01025 
01026 const char* vtkMedReader::GetFrequencyArrayName(int index)
01027 {
01028   return this->Frequencies->GetKey(index);
01029 }
01030 
01031 struct compTimes
01032 {
01033   bool operator()(pair<double, int> i, pair<double, int> j)
01034   {
01035     if(i.first!=j.first)
01036       return (i.first<j.first);
01037     return i.second<j.second;
01038   }
01039 };
01040 
01041 vtkDoubleArray* vtkMedReader::GetAvailableTimes()
01042 {
01043   this->AvailableTimes->Initialize();
01044   this->AvailableTimes->SetNumberOfComponents(1);
01045 
01046   std::set<std::string> newFrequencies;
01047 
01048   int tid = 0;
01049   std::map<med_float, std::set<med_int> >::iterator it =
01050       this->Internal->GlobalComputeStep.begin();
01051   while(it != this->Internal->GlobalComputeStep.end())
01052     {
01053     double time = it->first;
01054     this->AvailableTimes->InsertNextValue(time);
01055     string name = vtkMedUtilities::GetModeKey(tid, time, this->Internal->GlobalComputeStep.size()-1);
01056     this->Frequencies->AddKey(name.c_str());
01057     newFrequencies.insert(name);
01058     tid++;
01059     it++;
01060     }
01061 
01062   // now check if old frequencies have been removed
01063   for(int f = 0; f < this->Frequencies->GetNumberOfKey(); f++)
01064     {
01065     const char* name = this->Frequencies->GetKey(f);
01066     if(newFrequencies.find(name) == newFrequencies.end())
01067       {
01068       this->Frequencies->RemoveKeyByIndex(f);
01069       f--;
01070       }
01071     }
01072 
01073   return this->AvailableTimes;
01074 }
01075 
01076 void vtkMedReader::ChooseRealAnimationMode()
01077 {
01078   if(this->AnimationMode!=Default)
01079     {
01080     this->Internal->RealAnimationMode=this->AnimationMode;
01081     return;
01082     }
01083 
01084   // if there is exactly one physical time and more than one iteration
01085   // set the animation mode to iteration, else default to physical time.
01086   if (this->Internal->GlobalComputeStep.size() == 1 &&
01087       this->Internal->GlobalComputeStep[0].size() > 1)
01088     {
01089     this->Internal->RealAnimationMode=Iteration;
01090     return;
01091     }
01092 
01093   this->Internal->RealAnimationMode=PhysicalTime;
01094 }
01095 
01096 int vtkMedReader::GetEntityStatus(const vtkMedEntity& entity)
01097 {
01098   if (entity.EntityType==MED_NODE)
01099     return 1;
01100   if(entity.EntityType == MED_DESCENDING_FACE
01101      || entity.EntityType == MED_DESCENDING_EDGE)
01102     return 0;
01103 
01104   return this->Entities->GetKeyStatus(vtkMedUtilities::EntityKey(entity).c_str());
01105 }
01106 
01107 int vtkMedReader::GetFamilyStatus(vtkMedMesh* mesh, vtkMedFamily* family)
01108 {
01109   if(!mesh||!family)
01110     return 0;
01111 
01112   if(this->Internal->GroupSelectionMTime > this->Internal->FamilySelectionMTime)
01113     {
01114     this->SelectFamiliesFromGroups();
01115     }
01116 
01117   int status =  this->Internal->Families->GetKeyStatus(vtkMedUtilities::FamilyKey(
01118       mesh->GetName(), family->GetPointOrCell(),
01119       family->GetName()).c_str());
01120 
01121   return status;
01122 }
01123 
01124 int vtkMedReader::IsMeshSelected(vtkMedMesh* mesh)
01125 {
01126   for(int fam=0; fam<mesh->GetNumberOfPointFamily(); fam++)
01127     {
01128     if(this->GetFamilyStatus(mesh, mesh->GetPointFamily(fam))!=0)
01129       return 1;
01130     }
01131 
01132   for(int fam=0; fam<mesh->GetNumberOfCellFamily(); fam++)
01133     {
01134     if(this->GetFamilyStatus(mesh, mesh->GetCellFamily(fam))!=0)
01135       return 1;
01136     }
01137   return 0;
01138 }
01139 
01140 void vtkMedReader::GatherComputeSteps()
01141 {
01142   this->Internal->GlobalComputeStep.clear();
01143 
01144   std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
01145       fieldfileit = this->Internal->MedFiles.begin();
01146   while(fieldfileit != this->Internal->MedFiles.end())
01147     {
01148     vtkMedFile* file = fieldfileit->second;
01149     fieldfileit++;
01150 
01151     // first loop over all fields to gather their compute steps
01152     for(int fieldId=0; fieldId<file->GetNumberOfField(); fieldId++)
01153       {
01154       vtkMedField* field=file->GetField(fieldId);
01155 
01156       for(int stepId=0; stepId<field->GetNumberOfFieldStep(); stepId++)
01157         {
01158         vtkMedFieldStep* step=field->GetFieldStep(stepId);
01159         const vtkMedComputeStep& cs = step->GetComputeStep();
01160         this->Internal->GlobalComputeStep[cs.TimeOrFrequency].insert(cs.IterationIt);
01161         }
01162       }//fields
01163 
01164     // then loop over all meshes to gather their grid steps too.
01165     // for meshes, do not add the MED_UNDEF_DT time
01166     for(int meshId=0; meshId<file->GetNumberOfMesh(); meshId++)
01167       {
01168       vtkMedMesh* mesh=file->GetMesh(meshId);
01169 
01170       for(int stepId=0; stepId<mesh->GetNumberOfGridStep(); stepId++)
01171         {
01172         vtkMedGrid* grid=mesh->GetGridStep(stepId);
01173         const vtkMedComputeStep& cs = grid->GetComputeStep();
01174         if(cs.TimeOrFrequency != MED_UNDEF_DT || cs.TimeIt != MED_NO_DT)
01175           {
01176           this->Internal->GlobalComputeStep[cs.TimeOrFrequency].insert(cs.IterationIt);
01177           }
01178         }
01179       }//mesh
01180     }
01181   if(this->Internal->GlobalComputeStep.size() == 0)
01182     {
01183     this->Internal->GlobalComputeStep[MED_UNDEF_DT].insert(MED_NO_IT);
01184     }
01185 }
01186 
01187 int vtkMedReader::IsFieldSelected(vtkMedField* field)
01188 {
01189   return this->IsPointFieldSelected(field)||this->IsCellFieldSelected(field)
01190       ||this->IsQuadratureFieldSelected(field) || this->IsElnoFieldSelected(field);
01191 }
01192 
01193 int vtkMedReader::IsPointFieldSelected(vtkMedField* field)
01194 {
01195   return field->GetFieldType()==vtkMedField::PointField
01196       &&this->GetPointFieldArrayStatus(vtkMedUtilities::SimplifyName(
01197           field->GetName()).c_str());
01198 }
01199 
01200 int vtkMedReader::IsCellFieldSelected(vtkMedField* field)
01201 {
01202   return field->GetFieldType()==vtkMedField::CellField
01203       &&this->GetCellFieldArrayStatus(vtkMedUtilities::SimplifyName(
01204           field->GetName()).c_str());
01205 }
01206 
01207 vtkMedProfile* vtkMedReader::GetProfile(const char* pname)
01208 {
01209   std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
01210       fileit = this->Internal->MedFiles.begin();
01211   while(fileit != this->Internal->MedFiles.end())
01212     {
01213     vtkMedFile* file = fileit->second;
01214     fileit++;
01215     vtkMedProfile* profile = file->GetProfile(pname);
01216     if(profile != NULL)
01217       return profile;
01218     }
01219   return NULL;
01220 }
01221 
01222 vtkMedLocalization* vtkMedReader::GetLocalization(const char* lname)
01223 {
01224   std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
01225       fileit = this->Internal->MedFiles.begin();
01226   while(fileit != this->Internal->MedFiles.end())
01227     {
01228     vtkMedFile* file = fileit->second;
01229     fileit++;
01230     vtkMedLocalization* loc = file->GetLocalization(lname);
01231     if(loc != NULL)
01232       return loc;
01233     }
01234   return NULL;
01235 
01236 }
01237 
01238 int vtkMedReader::GetLocalizationKey(vtkMedFieldOnProfile* fop)
01239 {
01240   vtkMedLocalization* def=this->GetLocalization(fop->GetLocalizationName());
01241 
01242   // This is not a quadrature field with explicit definition.
01243   // There are two possible cases : either the intergration point is
01244   // at the center of the cell
01245   //1 quadrature point at the cell center
01246   int nloc = 0;
01247   std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
01248       fileit = this->Internal->MedFiles.begin();
01249   while(fileit != this->Internal->MedFiles.end())
01250     {
01251     vtkMedFile* file = fileit->second;
01252     fileit++;
01253 
01254     if(def && def->GetParentFile() == file)
01255       return nloc + def->GetMedIterator() - 1;
01256 
01257     nloc += file->GetNumberOfLocalization();
01258     }
01259 
01260   // center of a cell
01261   if(fop->GetNumberOfIntegrationPoint()==1)
01262     return nloc + 1 + fop->GetParentFieldOverEntity()->GetEntity().GeometryType;
01263 
01264   // or it is an elno field (field stored on nodes of the cells,
01265   // but with discontinuities at the vertices)
01266   return -fop->GetParentFieldOverEntity()->GetEntity().GeometryType;//ELNO
01267 }
01268 
01269 int vtkMedReader::IsQuadratureFieldSelected(vtkMedField* field)
01270 {
01271   return field->GetFieldType()==vtkMedField::QuadratureField
01272       &&this->GetQuadratureFieldArrayStatus(vtkMedUtilities::SimplifyName(
01273           field->GetName()).c_str());
01274 }
01275 
01276 int vtkMedReader::IsElnoFieldSelected(vtkMedField* field)
01277 {
01278   return field->GetFieldType()==vtkMedField::ElnoField
01279       &&this->GetElnoFieldArrayStatus(vtkMedUtilities::SimplifyName(
01280           field->GetName()).c_str());
01281 }
01282 
01283 // Description:
01284 // Give the animation steps to the pipeline
01285 void vtkMedReader::AdvertiseTime(vtkInformation* info)
01286 {
01287   this->ChooseRealAnimationMode();
01288 
01289   if(this->Internal->RealAnimationMode==PhysicalTime)
01290     {
01291     // I advertise the union of all times available
01292     // in all selected fields and meshes
01293     set<double> timeset;
01294 
01295     std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
01296         fieldfileit = this->Internal->MedFiles.begin();
01297     while(fieldfileit != this->Internal->MedFiles.end())
01298       {
01299       vtkMedFile* file = fieldfileit->second;
01300       fieldfileit++;
01301 
01302       // first loop over all fields to gather their compute steps
01303       for(int fieldId=0; fieldId<file->GetNumberOfField(); fieldId++)
01304         {
01305         vtkMedField* field=file->GetField(fieldId);
01306 
01307         if(!this->IsFieldSelected(field))
01308           continue;
01309 
01310         field->GatherFieldTimes(timeset);
01311         }//fields
01312 
01313       // then loop over all meshes to gather their grid steps too.
01314       for(int meshId=0; meshId<file->GetNumberOfMesh(); meshId++)
01315         {
01316         vtkMedMesh* mesh=file->GetMesh(meshId);
01317 
01318         if(!this->IsMeshSelected(mesh))
01319           continue;
01320 
01321         mesh->GatherGridTimes(timeset);
01322         }//meshes
01323       }
01324 
01325     if(timeset.size() > 0)
01326       {
01327       // remove MED_UNDEF_DT if there are other time step
01328       if(timeset.size() > 1)
01329         timeset.erase(MED_UNDEF_DT);
01330 
01331       vector<double> times;
01332       set<double>::iterator it = timeset.begin();
01333       while(it != timeset.end())
01334         {
01335         times.push_back(*it);
01336         it++;
01337         }
01338       sort(times.begin(), times.end());
01339 
01340       info->Set(vtkStreamingDemandDrivenPipeline::TIME_STEPS(), &times[0],
01341           times.size());
01342       double timeRange[2];
01343       timeRange[0]=times[0];
01344       timeRange[1]=times[times.size()-1];
01345       info->Set(vtkStreamingDemandDrivenPipeline::TIME_RANGE(), &timeRange[0],
01346           2);
01347       }
01348     else
01349       {
01350       info->Remove(vtkStreamingDemandDrivenPipeline::TIME_STEPS());
01351       info->Remove(vtkStreamingDemandDrivenPipeline::TIME_RANGE());
01352       }
01353     }
01354   else if(this->Internal->RealAnimationMode==Iteration)
01355     {
01356     // I advertise the union of all iterations available at the given
01357     // Time for all selected fields.
01358     set<med_int> iterationsets;
01359     med_float time = MED_UNDEF_DT;
01360     if(this->TimeIndexForIterations >= 0 &&
01361        this->TimeIndexForIterations <
01362        this->AvailableTimes->GetNumberOfTuples())
01363       {
01364       time = this->AvailableTimes->
01365                      GetValue((vtkIdType)this->TimeIndexForIterations);
01366       }
01367 
01368     std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
01369         fieldfileit = this->Internal->MedFiles.begin();
01370     while(fieldfileit != this->Internal->MedFiles.end())
01371       {
01372       vtkMedFile* file = fieldfileit->second;
01373       fieldfileit++;
01374 
01375       for(int fieldId=0; fieldId<file->GetNumberOfField(); fieldId++)
01376         {
01377         vtkMedField* field=file->GetField(fieldId);
01378         if(!this->IsFieldSelected(field))
01379           continue;
01380 
01381         field->GatherFieldIterations(time, iterationsets);
01382         }
01383       // then loop over all meshes to gather their grid steps too.
01384       for(int meshId=0; meshId<file->GetNumberOfMesh(); meshId++)
01385         {
01386         vtkMedMesh* mesh=file->GetMesh(meshId);
01387 
01388         if(!this->IsMeshSelected(mesh))
01389           continue;
01390 
01391         mesh->GatherGridIterations(time, iterationsets);
01392         }//meshes
01393       }
01394 
01395     if(iterationsets.size()>0)
01396       {
01397       // remove MED_NO_IT if there are other available iterations.
01398       if(iterationsets.size()>1)
01399         iterationsets.erase(MED_NO_IT);
01400 
01401       vector<double> iterations;
01402       set<med_int>::iterator it=iterationsets.begin();
01403       while(it!=iterationsets.end())
01404         {
01405         iterations.push_back((double)*it);
01406         it++;
01407         }
01408       sort(iterations.begin(), iterations.end());
01409       info->Set(vtkStreamingDemandDrivenPipeline::TIME_STEPS(), &iterations[0],
01410           iterations.size());
01411       double timeRange[2];
01412       timeRange[0]=iterations[0];
01413       timeRange[1]=iterations[iterations.size()-1];
01414       info->Set(vtkStreamingDemandDrivenPipeline::TIME_RANGE(), &timeRange[0],
01415           2);
01416       }
01417     else
01418       {
01419       info->Remove(vtkStreamingDemandDrivenPipeline::TIME_STEPS());
01420       info->Remove(vtkStreamingDemandDrivenPipeline::TIME_RANGE());
01421       }
01422     }
01423   else
01424     {
01425     info->Remove(vtkStreamingDemandDrivenPipeline::TIME_STEPS());
01426     info->Remove(vtkStreamingDemandDrivenPipeline::TIME_RANGE());
01427     }
01428 }
01429 
01430 vtkIdType vtkMedReader::GetFrequencyIndex(double freq)
01431 {
01432   return this->AvailableTimes->LookupValue(freq);
01433 }
01434 
01435 int vtkMedReader::RequestDataObject(vtkInformation* request,
01436     vtkInformationVector** vtkNotUsed(inputVector), vtkInformationVector* outputVector)
01437 {
01438   vtkInformation *info = outputVector->GetInformationObject(0);
01439   if (vtkMultiBlockDataSet::SafeDownCast(
01440       info->Get(vtkDataObject::DATA_OBJECT())))
01441     {
01442     // The output is already created
01443     return 1;
01444     }
01445   else
01446     {
01447     vtkMultiBlockDataSet* output=vtkMultiBlockDataSet::New();
01448     this->GetExecutive()->SetOutputData(0, output);
01449     output->Delete();
01450     this->GetOutputPortInformation(0)->Set(vtkDataObject::DATA_EXTENT_TYPE(),
01451         output->GetExtentType());
01452     }
01453   return 1;
01454 }
01455 
01456 void vtkMedReader::ClearSelections()
01457 {
01458   this->PointFields->Initialize();
01459   this->CellFields->Initialize();
01460   this->QuadratureFields->Initialize();
01461   this->ElnoFields->Initialize();
01462 
01463   std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
01464       fieldfileit = this->Internal->MedFiles.begin();
01465   while(fieldfileit != this->Internal->MedFiles.end())
01466     {
01467     vtkMedFile* file = fieldfileit->second;
01468     fieldfileit++;
01469 
01470     for(int index=0; index < file->GetNumberOfField(); index++)
01471       {
01472       vtkMedField* field = file->GetField(index);
01473       switch(field->GetFieldType())
01474         {
01475         case vtkMedField::PointField :
01476         this->PointFields->AddKey(vtkMedUtilities::SimplifyName(
01477               field->GetName()).c_str());
01478         break;
01479         case vtkMedField::CellField :
01480         this->CellFields->AddKey(vtkMedUtilities::SimplifyName(
01481               field->GetName()).c_str());
01482         break;
01483         case vtkMedField::QuadratureField :
01484         this->QuadratureFields->AddKey(vtkMedUtilities::SimplifyName(
01485               field->GetName()).c_str());
01486         break;
01487         case vtkMedField::ElnoField :
01488         this->ElnoFields->AddKey(vtkMedUtilities::SimplifyName(
01489               field->GetName()).c_str());
01490         break;
01491         }
01492       }
01493 
01494     this->Internal->Families->Initialize();
01495     this->Groups->Initialize();
01496     for(int meshIndex=0; meshIndex < file->GetNumberOfMesh(); meshIndex++)
01497       {
01498       vtkMedMesh* mesh = file->GetMesh(meshIndex);
01499       for(int famIndex=0; famIndex<mesh->GetNumberOfPointFamily(); famIndex++)
01500         {
01501         vtkMedFamily* fam=mesh->GetPointFamily(famIndex);
01502 
01503         int ng=fam->GetNumberOfGroup();
01504         for(int gindex=0; gindex<ng; gindex++)
01505           {
01506           vtkMedGroup* group=fam->GetGroup(gindex);
01507           string gname=vtkMedUtilities::GroupKey(mesh->GetName(),
01508               fam->GetPointOrCell(), group->GetName());
01509 
01510           this->Groups->AddKey(gname.c_str());
01511           this->Groups->SetKeyStatus(gname.c_str(), 0);
01512           }
01513         }
01514       for(int famIndex=0; famIndex<mesh->GetNumberOfCellFamily(); famIndex++)
01515         {
01516         vtkMedFamily* fam=mesh->GetCellFamily(famIndex);
01517 
01518         int ng=fam->GetNumberOfGroup();
01519         for(int gindex=0; gindex<ng; gindex++)
01520           {
01521           vtkMedGroup* group=fam->GetGroup(gindex);
01522           string gname=vtkMedUtilities::GroupKey(mesh->GetName(),
01523               fam->GetPointOrCell(), group->GetName());
01524 
01525           this->Groups->AddKey(gname.c_str());
01526           this->Groups->SetKeyStatus(gname.c_str(), 1);
01527           }
01528         }
01529       }
01530     this->Internal->GroupSelectionMTime.Modified();
01531 
01532     for(int meshIndex=0; meshIndex< file->GetNumberOfMesh(); meshIndex++)
01533       {
01534       if(file->GetMesh(meshIndex)->GetNumberOfGridStep() == 0)
01535         continue;
01536 
01537       vtkMedGrid* grid=file->GetMesh(meshIndex)->GetGridStep(0);
01538 
01539       for(int entityIndex=0; entityIndex<grid->GetNumberOfEntityArray();
01540         entityIndex++)
01541         {
01542         vtkMedEntityArray* array=grid->GetEntityArray(entityIndex);
01543         string name=vtkMedUtilities::EntityKey(array->GetEntity());
01544         this->Entities->AddKey(name.c_str());
01545         }
01546       }
01547     }
01548   this->Modified();
01549 }
01550 
01551 void vtkMedReader::SelectFamiliesFromGroups()
01552 {
01553   this->Internal->Families->Initialize();
01554 
01555   std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
01556       meshfileit = this->Internal->MedFiles.begin();
01557   while(meshfileit != this->Internal->MedFiles.end())
01558     {
01559     vtkMedFile* file = meshfileit->second;
01560     meshfileit++;
01561 
01562     for(int meshIndex=0; meshIndex < file->GetNumberOfMesh(); meshIndex++)
01563       {
01564       vtkMedMesh* mesh = file->GetMesh(meshIndex);
01565       for(int famIndex=0; famIndex<mesh->GetNumberOfFamily(); famIndex++)
01566         {
01567         vtkMedFamily* fam=mesh->GetFamily(famIndex);
01568         string name=vtkMedUtilities::FamilyKey(mesh->GetName(),
01569             fam->GetPointOrCell(), fam->GetName());
01570 
01571         this->Internal->Families->SetKeyStatus(name.c_str(), 0);
01572 
01573         for(int gindex=0; gindex<fam->GetNumberOfGroup(); gindex++)
01574           {
01575           vtkMedGroup* group=fam->GetGroup(gindex);
01576           string gname=vtkMedUtilities::GroupKey(mesh->GetName(),
01577               fam->GetPointOrCell(), group->GetName());
01578           int state=this->Groups->GetKeyStatus(gname.c_str());
01579 
01580           if(state)
01581             {
01582             this->SetFamilyStatus(name.c_str(), 1);
01583             }
01584           }
01585         }
01586       }
01587   }
01588 
01589   this->Internal->FamilySelectionMTime.Modified();
01590 }
01591 
01592 int vtkMedReader::GetNumberOfPointFieldArrays()
01593 {
01594   return this->PointFields->GetNumberOfKey();
01595 }
01596 
01597 const char*
01598 vtkMedReader::GetPointFieldArrayName(int index)
01599 {
01600   return this->PointFields->GetKey(index);
01601 }
01602 
01603 int vtkMedReader::GetPointFieldArrayStatus(const char* name)
01604 {
01605   return this->PointFields->GetKeyStatus(name);
01606 }
01607 
01608 void vtkMedReader::SetPointFieldArrayStatus(const char* name, int status)
01609 {
01610   if(this->PointFields->KeyExists(name)&&this->PointFields->GetKeyStatus(
01611       name)==status)
01612     return;
01613 
01614   this->PointFields->SetKeyStatus(name, status);
01615 
01616   this->Modified();
01617 }
01618 
01619 int vtkMedReader::GetNumberOfCellFieldArrays()
01620 {
01621   return this->CellFields->GetNumberOfKey();
01622 }
01623 
01624 const char*
01625 vtkMedReader::GetCellFieldArrayName(int index)
01626 {
01627   return this->CellFields->GetKey(index);
01628 }
01629 
01630 int vtkMedReader::GetCellFieldArrayStatus(const char* name)
01631 {
01632   return this->CellFields->GetKeyStatus(name);
01633 }
01634 
01635 void vtkMedReader::SetCellFieldArrayStatus(const char* name, int status)
01636 {
01637   if(this->CellFields->KeyExists(name)&&this->CellFields->GetKeyStatus(
01638       name)==status)
01639     return;
01640 
01641   this->CellFields->SetKeyStatus(name, status);
01642 
01643   this->Modified();
01644 }
01645 
01646 int vtkMedReader::GetNumberOfQuadratureFieldArrays()
01647 {
01648   return this->QuadratureFields->GetNumberOfKey();
01649 }
01650 
01651 const char* vtkMedReader::GetQuadratureFieldArrayName(int index)
01652 {
01653   return this->QuadratureFields->GetKey(index);
01654 }
01655 
01656 int vtkMedReader::GetQuadratureFieldArrayStatus(const char* name)
01657 {
01658   return this->QuadratureFields->GetKeyStatus(name);
01659 }
01660 
01661 void vtkMedReader::SetQuadratureFieldArrayStatus(const char* name, int status)
01662 {
01663   if(this->QuadratureFields->KeyExists(name)
01664       &&this->QuadratureFields->GetKeyStatus(name)==status)
01665     return;
01666 
01667   this->QuadratureFields->SetKeyStatus(name, status);
01668 
01669   this->Modified();
01670 }
01671 
01672 int vtkMedReader::GetNumberOfElnoFieldArrays()
01673 {
01674   return this->ElnoFields->GetNumberOfKey();
01675 }
01676 
01677 const char* vtkMedReader::GetElnoFieldArrayName(int index)
01678 {
01679   return this->ElnoFields->GetKey(index);
01680 }
01681 
01682 int vtkMedReader::GetElnoFieldArrayStatus(const char* name)
01683 {
01684   return this->ElnoFields->GetKeyStatus(name);
01685 }
01686 
01687 void vtkMedReader::SetElnoFieldArrayStatus(const char* name, int status)
01688 {
01689   if(this->ElnoFields->KeyExists(name)
01690       &&this->ElnoFields->GetKeyStatus(name)==status)
01691     return;
01692 
01693   this->ElnoFields->SetKeyStatus(name, status);
01694 
01695   this->Modified();
01696 }
01697 
01698 void vtkMedReader::SetEntityStatus(const char* name, int status)
01699 {
01700   if(this->Entities->KeyExists(name)&&this->Entities->GetKeyStatus(name)
01701       ==status)
01702     return;
01703 
01704   this->Entities->SetKeyStatus(name, status);
01705 
01706   this->Modified();
01707 }
01708 
01709 void vtkMedReader::SetFamilyStatus(const char* name, int status)
01710 {
01711   if(this->Internal->Families->KeyExists(name)
01712       &&this->Internal->Families->GetKeyStatus(name)==status)
01713     return;
01714 
01715   this->Internal->Families->SetKeyStatus(name, status);
01716 }
01717 
01718 void vtkMedReader::SetGroupStatus(const char* name, int status)
01719 {
01720 
01721   if(this->Groups->KeyExists(name)&&this->Groups->GetKeyStatus(name)
01722       ==status)
01723     return;
01724 
01725   this->Groups->SetKeyStatus(name, status);
01726 
01727   this->Modified();
01728 
01729   this->Internal->GroupSelectionMTime.Modified();
01730 }
01731 
01732 int vtkMedReader::GetGroupStatus(const char* key)
01733 {
01734   return this->Groups->GetKeyStatus(key);
01735 }
01736 
01737 void vtkMedReader::AddQuadratureSchemeDefinition(vtkInformation* info,
01738     vtkMedLocalization* loc)
01739 {
01740   if(info==NULL||loc==NULL)
01741     return;
01742 
01743   vtkInformationQuadratureSchemeDefinitionVectorKey *key=
01744       vtkQuadratureSchemeDefinition::DICTIONARY();
01745 
01746   vtkQuadratureSchemeDefinition* def=vtkQuadratureSchemeDefinition::New();
01747   int cellType=vtkMedUtilities::GetVTKCellType(loc->GetGeometryType());
01748   def->Initialize(cellType, vtkMedUtilities::GetNumberOfPoint(
01749       loc->GetGeometryType()), loc->GetNumberOfQuadraturePoint(),
01750       (double*)loc->GetShapeFunction()->GetVoidPointer(0),
01751       (double*)loc->GetWeights()->GetVoidPointer(0));
01752   key->Set(info, def, cellType);
01753   def->Delete();
01754 
01755 }
01756 
01757 void vtkMedReader::LoadConnectivity(vtkMedEntityArray* array)
01758 {
01759   vtkMedGrid* grid = array->GetParentGrid();
01760   array->LoadConnectivity();
01761   if (array->GetConnectivity()==MED_NODAL||vtkMedUtilities::GetDimension(
01762       array->GetEntity().GeometryType)<2
01763       || grid->GetParentMesh()->GetMeshType() == MED_STRUCTURED_MESH)
01764     return;
01765 
01766   vtkMedEntity subentity;
01767   subentity.EntityType = vtkMedUtilities::GetSubType(array->GetEntity().EntityType);
01768 
01769   vtkMedUnstructuredGrid* ugrid = vtkMedUnstructuredGrid::SafeDownCast(grid);
01770   if(ugrid == NULL)
01771     return;
01772 
01773   for(int index=0; index<vtkMedUtilities::GetNumberOfSubEntity(
01774       array->GetEntity().GeometryType); index++)
01775     {
01776     subentity.GeometryType = vtkMedUtilities::GetSubGeometry(
01777         array->GetEntity().GeometryType, index);
01778     vtkMedEntityArray* subarray=ugrid->GetEntityArray(subentity);
01779 
01780     if(subarray==NULL)
01781       {
01782       vtkErrorMacro("DESC connectivity used, but sub types do not exist in file.");
01783       return;
01784       }
01785     subarray->LoadConnectivity();
01786     }
01787 }
01788 
01789 vtkDataSet* vtkMedReader::CreateUnstructuredGridForPointSupport(
01790     vtkMedFamilyOnEntityOnProfile* foep)
01791 {
01792   vtkUnstructuredGrid* vtkgrid = vtkUnstructuredGrid::New();
01793   foep->ComputeIntersection(NULL);
01794   vtkMedFamilyOnEntity* foe = foep->GetFamilyOnEntity();
01795   vtkMedGrid* medgrid=foe->GetParentGrid();
01796   vtkMedUnstructuredGrid* medugrid=vtkMedUnstructuredGrid::SafeDownCast(
01797       medgrid);
01798   vtkMedCurvilinearGrid* medcgrid=vtkMedCurvilinearGrid::SafeDownCast(
01799       medgrid);
01800 
01801   medgrid->LoadCoordinates();
01802 
01803   vtkIdType npts=medgrid->GetNumberOfPoints();
01804 
01805   bool shallowCopy= (medugrid != NULL || medcgrid!=NULL);
01806   if(medgrid->GetParentMesh()->GetSpaceDimension()!=3)
01807     {
01808     shallowCopy=false;
01809     }
01810   else
01811     {
01812     shallowCopy = foep->CanShallowCopyPointField(NULL);
01813     }
01814 
01815   vtkDataArray* coords = NULL;
01816 
01817   if(medugrid != NULL)
01818     coords = medugrid->GetCoordinates();
01819   if(medcgrid != NULL)
01820     coords = medcgrid->GetCoordinates();
01821 
01822 
01823   vtkIdType numberOfPoints;
01824   vtkPoints* points=vtkPoints::New(coords->GetDataType());
01825   vtkgrid->SetPoints(points);
01826   points->Delete();
01827 
01828   vtkIdTypeArray* pointGlobalIds=vtkIdTypeArray::New();
01829   pointGlobalIds->SetName("MED_POINT_ID");
01830   pointGlobalIds->SetNumberOfComponents(1);
01831   vtkgrid->GetPointData()->SetGlobalIds(pointGlobalIds);
01832   pointGlobalIds->Delete();
01833 
01834   if (shallowCopy)
01835     {
01836     vtkgrid->GetPoints()->SetData(coords);
01837     numberOfPoints=npts;
01838 
01839     pointGlobalIds->SetNumberOfTuples(numberOfPoints);
01840     vtkIdType* ptr=pointGlobalIds->GetPointer(0);
01841     for(int pid=0; pid<numberOfPoints; pid++)
01842       ptr[pid]=pid+1;
01843     }
01844   if(!shallowCopy)
01845     {
01846     vtkIdType currentIndex=0;
01847 
01848     for(vtkIdType index=0; index<medgrid->GetNumberOfPoints(); index++)
01849       {
01850       if (!foep->KeepPoint(index))
01851         {
01852         continue;
01853         }
01854 
01855       double coord[3]={0.0, 0.0, 0.0};
01856       double * tuple=medgrid->GetCoordTuple(index);
01857       for(int dim=0; dim<medgrid->GetParentMesh()->GetSpaceDimension()&&dim<3; dim++)
01858         {
01859         coord[dim]=tuple[dim];
01860         }
01861       vtkgrid->GetPoints()->InsertPoint(currentIndex, coord);
01862       pointGlobalIds->InsertNextValue(index+1);
01863       currentIndex++;
01864       }
01865     vtkgrid->GetPoints()->Squeeze();
01866     pointGlobalIds->Squeeze();
01867     numberOfPoints=currentIndex;
01868     }
01869 
01870   // now create the VTK_VERTEX cells
01871   for(vtkIdType id=0; id<numberOfPoints; id++)
01872     {
01873     vtkgrid->InsertNextCell(VTK_VERTEX, 1, &id);
01874     }
01875   vtkgrid->Squeeze();
01876 
01877   // in this particular case, the global ids of the cells is the same as the global ids of the points.
01878   vtkgrid->GetCellData()->SetGlobalIds(vtkgrid->GetPointData()->GetGlobalIds());
01879 
01880   return vtkgrid;
01881 }
01882 
01883 vtkMedGrid* vtkMedReader::FindGridStep(vtkMedMesh* mesh)
01884 {
01885   if(this->Internal->RealAnimationMode == vtkMedReader::PhysicalTime)
01886     {
01887     vtkMedComputeStep cs;
01888     cs.TimeOrFrequency = this->Internal->UpdateTimeStep;
01889     return mesh->FindGridStep(cs, vtkMedReader::PhysicalTime);
01890     }
01891   else if(this->Internal->RealAnimationMode == vtkMedReader::Modes)
01892     {
01893     vtkMedComputeStep cs;
01894     cs.IterationIt = MED_NO_IT;
01895     cs.TimeIt = MED_NO_DT;
01896     cs.TimeOrFrequency = MED_NO_DT;
01897     return mesh->FindGridStep(cs, vtkMedReader::Modes);
01898     }
01899   else // Iterations
01900     {
01901     vtkMedComputeStep cs;
01902     // the time is set by choosing its index in the global
01903     // array giving the available times : this->TimeIndexForIterations
01904     cs.TimeOrFrequency = (med_int)this->AvailableTimes->GetValue(
01905         (vtkIdType)this->TimeIndexForIterations);
01906     // the iteration is asked by the pipeline
01907     cs.IterationIt = (med_int)this->Internal->UpdateTimeStep;
01908     return mesh->FindGridStep(cs, vtkMedReader::Iteration);
01909     }
01910   return NULL;
01911 }
01912 
01913 void vtkMedReader::CreateMedSupports()
01914 {
01915   this->Internal->UsedSupports.clear();
01916 
01917   std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
01918       meshfileit = this->Internal->MedFiles.begin();
01919   while(meshfileit != this->Internal->MedFiles.end())
01920     {
01921     vtkMedFile* file = meshfileit->second;
01922     meshfileit++;
01923 
01924     for(int meshIndex=0; meshIndex<file->GetNumberOfMesh(); meshIndex++)
01925       {
01926       vtkMedMesh* mesh = file->GetMesh(meshIndex);
01927       vtkMedGrid* grid = this->FindGridStep(mesh);
01928       if(grid == NULL)
01929         continue;
01930 
01931       for(int entityIndex=0; entityIndex<grid->GetNumberOfEntityArray();
01932         entityIndex++)
01933         {
01934         vtkMedEntityArray* array=grid->GetEntityArray(entityIndex);
01935         if(this->GetEntityStatus(array->GetEntity())==0)
01936           {
01937           continue;
01938           }
01939 
01940         file->GetMedDriver()->LoadFamilyIds(array);
01941         for(int foeIndex=0; foeIndex<array->GetNumberOfFamilyOnEntity();
01942           foeIndex++)
01943           {
01944           vtkMedFamilyOnEntity* foe=array->GetFamilyOnEntity(foeIndex);
01945           vtkMedFamily* family=foe->GetFamily();
01946           if(this->GetFamilyStatus(mesh, family)==0)
01947             continue;
01948 
01949           // now, I look over all non-point fields to see which profiles
01950           // have to be used on points.
01951           bool selectedSupport = false;
01952 
01953           std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
01954               fieldfileit = this->Internal->MedFiles.begin();
01955           while(fieldfileit != this->Internal->MedFiles.end())
01956             {
01957             vtkMedFile* fieldfile = fieldfileit->second;
01958             fieldfileit++;
01959 
01960             for(int fieldId=0; fieldId<fieldfile->GetNumberOfField(); fieldId++)
01961               {
01962               vtkMedField* field=fieldfile->GetField(fieldId);
01963 
01964               if (!this->IsFieldSelected(field))
01965                 continue;
01966 
01967               vtkMedListOfFieldSteps steps;
01968 
01969               this->GatherFieldSteps(field, steps);
01970 
01971               vtkMedListOfFieldSteps::iterator it=steps.begin();
01972               while(it!=steps.end())
01973                 {
01974                 vtkMedFieldStep *step = *it;
01975                 step->LoadInformation();
01976                 it++;
01977 
01978                 for(int eid=0; eid<step->GetNumberOfFieldOverEntity(); eid++)
01979                   {
01980                   vtkMedFieldOverEntity* fieldOverEntity =
01981                       step->GetFieldOverEntity(eid);
01982 
01983                   for(int pid = 0; pid < fieldOverEntity->GetNumberOfFieldOnProfile(); pid++)
01984                     {
01985                     vtkMedFieldOnProfile* fop =
01986                         fieldOverEntity->GetFieldOnProfile(pid);
01987                     vtkMedProfile* prof = fop->GetProfile();
01988                     vtkMedFamilyOnEntityOnProfile* foep =
01989                         foe->GetFamilyOnEntityOnProfile(prof);
01990                     if(foep != NULL)
01991                       {
01992                       this->Internal->UsedSupports.insert(foep);
01993                       selectedSupport = true;
01994                       }
01995                     }
01996                   }
01997                 }
01998               }
01999             }
02000           // If no field use this family on entity, I nevertheless create the
02001           // support, with an empty profile.
02002           if(!selectedSupport)
02003             {
02004             vtkMedFamilyOnEntityOnProfile* foep =
02005                 foe->GetFamilyOnEntityOnProfile((vtkMedProfile*)NULL);
02006             if(foep == NULL)
02007               {
02008               foep = vtkMedFamilyOnEntityOnProfile::New();
02009               foep->SetFamilyOnEntity(foe);
02010               foep->SetProfile(NULL);
02011               foe->AddFamilyOnEntityOnProfile(foep);
02012               foep->Delete();
02013               }
02014             this->Internal->UsedSupports.insert(foep);
02015             }
02016           }
02017         }
02018       }
02019   }
02020 }
02021 
02022 bool vtkMedReader::BuildVTKSupport(
02023     vtkMedFamilyOnEntityOnProfile* foep,
02024     int doBuildSupport)
02025 {
02026 
02027   vtkMedFamilyOnEntity* foe = foep->GetFamilyOnEntity();
02028 
02029   int numProc = 1;
02030   vtkMultiProcessController* controller = vtkMultiProcessController::GetGlobalController();
02031   if (controller != NULL)
02032     {
02033     numProc = controller->GetNumberOfProcesses();
02034     }
02035 
02036   if ((foep->GetValid() == 0) && numProc == 1)
02037     {
02038     return false;
02039     }
02040 
02041   vtkMedGrid* grid=foe->GetParentGrid();
02042 
02043   vtkMedEntityArray* array=foe->GetEntityArray();
02044   vtkMedMesh* mesh=grid->GetParentMesh();
02045   vtkSmartPointer<vtkStringArray> path = vtkSmartPointer<vtkStringArray>::New();
02046   string meshName=vtkMedUtilities::SimplifyName(mesh->GetName());
02047   path->InsertNextValue(meshName);
02048   std::string finalName;
02049 
02050   if(foe->GetPointOrCell()==vtkMedUtilities::OnPoint)
02051     {
02052     path->InsertNextValue(vtkMedUtilities::OnPointName);
02053     finalName=vtkMedUtilities::SimplifyName(foe->GetFamily()->GetName());
02054     }
02055   else
02056     {
02057     path->InsertNextValue(vtkMedUtilities::OnCellName);
02058     path->InsertNextValue(vtkMedUtilities::SimplifyName(foe->GetFamily()->GetName()));
02059     finalName=vtkMedUtilities::EntityKey(array->GetEntity());
02060     }
02061 
02062   if(foep->GetProfile() != NULL)
02063     {
02064     path->InsertNextValue(finalName);
02065     finalName = foep->GetProfile()->GetName();
02066     }
02067 
02068   ostringstream progressBarTxt;
02069   for(int depth=0; depth<path->GetNumberOfValues(); depth++)
02070     {
02071     progressBarTxt<<path->GetValue(depth)<<" ";
02072     }
02073   progressBarTxt<<finalName;
02074   SetProgressText(progressBarTxt.str().c_str());
02075 
02076   vtkDataSet* cachedDataSet = NULL;
02077   if(this->Internal->DataSetCache.find(foep)!=this->Internal->DataSetCache.end())
02078     {
02079     cachedDataSet = this->Internal->DataSetCache[foep];
02080     }
02081   else
02082     {
02083     vtkDataSet* dataset = NULL;
02084     if(doBuildSupport)
02085       {
02086       if(foe->GetPointOrCell()==vtkMedUtilities::OnPoint)
02087         {
02088         dataset = this->CreateUnstructuredGridForPointSupport(foep);
02089         }
02090       else
02091         {
02092         dataset = foep->GetFamilyOnEntity()->GetParentGrid()->
02093                   CreateVTKDataSet(foep);
02094         }
02095       }
02096 
02097     if(dataset == NULL)
02098       {
02099       return false;
02100       }
02101 
02102     this->Internal->DataSetCache[foep]=dataset;
02103     cachedDataSet = dataset;
02104     if(dataset != NULL)
02105       dataset->Delete();
02106   }
02107 
02108   vtkMultiBlockDataSet* root=vtkMedUtilities::GetParent(this->GetOutput(), path);
02109   int nb=root->GetNumberOfBlocks();
02110 
02111   if(cachedDataSet != NULL)
02112     {
02113     vtkDataSet* realDataSet=cachedDataSet->NewInstance();
02114     root->SetBlock(nb, realDataSet);
02115     realDataSet->Delete();
02116 
02117     root->GetMetaData(nb)->Set(vtkCompositeDataSet::NAME(), finalName.c_str());
02118     realDataSet->ShallowCopy(cachedDataSet);
02119 
02120     this->Internal->DataSetCache[foep]=cachedDataSet;
02121     this->Internal->CurrentDataSet[foep]=realDataSet;
02122 
02123     path->InsertNextValue(finalName);
02124     path->SetName("BLOCK_NAME");
02125     realDataSet->GetFieldData()->AddArray(path);
02126     realDataSet->GetInformation()->Remove(vtkMedUtilities::BLOCK_NAME());
02127     for(int depth=0; depth<path->GetNumberOfValues(); depth++)
02128       {
02129       realDataSet->GetInformation()->Set(vtkMedUtilities::BLOCK_NAME(),
02130                                          path->GetValue(depth), depth);
02131       }
02132     }
02133 }
02134 
02135 void vtkMedReader::MapFieldOnSupport(vtkMedFieldOnProfile* fop,
02136                                      vtkMedFamilyOnEntityOnProfile* foep,
02137                                      int doCreateField)
02138 {
02139   bool cached = false;
02140 
02141   if(this->Internal->FieldCache.find(foep)
02142       !=this->Internal->FieldCache.end())
02143     {
02144     map<vtkMedFieldOnProfile*, VTKField>& fieldCache =
02145         this->Internal->FieldCache[foep];
02146     if(fieldCache.find(fop)!=fieldCache.end())
02147       {
02148       cached=true;
02149       }
02150     }
02151 
02152   if(!cached)
02153     {
02154     this->CreateVTKFieldOnSupport(fop, foep, doCreateField);
02155     }
02156   this->SetVTKFieldOnSupport(fop, foep);
02157 }
02158 
02159 void vtkMedReader::MapFieldsOnSupport(vtkMedFamilyOnEntityOnProfile* foep,
02160                                       int doCreateField)
02161 {
02162   // now loop over all fields to map it to the created grids
02163   std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
02164       fieldfileit = this->Internal->MedFiles.begin();
02165   while(fieldfileit != this->Internal->MedFiles.end())
02166     {
02167     vtkMedFile* file = fieldfileit->second;
02168     fieldfileit++;
02169 
02170     for(int fieldId=0; fieldId<file->GetNumberOfField(); fieldId++)
02171       {
02172       vtkMedField* field=file->GetField(fieldId);
02173     
02174       if(strcmp(foep->GetFamilyOnEntity()->
02175                 GetParentGrid()->GetParentMesh()->GetName(),
02176                 field->GetMeshName()) != 0)
02177         continue;
02178 
02179       if(strcmp(foep->GetFamilyOnEntity()->
02180                 GetParentGrid()->GetParentMesh()->GetName(),
02181                 field->GetMeshName()) != 0)
02182         continue;
02183 
02184       if (!this->IsFieldSelected(field))
02185         continue;
02186       
02187       vtkMedListOfFieldSteps steps;
02188 
02189       this->GatherFieldSteps(field, steps);
02190       
02191       vtkMedListOfFieldSteps::iterator it=steps.begin();
02192       while(it!=steps.end())
02193         {
02194         vtkMedFieldStep *step = *it;
02195         step->LoadInformation();
02196         it++;
02197       
02198         for(int eid=0; eid<step->GetNumberOfFieldOverEntity(); eid++)
02199           {
02200           vtkMedFieldOverEntity* fieldOverEntity = step->GetFieldOverEntity(eid);
02201           for(int pid = 0; pid < fieldOverEntity->GetNumberOfFieldOnProfile(); pid++)
02202             {
02203             vtkMedFieldOnProfile* fop = fieldOverEntity->GetFieldOnProfile(pid);
02204             if(foep->CanMapField(fop))
02205               {
02206               this->MapFieldOnSupport(fop, foep, doCreateField);
02207               }
02208             }
02209           }
02210         }
02211       }
02212     }
02213 }
02214 
02215 void vtkMedReader::GatherFieldSteps(vtkMedField* field,
02216                                     vtkMedListOfFieldSteps& steps)
02217 {
02218   if(this->Internal->RealAnimationMode == vtkMedReader::PhysicalTime)
02219     {
02220     vtkMedComputeStep cs;
02221     cs.TimeOrFrequency = this->Internal->UpdateTimeStep;
02222     vtkMedFieldStep* fs =
02223         field->FindFieldStep(cs, vtkMedReader::PhysicalTime);
02224     if(fs != NULL)
02225       steps.push_back(fs);
02226     }
02227   else if(this->Internal->RealAnimationMode == vtkMedReader::Iteration)
02228     {
02229     vtkMedComputeStep cs;
02230     cs.IterationIt = (med_int)this->Internal->UpdateTimeStep;
02231     cs.TimeOrFrequency = (med_int)this->AvailableTimes->GetValue(
02232         (vtkIdType)this->TimeIndexForIterations);
02233     vtkMedFieldStep* fs =
02234         field->FindFieldStep(cs, vtkMedReader::Iteration);
02235     if(fs != NULL)
02236       steps.push_back(fs);
02237     }
02238   else // modes
02239     {
02240     for(int modeid = 0; modeid < this->Frequencies->GetNumberOfKey(); modeid++)
02241       {
02242       if(this->Frequencies->GetKeyStatus(
02243           this->Frequencies->GetKey(modeid)) == 0)
02244         {
02245         continue;
02246         }
02247 
02248       vtkMedComputeStep cs;
02249       cs.TimeOrFrequency = this->AvailableTimes->GetValue(modeid);
02250       vtkMedFieldStep* fs =
02251           field->FindFieldStep(cs, vtkMedReader::PhysicalTime);
02252       if(fs != NULL)
02253         steps.push_back(fs);
02254       }
02255     }
02256 }
02257 
02258 void vtkMedReader::SetVTKFieldOnSupport(vtkMedFieldOnProfile* fop,
02259     vtkMedFamilyOnEntityOnProfile* foep)
02260 {
02261   vtkMedFamilyOnEntity* foe = foep->GetFamilyOnEntity();
02262   vtkMedFieldOverEntity* fieldOverEntity = fop->GetParentFieldOverEntity();
02263   vtkMedFieldStep* step = fieldOverEntity->GetParentStep();
02264   vtkMedField* field = step->GetParentField();
02265   
02266   const vtkMedComputeStep& cs = step->GetComputeStep();
02267 
02268   vtkDataSet* ds=this->Internal->CurrentDataSet[foep];
02269   if (ds == NULL)
02270     {
02271   // ds == NULL could arrive is some cases when working in parallel
02272   vtkWarningMacro( "--- vtkMedReader::SetVTKFieldOnSupport: ds is NULL !!");
02273   return;
02274     }
02275 
02276   VTKField& vtkfield=this->Internal->FieldCache[foep][fop];
02277 
02278   std::string name=vtkMedUtilities::SimplifyName(field->GetName());
02279   std::string vectorname = name+"_Vector";
02280   if(this->AnimationMode==Modes)
02281     {
02282     double freq = cs.TimeOrFrequency;
02283     int index = this->GetFrequencyIndex(freq);
02284     name += string(" ") + vtkMedUtilities::GetModeKey(index, freq,
02285       this->AvailableTimes->GetNumberOfTuples()-1);
02286     vectorname += string(" ") + vtkMedUtilities::GetModeKey(index, freq,
02287       this->AvailableTimes->GetNumberOfTuples()-1);
02288     }
02289 
02290   vtkfield.DataArray->SetName(name.c_str());
02291   vtkfield.DataArray->Squeeze();
02292   if(vtkfield.Vectors != NULL)
02293     {
02294     vtkfield.Vectors->SetName(vectorname.c_str());
02295     vtkfield.Vectors->Squeeze();
02296     }
02297   if(vtkfield.QuadratureIndexArray!=NULL)
02298     {
02299     vtkfield.QuadratureIndexArray->Squeeze();
02300     }
02301 
02302   if(foe->GetPointOrCell()==vtkMedUtilities::OnCell)
02303     {
02304     if(field->GetFieldType()==vtkMedField::PointField)
02305       {
02306       if(vtkfield.DataArray->GetNumberOfTuples()!=ds->GetNumberOfPoints())
02307         {
02308          vtkDebugMacro("the data array " << vtkfield.DataArray->GetName()
02309                       << " do not have the good number of tuples");
02310         return;
02311         }
02312       ds->GetPointData()->AddArray(vtkfield.DataArray);
02313       if(vtkfield.Vectors != NULL && this->GenerateVectors)
02314         {
02315         ds->GetPointData()->AddArray(vtkfield.Vectors);
02316         ds->GetPointData()->SetActiveVectors(vtkfield.Vectors->GetName());
02317         }
02318       switch (vtkfield.DataArray->GetNumberOfComponents())
02319         {
02320         case 1:
02321         ds->GetPointData()->SetActiveScalars(vtkfield.DataArray->GetName());
02322         break;
02323         case 3:
02324         ds->GetPointData()->SetActiveVectors(vtkfield.DataArray->GetName());
02325         break;
02326         }
02327       // if the data set is only composed of VTK_VERTEX cells,
02328       // and no field called with the same name exist on cells,
02329       // map this field to cells too
02330       if(foe->GetVertexOnly()==1 && ds->GetCellData()->GetArray(
02331               vtkfield.DataArray->GetName())==NULL)
02332         {
02333         ds->GetCellData()->AddArray(vtkfield.DataArray);
02334         if(vtkfield.Vectors != NULL && this->GenerateVectors)
02335           {
02336           ds->GetCellData()->AddArray(vtkfield.Vectors);
02337           ds->GetCellData()->SetActiveVectors(vtkfield.Vectors->GetName());
02338           }
02339         switch (vtkfield.DataArray->GetNumberOfComponents())
02340           {
02341           case 1:
02342           ds->GetCellData()->SetActiveScalars(vtkfield.DataArray->GetName());
02343           break;
02344           case 3:
02345           ds->GetCellData()->SetActiveVectors(vtkfield.DataArray->GetName());
02346           break;
02347           }
02348         }
02349       }
02350     if(field->GetFieldType()==vtkMedField::CellField)
02351       {
02352       if((this->Internal->NumberOfPieces == 1) && vtkfield.DataArray->GetNumberOfTuples()!=ds->GetNumberOfCells()  )
02353         {
02354         vtkDebugMacro("the data array " << vtkfield.DataArray->GetName()
02355                       << " do not have the good number of tuples"
02356                       << " ncell=" << ds->GetNumberOfCells()
02357                       << " ntuple=" << vtkfield.DataArray->GetNumberOfTuples());
02358         return;
02359         }
02360       // In case we are in parallel and our process does not contain the data
02361       if(ds->GetNumberOfCells()==0)
02362         {
02363         return;
02364         }
02365       ds->GetCellData()->AddArray(vtkfield.DataArray);
02366 
02367       if(vtkfield.Vectors != NULL && this->GenerateVectors)
02368         {
02369         ds->GetCellData()->AddArray(vtkfield.Vectors);
02370         ds->GetCellData()->SetActiveVectors(vtkfield.Vectors->GetName());
02371         }
02372       switch (vtkfield.DataArray->GetNumberOfComponents())
02373         {
02374         case 1:
02375         ds->GetCellData()->SetActiveScalars(vtkfield.DataArray->GetName());
02376         break;
02377         case 3:
02378         ds->GetCellData()->SetActiveVectors(vtkfield.DataArray->GetName());
02379         break;
02380         }
02381       // if the data set is only composed of VTK_VERTEX cells,
02382       // and no field called with the same name exist on points,
02383       // map this field to points too
02384       if(foe->GetVertexOnly()==1 && ds->GetPointData()->GetArray(
02385               vtkfield.DataArray->GetName())==NULL)
02386         {
02387         ds->GetPointData()->AddArray(vtkfield.DataArray);
02388         if(vtkfield.Vectors != NULL && this->GenerateVectors)
02389           {
02390           ds->GetPointData()->AddArray(vtkfield.Vectors);
02391           ds->GetPointData()->SetActiveVectors(vtkfield.Vectors->GetName());
02392           }
02393         switch (vtkfield.DataArray->GetNumberOfComponents())
02394           {
02395           case 1:
02396           ds->GetPointData()->SetActiveScalars(vtkfield.DataArray->GetName());
02397           break;
02398           case 3:
02399           ds->GetPointData()->SetActiveVectors(vtkfield.DataArray->GetName());
02400           break;
02401           }
02402         }
02403       }
02404     if(field->GetFieldType()==vtkMedField::QuadratureField ||
02405        field->GetFieldType()==vtkMedField::ElnoField )
02406       {
02407       vtkIdType ncells=ds->GetNumberOfCells();
02408       vtkIdType nid=vtkfield.QuadratureIndexArray->GetNumberOfTuples();
02409       vtkIdType nda=vtkfield.DataArray->GetNumberOfTuples();
02410       if((nid!=ncells) && (this->Internal->NumberOfPieces == 1))
02411         {
02412         vtkDebugMacro(
02413             "There should be as many quadrature index values as there are cells");
02414         return;
02415         }
02416       else
02417         {
02418       if (ncells == 0)
02419       {
02420         vtkfield.DataArray->SetNumberOfTuples( 0 );
02421         vtkfield.DataArray->Squeeze();
02422       }
02423       if (nid>ncells)  // PROBABLY NOT NECESSARY
02424       {
02425         vtkfield.QuadratureIndexArray->SetNumberOfTuples(ncells);
02426         int nquad=fop->GetNumberOfIntegrationPoint();
02427         vtkfield.DataArray->SetNumberOfTuples( nquad * ds->GetNumberOfCells() );
02428         vtkfield.DataArray->Squeeze();
02429       }
02430         ds->GetFieldData()->AddArray(vtkfield.DataArray);
02431         ds->GetCellData()->AddArray(vtkfield.QuadratureIndexArray);
02432 
02433         nid=vtkfield.QuadratureIndexArray->GetNumberOfTuples();
02434         nda=vtkfield.DataArray->GetNumberOfTuples();
02435 
02436         if(vtkfield.Vectors != NULL && this->GenerateVectors)
02437           {
02438           ds->GetFieldData()->AddArray(vtkfield.Vectors);
02439           }
02440         }
02441       }
02442     }//support OnCell
02443   else
02444     {//support OnPoint
02445     if(vtkfield.DataArray->GetNumberOfTuples()!=ds->GetNumberOfPoints())
02446       {
02447       vtkDebugMacro("the data array " << vtkfield.DataArray->GetName() << " do not have the good number of tuples");
02448       return;
02449       }
02450     ds->GetPointData()->AddArray(vtkfield.DataArray);
02451     if(vtkfield.Vectors != NULL && this->GenerateVectors)
02452       {
02453       ds->GetPointData()->AddArray(vtkfield.Vectors);
02454       ds->GetPointData()->SetActiveVectors(vtkfield.Vectors->GetName());
02455       }
02456     switch (vtkfield.DataArray->GetNumberOfComponents())
02457       {
02458       case 1:
02459       ds->GetPointData()->SetActiveScalars(vtkfield.DataArray->GetName());
02460       break;
02461       case 3:
02462       ds->GetPointData()->SetActiveVectors(vtkfield.DataArray->GetName());
02463       break;
02464       }
02465     // all the VTK_VERTEX created cells have the same order than the points,
02466     // I can safely map the point array to the cells in this particular case.
02467     ds->GetCellData()->AddArray(vtkfield.DataArray);
02468     if(vtkfield.Vectors != NULL && this->GenerateVectors)
02469       {
02470       ds->GetCellData()->AddArray(vtkfield.Vectors);
02471       ds->GetCellData()->SetActiveVectors(vtkfield.Vectors->GetName());
02472       }
02473     switch (vtkfield.DataArray->GetNumberOfComponents())
02474       {
02475       case 1:
02476       ds->GetCellData()->SetActiveScalars(vtkfield.DataArray->GetName());
02477       break;
02478       case 3:
02479       ds->GetCellData()->SetActiveVectors(vtkfield.DataArray->GetName());
02480       break;
02481       }
02482     }
02483 }
02484 
02485 void vtkMedReader::InitializeQuadratureOffsets(vtkMedFieldOnProfile* fop,
02486     vtkMedFamilyOnEntityOnProfile* foep)
02487 {
02488   vtkMedFamilyOnEntity* foe = foep->GetFamilyOnEntity();
02489 
02490   vtkMedFieldOverEntity* fieldOverEntity = fop->GetParentFieldOverEntity();
02491   vtkMedFieldStep* step = fieldOverEntity->GetParentStep();
02492   vtkMedField *field= step->GetParentField();
02493 
02494   // then I compute the quadrature key if needed, and try to find the quadrature offsets.
02495   if(this->Internal->QuadratureOffsetCache.find(foep)
02496       ==this->Internal->QuadratureOffsetCache.end())
02497     this->Internal->QuadratureOffsetCache[foep]=map<LocalizationKey,
02498         vtkSmartPointer<vtkIdTypeArray> > ();
02499 
02500   map<LocalizationKey, vtkSmartPointer<vtkIdTypeArray> >& quadOffsets=
02501       this->Internal->QuadratureOffsetCache[foep];
02502 
02503   LocalizationKey quadkey=this->GetLocalizationKey(fop);
02504 
02505   if(quadOffsets.find(quadkey)!=quadOffsets.end())
02506     {// the quadrature offset array has already been created
02507     return;
02508     }
02509 
02510   vtkIdTypeArray* qoffsets=vtkIdTypeArray::New();
02511   quadOffsets[quadkey]=qoffsets;
02512   qoffsets->Delete();
02513 
02514   ostringstream sstr;
02515   if(field->GetFieldType() == vtkMedField::ElnoField)
02516     {
02517     qoffsets->GetInformation()->Set(vtkMedUtilities::ELNO(), 1);
02518     sstr<<"ELNO";
02519     }
02520   else if(field->GetFieldType() == vtkMedField::QuadratureField)
02521     {
02522     qoffsets->GetInformation()->Set(vtkMedUtilities::ELGA(), 1);
02523     sstr<<"ELGA";
02524     }
02525   else
02526     {
02527     sstr<<"QuadraturePointOffset";
02528     }
02529 
02530   qoffsets->SetName(sstr.str().c_str());
02531 
02532   vtkSmartPointer<vtkMedLocalization> loc=
02533       this->GetLocalization(fop->GetLocalizationName());
02534 
02535   if(loc == NULL)
02536     {
02537     if(fop->GetNumberOfIntegrationPoint()==1)
02538       {// cell-centered fields can be seen as quadrature fields with 1
02539       // quadrature point at the center
02540       vtkMedLocalization* center=vtkMedLocalization::New();
02541       loc=center;
02542       center->Delete();
02543       center->BuildCenter(fieldOverEntity->GetEntity().GeometryType);
02544       }
02545     else if(loc == NULL && field->GetFieldType() == vtkMedField::ElnoField)
02546       {// ELNO fields have no vtkMedLocalization,
02547       // I need to create a dummy one
02548       vtkMedLocalization* elnodef=vtkMedLocalization::New();
02549       loc=elnodef;
02550       elnodef->Delete();
02551       elnodef->BuildELNO(fieldOverEntity->GetEntity().GeometryType);
02552       }
02553     else
02554       {
02555       vtkErrorMacro("Cannot find localization of quadrature field "
02556                     << field->GetName());
02557       }
02558     }
02559   this->AddQuadratureSchemeDefinition(qoffsets->GetInformation(), loc);
02560 }
02561 
02562 void vtkMedReader::CreateVTKFieldOnSupport(vtkMedFieldOnProfile* fop,
02563     vtkMedFamilyOnEntityOnProfile* foep, int doCreateField)
02564 {
02565   vtkMedFamilyOnEntity* foe = foep->GetFamilyOnEntity();
02566   vtkMedFieldOverEntity* fieldOverEntity = fop->GetParentFieldOverEntity();
02567   vtkMedFieldStep* step = fieldOverEntity->GetParentStep();
02568   vtkMedField* field = step->GetParentField();
02569 
02570   if(vtkMedUnstructuredGrid::SafeDownCast(
02571       foep->GetFamilyOnEntity()->GetParentGrid()) != NULL)
02572     {
02573     fop->Load(MED_COMPACT_STMODE);
02574     }
02575   else
02576     {
02577     fop->Load(MED_GLOBAL_STMODE);
02578     }
02579 
02580   vtkMedIntArray* profids=NULL;
02581 
02582   vtkMedProfile* profile=fop->GetProfile();
02583 
02584   if(profile)
02585     {
02586     profile->Load();
02587     profids=profile->GetIds();
02588     }//has profile
02589 
02590   VTKField& vtkfield=this->Internal->FieldCache[foep][fop];
02591 
02592   bool shallowok = true;
02593 
02594   // for structured grid, the values are loaded globally, and cells which are
02595   // not on the profile or not on the family are blanked out.
02596   // shallow copy is therefore always possible
02597   if(vtkMedUnstructuredGrid::SafeDownCast(
02598       foep->GetFamilyOnEntity()->GetParentGrid()) != NULL)
02599     {
02600     shallowok = foep->CanShallowCopy(fop);
02601     }
02602 
02603   if(shallowok)
02604     {
02605     vtkfield.DataArray = fop->GetData();
02606     }
02607   else
02608     {
02609     vtkDataArray* data=vtkMedUtilities::NewArray(field->GetDataType());
02610     vtkfield.DataArray=data;
02611     data->Delete();
02612     vtkfield.DataArray->SetNumberOfComponents(field->GetNumberOfComponent());
02613     }
02614 
02615   for(int comp=0; comp<field->GetNumberOfComponent(); comp++)
02616     {
02617     vtkfield.DataArray->SetComponentName(comp, vtkMedUtilities::SimplifyName(
02618         field->GetComponentName()->GetValue(comp)).c_str());
02619     }
02620 
02621   bool createOffsets=false;
02622   if(field->GetFieldType()==vtkMedField::QuadratureField ||
02623      field->GetFieldType()==vtkMedField::ElnoField )
02624     {
02625     this->InitializeQuadratureOffsets(fop, foep);
02626 
02627     LocalizationKey quadKey = this->GetLocalizationKey(fop);
02628     vtkfield.QuadratureIndexArray
02629         =this->Internal->QuadratureOffsetCache[foep][quadKey];
02630     vtkDataSet* ds = this->Internal->CurrentDataSet[foep];
02631 
02632     vtkfield.DataArray->GetInformation()->Set(
02633         vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME(),
02634         vtkfield.QuadratureIndexArray->GetName());
02635     vtkfield.QuadratureIndexArray->GetInformation()->Set(
02636         vtkAbstractArray::GUI_HIDE(), 1);
02637 
02638     if(vtkfield.QuadratureIndexArray->GetNumberOfTuples()
02639         !=ds->GetNumberOfCells())
02640       {
02641       vtkfield.QuadratureIndexArray->SetNumberOfTuples(0);
02642       createOffsets=true;
02643       }
02644     }
02645 
02646   if(!doCreateField)
02647     return;
02648 
02649   if(shallowok)
02650     {
02651     // build the quadrature offset array if needed
02652     if(createOffsets)
02653       {
02654       vtkIdType noffsets;
02655       int nquad=fop->GetNumberOfIntegrationPoint();
02656       noffsets=fop->GetData()->GetNumberOfTuples()/nquad;
02657 
02658       vtkDataSet* ds=this->Internal->CurrentDataSet[foep];
02659       if(noffsets!=ds->GetNumberOfCells())
02660         {
02661         vtkErrorMacro(
02662             "number of quadrature offsets (" << noffsets << ") must "
02663             << "match the number of cells (" << ds->GetNumberOfCells() << ")!");
02664         }
02665 
02666       vtkfield.QuadratureIndexArray->SetNumberOfTuples(noffsets);
02667       for(vtkIdType id=0; id<noffsets; id++)
02668         {
02669         vtkfield.QuadratureIndexArray->SetValue(id, nquad*id);
02670         }
02671       }
02672 
02673     }
02674   else if(foe->GetPointOrCell() == vtkMedUtilities::OnCell
02675      && field->GetFieldType() != vtkMedField::PointField)
02676     {
02677     // Cell-centered field on cell support
02678     int nquad = 1;
02679     if (field->GetFieldType()==vtkMedField::QuadratureField ||
02680         field->GetFieldType()==vtkMedField::ElnoField)
02681       {
02682       nquad=fop->GetNumberOfIntegrationPoint();
02683       }
02684     vtkMedIntArray* profids=NULL;
02685     if (profile)
02686       {
02687       profids=profile->GetIds();
02688       }
02689     vtkIdType maxIndex=fop->GetData()->GetNumberOfTuples();
02690     maxIndex/=nquad;
02691     vtkIdType quadIndex = 0;
02692 
02693     for (vtkIdType id = 0; id<maxIndex; id++)
02694       {
02695       vtkIdType realIndex = (profids!=NULL ? profids->GetValue(id)-1 : id);
02696       if (!foep->KeepCell(realIndex))
02697         continue;
02698 
02699       if (field->GetFieldType()==vtkMedField::QuadratureField ||
02700           field->GetFieldType()==vtkMedField::ElnoField)
02701         {
02702         for (int q = 0; q<nquad; q++)
02703           {
02704           vtkfield.DataArray->InsertNextTuple(nquad*realIndex+q,
02705               fop->GetData());
02706           }
02707         if (createOffsets)
02708           {
02709           vtkfield.QuadratureIndexArray->InsertNextValue(quadIndex);
02710           quadIndex+=nquad;
02711           }
02712         }
02713       else
02714         {
02715         vtkfield.DataArray->InsertNextTuple(id, fop->GetData());
02716         }
02717       }//copy all tuples
02718     vtkfield.DataArray->Squeeze();
02719     vtkDataSet* ds = this->Internal->CurrentDataSet[foep];
02720     }
02721   else if(foe->GetPointOrCell() == vtkMedUtilities::OnCell
02722      && field->GetFieldType() == vtkMedField::PointField)
02723     {// point field on cell support
02724     vtkMedIntArray* profids=NULL;
02725 
02726     vtkMedProfile* profile=fop->GetProfile();
02727 
02728     if(profile)
02729       {
02730       profile->Load();
02731 
02732       profids=profile->GetIds();
02733       }//has profile
02734 
02735     vtkIdType maxId=fop->GetData()->GetNumberOfTuples();
02736 
02737     for(vtkIdType id=0; id<maxId; id++)
02738       {
02739       // if I have a profile, then I should insert the value at the position given in the profile.
02740       vtkIdType destIndex;
02741       if(profids!=NULL)
02742         {
02743         destIndex=profids->GetValue(id)-1; // -1 because med indices start at 1
02744         }
02745       else
02746         {
02747         destIndex=id;
02748         }
02749 
02750       if(!foep->KeepPoint(destIndex))
02751         continue;
02752       // if I use the med2VTKIndex, then the index to insert
02753       // this value is given by the map.
02754       destIndex = foep->GetVTKPointIndex(destIndex);
02755       vtkfield.DataArray->InsertTuple(destIndex, id, fop->GetData());
02756       }
02757     vtkfield.DataArray->Squeeze();
02758     }// point field on cell support
02759   else if(foe->GetPointOrCell() == vtkMedUtilities::OnPoint &&
02760      field->GetFieldType() == vtkMedField::PointField)
02761     {//support OnPoint
02762 
02763     vtkIdType maxId = fop->GetData()->GetNumberOfTuples();
02764 
02765     for(vtkIdType id=0; id<maxId; id++)
02766       {
02767       vtkIdType realIndex=id;
02768       if(profids!=NULL)
02769         {
02770         realIndex=profids->GetValue(id)-1; // -1 because med indices start at 1
02771         }
02772 
02773       if(!foep->KeepPoint(realIndex))
02774         continue;
02775 
02776       vtkfield.DataArray->InsertNextTuple(fop->GetData()->GetTuple(id));
02777       }
02778     vtkfield.DataArray->Squeeze();
02779     }// support on point
02780 
02781   // now generate the vector field if asked for
02782   if(this->GenerateVectors)
02783     {
02784     int ncomp = vtkfield.DataArray->GetNumberOfComponents();
02785     if(ncomp > 1 && ncomp != 3)
02786       {
02787       vtkDataArray* vectors = vtkfield.DataArray->NewInstance();
02788       vectors->SetNumberOfComponents(3);
02789       vectors->SetNumberOfTuples(vtkfield.DataArray->GetNumberOfTuples());
02790       vtkfield.Vectors = vectors;
02791       vectors->Delete();
02792 
02793       vectors->CopyInformation(vtkfield.DataArray->GetInformation());
02794 
02795       if(ncomp < 3)
02796         vectors->SetComponentName(2, "Z");
02797       else
02798         vectors->SetComponentName(2, vtkfield.DataArray->GetComponentName(2));
02799 
02800       vectors->SetComponentName(1, vtkfield.DataArray->GetComponentName(1));
02801       vectors->SetComponentName(0, vtkfield.DataArray->GetComponentName(0));
02802 
02803       int tuplesize = (ncomp > 3? ncomp: 3);
02804       double* tuple = new double[tuplesize];
02805       for(int tid=0; tid<tuplesize; tid++)
02806         {
02807         tuple[tid] = 0.0;
02808         }
02809 
02810       for(vtkIdType id=0; id < vtkfield.DataArray->GetNumberOfTuples(); id++)
02811         {
02812         vtkfield.DataArray->GetTuple(id, tuple);
02813         vectors->SetTuple(id, tuple);
02814         }
02815       }
02816     }
02817 }
02818 
02819 int vtkMedReader::HasMeshAnyCellSelectedFamily(vtkMedMesh* mesh)
02820 {
02821   int nfam = mesh->GetNumberOfCellFamily();
02822   for (int famid = 0; famid<nfam; famid++)
02823     {
02824     vtkMedFamily* fam = mesh->GetFamily(famid);
02825     if (fam->GetPointOrCell()!=vtkMedUtilities::OnCell||!this->GetFamilyStatus(
02826         mesh, fam))
02827       continue;
02828     return true;
02829     }
02830   return false;
02831 }
02832 
02833 int vtkMedReader::HasMeshAnyPointSelectedFamily(vtkMedMesh* mesh)
02834 {
02835   int nfam = mesh->GetNumberOfCellFamily();
02836   for (int famid = 0; famid<nfam; famid++)
02837     {
02838     vtkMedFamily* fam = mesh->GetFamily(famid);
02839     if (fam->GetPointOrCell()!=vtkMedUtilities::OnPoint
02840         ||!this->GetFamilyStatus(mesh, fam))
02841       continue;
02842     return true;
02843     }
02844   return false;
02845 }
02846 
02847 void vtkMedReader::BuildSIL(vtkMutableDirectedGraph* sil)
02848 {
02849   if(sil==NULL)
02850     return;
02851 
02852   sil->Initialize();
02853 
02854   vtkSmartPointer<vtkVariantArray> childEdge=
02855       vtkSmartPointer<vtkVariantArray>::New();
02856   childEdge->InsertNextValue(0);
02857 
02858   vtkSmartPointer<vtkVariantArray> crossEdge=
02859       vtkSmartPointer<vtkVariantArray>::New();
02860   crossEdge->InsertNextValue(1);
02861 
02862   // CrossEdge is an edge linking hierarchies.
02863   vtkUnsignedCharArray* crossEdgesArray=vtkUnsignedCharArray::New();
02864   crossEdgesArray->SetName("CrossEdges");
02865   sil->GetEdgeData()->AddArray(crossEdgesArray);
02866   crossEdgesArray->Delete();
02867   vtkstd::deque<vtkstd::string> names;
02868 
02869   // Now build the hierarchy.
02870   vtkIdType rootId=sil->AddVertex();
02871   names.push_back("SIL");
02872 
02873   // Add a global entry to encode global names for the families
02874   vtkIdType globalFamilyRoot=sil->AddChild(rootId, childEdge);
02875   names.push_back("Families");
02876 
02877   // Add a global entry to encode global names for the families
02878   vtkIdType globalGroupRoot=sil->AddChild(rootId, childEdge);
02879   names.push_back("Groups");
02880 
02881   // Add the groups subtree
02882   vtkIdType groupsRoot=sil->AddChild(rootId, childEdge);
02883   names.push_back("GroupTree");
02884 
02885   // Add a global entry to encode names for the cell types
02886   vtkIdType globalEntityRoot=sil->AddChild(rootId, childEdge);
02887   names.push_back("Entity");
02888 
02889   // Add the cell types subtree
02890   vtkIdType entityTypesRoot=sil->AddChild(rootId, childEdge);
02891   names.push_back("EntityTree");
02892 
02893   // this is the map that keep added cell types
02894   map<vtkMedEntity, vtkIdType> entityMap;
02895   map<med_entity_type, vtkIdType> typeMap;
02896 
02897   std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
02898       meshfileit = this->Internal->MedFiles.begin();
02899   while(meshfileit != this->Internal->MedFiles.end())
02900     {
02901     vtkMedFile* file = meshfileit->second;
02902     meshfileit++;
02903 
02904     int numMeshes=file->GetNumberOfMesh();
02905     for(int meshIndex=0; meshIndex<numMeshes; meshIndex++)
02906       {
02907       vtkMedMesh* mesh = file->GetMesh(meshIndex);
02908       vtkMedGrid* grid = this->FindGridStep(mesh);
02909 
02910       if(grid == NULL)
02911         continue;
02912 
02913       // add all entities
02914       for(int entityIndex=0; entityIndex<grid->GetNumberOfEntityArray(); entityIndex++)
02915         {
02916         vtkMedEntityArray* array=grid->GetEntityArray(entityIndex);
02917         vtkMedEntity entity = array->GetEntity();
02918 
02919         if(entityMap.find(entity)==entityMap.end())
02920           {
02921           vtkIdType entityGlobalId=sil->AddChild(globalEntityRoot, childEdge);
02922           names.push_back(vtkMedUtilities::EntityKey(entity));
02923 
02924           vtkIdType typeId;
02925           if(typeMap.find(entity.EntityType)==typeMap.end())
02926             {
02927             typeId=sil->AddChild(entityTypesRoot, childEdge);
02928             names.push_back(vtkMedUtilities::EntityName(entity.EntityType));
02929             typeMap[entity.EntityType]=typeId;
02930             }
02931           else
02932             {
02933             typeId=typeMap[entity.EntityType];
02934             }
02935           vtkIdType entityId=sil->AddChild(typeId, childEdge);
02936           names.push_back(entity.GeometryName);
02937 
02938           sil->AddEdge(entityId, entityGlobalId, crossEdge);
02939 
02940           entityMap[entity]=entityId;
02941           }
02942         }
02943 
02944       vtkIdType meshGroup = sil->AddChild(groupsRoot, childEdge);
02945       names.push_back(vtkMedUtilities::SimplifyName(mesh->GetName()));
02946 
02947       // add the two OnPoint and OnCell entries, for groups and for families
02948       vtkIdType meshCellGroups = sil->AddChild(meshGroup, childEdge);
02949       names.push_back(vtkMedUtilities::OnCellName);
02950 
02951       vtkIdType meshPointGroups = sil->AddChild(meshGroup, childEdge);
02952       names.push_back(vtkMedUtilities::OnPointName);
02953 
02954       // this maps will keep all added groups on nodes/cells of this mesh
02955       map<string, vtkIdType> nodeGroupMap;
02956       map<string, vtkIdType> cellGroupMap;
02957 
02958       // add all families
02959       for(int famIndex=0; famIndex<mesh->GetNumberOfFamily(); famIndex++)
02960         {
02961         vtkMedFamily* family=mesh->GetFamily(famIndex);
02962 
02963         vtkIdType globalFamilyId = sil->AddChild(globalFamilyRoot, childEdge);
02964         names.push_back(vtkMedUtilities::FamilyKey(mesh->GetName(),
02965                                                    family->GetPointOrCell(),
02966                                                    family->GetName()));
02967 
02968         // family with Id 0 is both on cell and on points, so add it to both
02969         map<string, vtkIdType> & groupMap=(family->GetPointOrCell()
02970             ==vtkMedUtilities::OnPoint? nodeGroupMap: cellGroupMap);
02971 
02972         vtkIdType groupRootId =
02973             (family->GetPointOrCell()==vtkMedUtilities::OnPoint?
02974              meshPointGroups : meshCellGroups);
02975 
02976         // add all the groups of this family
02977         for(vtkIdType groupIndex=0; groupIndex<family->GetNumberOfGroup();
02978           groupIndex++)
02979           {
02980           vtkMedGroup* group=family->GetGroup(groupIndex);
02981 
02982           vtkIdType familyGroupId = sil->AddChild(globalFamilyId, childEdge);
02983           names.push_back(vtkMedUtilities::FamilyKey(
02984               mesh->GetName(), family->GetPointOrCell(),
02985               family->GetName()));
02986 
02987           vtkIdType groupGlobalId;
02988           if(groupMap.find(group->GetName())==groupMap.end())
02989             {
02990             vtkIdType groupLocalId;
02991             groupLocalId=sil->AddChild(groupRootId, childEdge);
02992             names.push_back(vtkMedUtilities::SimplifyName(group->GetName()));
02993 
02994             groupGlobalId=sil->AddChild(globalGroupRoot, childEdge);
02995             names.push_back(vtkMedUtilities::GroupKey(
02996                 mesh->GetName(), family->GetPointOrCell(),
02997                 group->GetName()));
02998             groupMap[group->GetName()]=groupGlobalId;
02999 
03000             sil->AddEdge(groupLocalId, groupGlobalId, crossEdge);
03001             }
03002           vtkIdType groupId = groupMap[group->GetName()];
03003           sil->AddEdge(familyGroupId, groupId, childEdge);
03004 
03005           }//groupIndex
03006         }//famIndex
03007       }//meshIndex
03008     }// file iterator
03009 
03010   // This array is used to assign names to nodes.
03011   vtkStringArray* namesArray=vtkStringArray::New();
03012   namesArray->SetName("Names");
03013   namesArray->SetNumberOfTuples(sil->GetNumberOfVertices());
03014   sil->GetVertexData()->AddArray(namesArray);
03015   namesArray->Delete();
03016   vtkstd::deque<vtkstd::string>::iterator iter;
03017   vtkIdType cc;
03018   for(cc=0, iter=names.begin(); iter!=names.end(); ++iter, ++cc)
03019     {
03020     namesArray->SetValue(cc, (*iter).c_str());
03021     }
03022 }
03023 
03024 void vtkMedReader::ClearMedSupports()
03025 {
03026   this->Internal->DataSetCache.clear();
03027   //this->Internal->Med2VTKPointIndex.clear();
03028 
03029   std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
03030       meshfileit = this->Internal->MedFiles.begin();
03031   while(meshfileit != this->Internal->MedFiles.end())
03032     {
03033     vtkMedFile* file = meshfileit->second;
03034     meshfileit++;
03035     int numMeshes=file->GetNumberOfMesh();
03036     for(int meshIndex=0; meshIndex<numMeshes; meshIndex++)
03037       {
03038       vtkMedMesh* mesh=file->GetMesh(meshIndex);
03039       mesh->ClearMedSupports();
03040       }
03041 
03042     int numProf = file->GetNumberOfProfile();
03043     for (int prof = 0; prof<numProf; prof++)
03044       {
03045       vtkMedProfile* profile = file->GetProfile(prof);
03046       if (profile->GetIds()!=NULL)
03047         profile->GetIds()->Initialize();
03048       }
03049     }
03050 }
03051 
03052 void vtkMedReader::ClearMedFields()
03053 {
03054   this->Internal->FieldCache.clear();
03055   this->Internal->QuadOffsetKey.clear();
03056   this->Internal->QuadratureOffsetCache.clear();
03057 
03058   std::map<std::string, vtkSmartPointer<vtkMedFile> >::iterator
03059       fieldfileit = this->Internal->MedFiles.begin();
03060   while(fieldfileit != this->Internal->MedFiles.end())
03061     {
03062     vtkMedFile* file = fieldfileit->second;
03063     fieldfileit++;
03064 
03065     int numFields=file->GetNumberOfField();
03066     for(int ff=0; ff<numFields; ff++)
03067       {
03068       vtkMedField* field=file->GetField(ff);
03069       int nstep=field->GetNumberOfFieldStep();
03070       for(int sid=0; sid<nstep; sid++)
03071         {
03072         vtkMedFieldStep* step = field->GetFieldStep(sid);
03073         for(int id=0; id<step->GetNumberOfFieldOverEntity(); id++)
03074           {
03075           vtkMedFieldOverEntity * fieldOverEntity=step->GetFieldOverEntity(id);
03076           for(int pid = 0; pid < fieldOverEntity->GetNumberOfFieldOnProfile(); pid++)
03077             {
03078             vtkMedFieldOnProfile* fop = fieldOverEntity->GetFieldOnProfile(pid);
03079             if(fop->GetData() != NULL)
03080               fop->SetData(NULL);
03081             }
03082           }
03083         }
03084       }
03085     }
03086 }
03087 
03088 void vtkMedReader::ClearCaches(int when)
03089 {
03090   switch(when)
03091   {
03092     case Initialize:
03093       this->Internal->CurrentDataSet.clear();
03094       this->Internal->DataSetCache.clear();
03095       this->Internal->FieldCache.clear();
03096       this->Internal->UsedSupports.clear();
03097       this->Internal->QuadratureOffsetCache.clear();
03098       this->Internal->QuadOffsetKey.clear();
03099       //this->Internal->Med2VTKPointIndex.clear();
03100       break;
03101     case StartRequest:
03102       this->Internal->CurrentDataSet.clear();
03103       this->Internal->UsedSupports.clear();
03104       if(this->CacheStrategy==CacheNothing)
03105         {
03106         this->ClearMedSupports();
03107         this->ClearMedFields();
03108         }
03109       else if(this->CacheStrategy==CacheGeometry)
03110         {
03111         this->ClearMedFields();
03112         }
03113       break;
03114     case AfterCreateMedSupports:
03115       // TODO : clear the unused supports and associated cached datasets and fields
03116       break;
03117     case EndBuildVTKSupports:
03118       break;
03119     case EndRequest:
03120       if(this->CacheStrategy==CacheNothing)
03121         {
03122         this->ClearMedSupports();
03123         this->ClearMedFields();
03124         }
03125       else if(this->CacheStrategy==CacheGeometry && this->AnimationMode != Modes)
03126         {
03127         this->ClearMedFields();
03128         }
03129       break;
03130   }
03131 }
03132 
03133 void vtkMedReader::PrintSelf(ostream& os, vtkIndent indent)
03134 {
03135   PRINT_STRING(os, indent, FileName);
03136   PRINT_IVAR(os, indent, AnimationMode);
03137   PRINT_IVAR(os, indent, TimeIndexForIterations);
03138   PRINT_OBJECT(os, indent, PointFields);
03139   PRINT_OBJECT(os, indent, CellFields);
03140   PRINT_OBJECT(os, indent, QuadratureFields);
03141   PRINT_OBJECT(os, indent, ElnoFields);
03142   PRINT_OBJECT(os, indent, Groups);
03143   PRINT_OBJECT(os, indent, Entities);
03144   PRINT_IVAR(os, indent, CacheStrategy);
03145   this->Superclass::PrintSelf(os, indent);
03146 }