VTK
|
00001 /*========================================================================= 00002 00003 Program: Visualization Toolkit 00004 Module: vtkExodusReader.h 00005 00006 Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen 00007 All rights reserved. 00008 See Copyright.txt or http://www.kitware.com/Copyright.htm for details. 00009 00010 This software is distributed WITHOUT ANY WARRANTY; without even 00011 the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR 00012 PURPOSE. See the above copyright notice for more information. 00013 00014 =========================================================================*/ 00015 /*---------------------------------------------------------------------------- 00016 Copyright (c) Sandia Corporation 00017 See Copyright.txt or http://www.paraview.org/HTML/Copyright.html for details. 00018 ----------------------------------------------------------------------------*/ 00019 00039 #ifndef __vtkExodusReader_h 00040 #define __vtkExodusReader_h 00041 00042 #define ARRAY_TYPE_NAMES_IN_CXX_FILE 00043 00044 #include "vtkUnstructuredGridAlgorithm.h" 00045 00046 class vtkIntArray; 00047 class vtkFloatArray; 00048 class vtkDataArray; 00049 class vtkDataSet; 00050 class vtkPoints; 00051 class vtkExodusMetadata; 00052 class vtkExodusModel; 00053 class vtkExodusXMLParser; 00054 00055 00056 #include "vtkDSPFilterGroup.h" //for USE_EXO_DSP_FILTERS 00057 00058 00059 class VTK_HYBRID_EXPORT vtkExodusReader : public vtkUnstructuredGridAlgorithm 00060 { 00061 public: 00062 static vtkExodusReader *New(); 00063 vtkTypeMacro(vtkExodusReader,vtkUnstructuredGridAlgorithm); 00064 void PrintSelf(ostream& os, vtkIndent indent); 00065 00067 int CanReadFile(const char* fname); 00068 00070 00071 vtkSetStringMacro(FileName); 00072 vtkGetStringMacro(FileName); 00074 00076 00077 vtkSetStringMacro(XMLFileName); 00078 vtkGetStringMacro(XMLFileName); 00080 00082 00083 vtkSetMacro(TimeStep, int); 00084 vtkGetMacro(TimeStep, int); 00086 00088 00091 vtkSetMacro(GenerateBlockIdCellArray, int); 00092 vtkGetMacro(GenerateBlockIdCellArray, int); 00093 vtkBooleanMacro(GenerateBlockIdCellArray, int); 00094 const char *GetBlockIdArrayName() { return "BlockId"; } 00096 00097 00099 00102 vtkSetMacro(GenerateGlobalElementIdArray, int); 00103 vtkGetMacro(GenerateGlobalElementIdArray, int); 00104 vtkBooleanMacro(GenerateGlobalElementIdArray, int); 00106 //BTX 00107 enum { 00108 SEARCH_TYPE_ELEMENT=0, 00109 SEARCH_TYPE_NODE, 00110 SEARCH_TYPE_ELEMENT_THEN_NODE, 00111 SEARCH_TYPE_NODE_THEN_ELEMENT, 00112 ID_NOT_FOUND=-234121312 00113 }; 00114 //ETX 00115 static const char *GetGlobalElementIdArrayName() { return "GlobalElementId"; } 00116 static const char *GetPedigreeElementIdArrayName() { return "PedigreeElementId"; } 00117 static int GetGlobalElementID( vtkDataSet *data, int localID ); 00118 static int GetGlobalElementID ( vtkDataSet *data, int localID, 00119 int searchType ); 00120 00122 00126 vtkSetMacro(GenerateGlobalNodeIdArray, int); 00127 vtkGetMacro(GenerateGlobalNodeIdArray, int); 00128 vtkBooleanMacro(GenerateGlobalNodeIdArray, int); 00129 static const char *GetGlobalNodeIdArrayName() { return "GlobalNodeId"; } 00130 static const char *GetPedigreeNodeIdArrayName() { return "PedigreeNodeId"; } 00131 static int GetGlobalNodeID( vtkDataSet *data, int localID ); 00132 static int GetGlobalNodeID( vtkDataSet *data, int localID, 00133 int searchType ); 00135 00137 00141 vtkSetMacro(ApplyDisplacements, int); 00142 vtkGetMacro(ApplyDisplacements, int); 00143 vtkBooleanMacro(ApplyDisplacements, int); 00144 vtkSetMacro(DisplacementMagnitude, float); 00145 vtkGetMacro(DisplacementMagnitude, float); 00147 00149 00150 vtkGetStringMacro(Title); 00151 vtkGetMacro(Dimensionality, int); 00152 vtkGetMacro(NumberOfTimeSteps, int); 00153 int GetNumberOfElements() { return this->NumberOfUsedElements; } 00154 vtkGetMacro(NumberOfNodeSets, int); 00155 vtkGetMacro(NumberOfSideSets, int); 00156 vtkGetMacro(NumberOfBlocks, int); 00157 vtkGetVector2Macro(TimeStepRange, int); 00158 vtkSetVector2Macro(TimeStepRange, int); 00159 int GetNumberOfNodes() { return this->NumberOfUsedNodes; } 00160 int GetNumberOfElementsInBlock(int block_idx); 00161 int GetBlockId(int block_idx); 00162 virtual int GetTotalNumberOfNodes() { return this->NumberOfNodesInFile; } 00164 00165 00167 00172 int GetNumberOfPointArrays(); 00173 const char *GetPointArrayName(int index); 00174 int GetPointArrayID( const char *name ); 00175 int GetPointArrayNumberOfComponents(int index); 00176 void SetPointArrayStatus(int index, int flag); 00177 void SetPointArrayStatus(const char*, int flag); 00178 int GetPointArrayStatus(int index); 00179 int GetPointArrayStatus(const char*); 00181 00182 int GetNumberOfCellArrays(); 00183 const char *GetCellArrayName(int index); 00184 int GetCellArrayID( const char *name ); 00185 int GetCellArrayNumberOfComponents(int index); 00186 void SetCellArrayStatus(int index, int flag); 00187 void SetCellArrayStatus(const char*, int flag); 00188 int GetCellArrayStatus(int index); 00189 int GetCellArrayStatus(const char*); 00190 virtual int GetTotalNumberOfElements() 00191 { return this->NumberOfElementsInFile; } 00192 00194 00198 int GetNumberOfBlockArrays(); 00199 const char *GetBlockArrayName(int index); 00200 int GetBlockArrayID( const char *name ); 00201 void SetBlockArrayStatus(int index, int flag); 00202 void SetBlockArrayStatus(const char*, int flag); 00203 int GetBlockArrayStatus(int index); 00204 int GetBlockArrayStatus(const char*); 00206 00207 00209 00216 int GetNumberOfNodeSetArrays(){return this->GetNumberOfNodeSets();} 00217 int GetNodeSetArrayStatus(int index); 00218 int GetNodeSetArrayStatus(const char* name); 00219 void SetNodeSetArrayStatus(int index, int flag); 00220 void SetNodeSetArrayStatus(const char* name, int flag); 00221 const char *GetNodeSetArrayName(int index); 00223 00224 int GetNumberOfSideSetArrays(){return this->GetNumberOfSideSets();} 00225 int GetSideSetArrayStatus(int index); 00226 int GetSideSetArrayStatus(const char* name); 00227 void SetSideSetArrayStatus(int index, int flag); 00228 void SetSideSetArrayStatus(const char* name, int flag); 00229 const char *GetSideSetArrayName(int index); 00230 00232 00236 int GetNumberOfPartArrays(); 00237 const char *GetPartArrayName(int arrayIdx); 00238 int GetPartArrayID( const char *name ); 00239 const char *GetPartBlockInfo(int arrayIdx); 00240 void SetPartArrayStatus(int index, int flag); 00241 void SetPartArrayStatus(const char*, int flag); 00242 int GetPartArrayStatus(int index); 00243 int GetPartArrayStatus(const char*); 00245 00246 00248 00252 int GetNumberOfMaterialArrays(); 00253 const char *GetMaterialArrayName(int arrayIdx); 00254 int GetMaterialArrayID( const char *name ); 00255 void SetMaterialArrayStatus(int index, int flag); 00256 void SetMaterialArrayStatus(const char*, int flag); 00257 int GetMaterialArrayStatus(int index); 00258 int GetMaterialArrayStatus(const char*); 00260 00262 00266 int GetNumberOfAssemblyArrays(); 00267 const char *GetAssemblyArrayName(int arrayIdx); 00268 int GetAssemblyArrayID( const char *name ); 00269 void SetAssemblyArrayStatus(int index, int flag); 00270 void SetAssemblyArrayStatus(const char*, int flag); 00271 int GetAssemblyArrayStatus(int index); 00272 int GetAssemblyArrayStatus(const char*); 00274 00276 00283 int GetNumberOfHierarchyArrays(); 00284 const char *GetHierarchyArrayName(int arrayIdx); 00285 void SetHierarchyArrayStatus(int index, int flag); 00286 void SetHierarchyArrayStatus(const char*, int flag); 00287 int GetHierarchyArrayStatus(int index); 00288 int GetHierarchyArrayStatus(const char*); 00290 00292 00299 vtkGetMacro(HasModeShapes, int); 00300 vtkSetMacro(HasModeShapes, int); 00301 vtkBooleanMacro(HasModeShapes, int); 00303 00304 vtkGetMacro(DisplayType,int); 00305 virtual void SetDisplayType(int type); 00306 00312 vtkBooleanMacro(ExodusModelMetadata, int); 00313 vtkSetMacro(ExodusModelMetadata, int); 00314 vtkGetMacro(ExodusModelMetadata, int); 00315 00318 vtkExodusModel *GetExodusModel(){return this->ExodusModel;} 00319 00327 vtkSetMacro(PackExodusModelOntoOutput, int); 00328 vtkGetMacro(PackExodusModelOntoOutput, int); 00329 vtkBooleanMacro(PackExodusModelOntoOutput, int); 00330 00331 //BTX 00333 00334 enum ArrayType { 00335 CELL=0, 00336 POINT, 00337 BLOCK, 00338 PART, 00339 MATERIAL, 00340 ASSEMBLY, 00341 HIERARCHY, 00342 NUM_ARRAY_TYPES, 00343 UNKNOWN_TYPE 00344 }; 00346 //ETX 00347 00349 int IsValidVariable( const char *type, const char *name ); 00350 00351 //BTX 00353 00354 int GetNumberOfArrays( vtkExodusReader::ArrayType type ); 00355 const char *GetArrayName( vtkExodusReader::ArrayType type, int id ); 00357 //ETX 00358 00360 int GetVariableID ( const char *type, const char *name ); 00361 00362 void SetAllAssemblyArrayStatus( int status ); 00363 void SetAllBlockArrayStatus( int status ); 00364 void SetAllCellArrayStatus( int status ); 00365 void SetAllHierarchyArrayStatus( int status ); 00366 void SetAllMaterialArrayStatus( int status ); 00367 void SetAllPartArrayStatus( int status ); 00368 void SetAllPointArrayStatus( int status ); 00369 //BTX 00370 void SetAllArrayStatus ( vtkExodusReader::ArrayType type, int flag ); 00371 void SetArrayStatus ( vtkExodusReader::ArrayType type, const char *name, 00372 int flag ); 00373 //ETX 00374 void SetArrayStatus ( const char *type, const char *name, int flag ) 00375 { 00376 this->SetArrayStatus( this->GetArrayTypeID(type), name, flag ); 00377 } 00378 //BTX 00379 int GetArrayStatus ( vtkExodusReader::ArrayType type, const char *name ); 00380 //ETX 00381 int GetArrayStatus ( const char *type, const char *name ) 00382 { 00383 return this->GetArrayStatus( this->GetArrayTypeID( type ), name ); 00384 } 00385 00386 // Helper functions 00387 static int StringsEqual(const char* s1, char* s2); 00388 static void StringUppercase(const char* str, char* upperstr); 00389 static char *StrDupWithNew(const char *s); 00390 00391 // time series query functions 00392 int GetTimeSeriesData( int ID, const char *vName, const char *vType, 00393 vtkFloatArray *result ); 00394 00395 00396 //begin USE_EXO_DSP_FILTERS 00397 int GetNumberOfVariableArrays(); 00398 const char *GetVariableArrayName(int a_which); 00399 void EnableDSPFiltering(); 00400 void AddFilter(vtkDSPFilterDefinition *a_filter); 00401 void StartAddingFilter(); 00402 void AddFilterInputVar(char *name); 00403 void AddFilterOutputVar(char *name); 00404 void AddFilterNumeratorWeight(double weight); 00405 void AddFilterForwardNumeratorWeight(double weight); 00406 void AddFilterDenominatorWeight(double weight); 00407 void FinishAddingFilter(); 00408 void RemoveFilter(char *a_outputVariableName); 00409 void GetDSPOutputArrays(int exoid, vtkUnstructuredGrid* output); 00410 //BTX 00411 vtkExodusReader::ArrayType GetArrayTypeID( const char *type ); 00412 00413 #ifdef ARRAY_TYPE_NAMES_IN_CXX_FILE 00414 static const char *GetArrayTypeName( vtkExodusReader::ArrayType type ); 00415 #else 00416 static const char *ArrayTypeNames[NUM_ARRAY_TYPES]; 00417 00418 static const char *GetArrayTypeName( vtkExodusReader::ArrayType type ) 00419 { 00420 return ArrayTypeNames[type]; 00421 } 00422 #endif 00423 //ETX 00424 00425 vtkDSPFilterDefinition *AddingFilter; 00426 int DSPFilteringIsEnabled; 00427 vtkDSPFilterGroup **DSPFilters; 00428 //end USE_EXO_DSP_FILTERS 00429 00430 00431 00432 protected: 00433 vtkExodusReader(); 00434 ~vtkExodusReader(); 00435 00436 void NewExodusModel(); 00437 00438 void ReadGeometry(int exoid, vtkUnstructuredGrid* output); 00439 void ReadCells(int exoid, vtkUnstructuredGrid* output); 00440 void ReadPoints(int exoid, vtkUnstructuredGrid* output); 00441 void ReadArrays(int exoid, vtkUnstructuredGrid* output); 00442 void ReadNodeAndSideSets(int exoid, vtkUnstructuredGrid* output); 00443 vtkDataArray *ReadPointArray(int exoid, int varIndex); 00444 vtkDataArray *ReadPointVector(int handle, int varIndex, int dim); 00445 vtkDataArray *ReadCellArray(int exoid, int varIndex); 00446 vtkDataArray *ReadCellVector(int handle, int varIndex, int dim); 00447 void ReadNodeSetMetadata(); 00448 void ReadSideSetMetadata(); 00449 00450 // helper for finding IDs 00451 static int GetIDHelper ( const char *arrayName, vtkDataSet *data, int localID, 00452 int searchType ); 00453 static int GetGlobalID( const char *arrayName, vtkDataSet *data, int localID, 00454 int searchType ); 00455 00456 // This method is a helper for determining the 00457 // number of additional cell scalar field 00458 // values needed to 'pad' for node and side sets 00459 int GetExtraCellCountForNodeSideSets(); 00460 00461 // This method generates arrays like blockid, global nodeid 00462 // and global element id 00463 void GenerateExtraArrays(vtkUnstructuredGrid* output); 00464 00465 // Parameters for controlling what is read in. 00466 char *FileName; 00467 char *XMLFileName; 00468 int TimeStep; 00469 int ActualTimeStep; 00470 double TimeValue; 00471 int GenerateBlockIdCellArray; 00472 int GenerateGlobalElementIdArray; 00473 int GenerateGlobalNodeIdArray; 00474 int ApplyDisplacements; 00475 double DisplacementMagnitude; 00476 00477 // Information specific for exodus files. 00478 vtkSetStringMacro(Title); 00479 char *Title; 00480 int Dimensionality; 00481 int NumberOfNodeSets; 00482 int NumberOfSideSets; 00483 int NumberOfBlocks; 00484 int NumberOfUsedNodes; 00485 int NumberOfNodesInFile; 00486 int NumberOfUsedElements; 00487 int NumberOfElementsInFile; 00488 int NumberOfTimeSteps; 00489 int ExodusCPUWordSize; 00490 int ExodusIOWordSize; 00491 float ExodusVersion; 00492 vtkIntArray *CellVarTruthTable; 00493 00494 //1=display Block names 00495 //2=display Part names 00496 //3=display Material names 00497 int DisplayType; 00498 00499 //Parser that understands the xml part and material file 00500 vtkExodusXMLParser *Parser; 00501 00502 // **KEN** By VTK convention, metaData should be Metadata. 00503 00505 // Scalar Array and Block Info 00507 vtkExodusMetadata *MetaData; 00508 00509 00511 00512 int CurrentHandle; 00513 char* CurrentFileName; 00514 char* CurrentXMLFileName; 00515 vtkSetStringMacro(CurrentFileName); 00516 vtkSetStringMacro(CurrentXMLFileName); 00518 00519 // Open the exodus file, and set some basic information 00520 int OpenCurrentFile(); 00521 00522 // Close the exodus file 00523 void CloseCurrentFile(); 00524 00525 00527 int TimeStepRange[2]; 00528 00529 // DataCache: this object keeps the points and cells 00530 // around so they don't need to be re-read when the 00531 // timestep changes or an scalar array is switched 00532 vtkUnstructuredGrid *DataCache; 00533 00534 // Should I re-read in the geometry and topology of the dataset 00535 int RemakeDataCacheFlag; 00536 00537 // vtkExodusModel needs to count changes in geometry, so it knows 00538 // if geometry has changed since it last updated model data. 00539 00540 int NewGeometryCount; 00541 00542 // PointMap keeps track of which points are actually 00543 // used by the cells that are read in (blocks) 00544 vtkIntArray *PointMap; 00545 vtkIntArray *ReversePointMap; 00546 void SetUpPointMap(int num_points); 00547 int GetPointMapIndex(int point_id); 00548 00549 // Global element ID cache 00550 int *GlobalElementIdCache; 00551 void SetGlobalElementIdCache(int *list); 00552 00553 // Time query function. Called by ExecuteInformation(). 00554 // Fills the TimestepValues array. 00555 void GetAllTimes(vtkInformationVector *); 00556 00557 int HasModeShapes; 00558 00559 vtkExodusModel *ExodusModel; 00560 int PackExodusModelOntoOutput; 00561 int ExodusModelMetadata; 00562 00563 double *TimeSteps; 00564 00565 int RequestInformation( 00566 vtkInformation *, vtkInformationVector **, vtkInformationVector *); 00567 int RequestData( 00568 vtkInformation *, vtkInformationVector **, vtkInformationVector *); 00569 00570 // Used to determine current progress. 00571 double ProgressOffset; 00572 double ProgressScale; 00573 00574 private: 00575 vtkExodusReader(const vtkExodusReader&); // Not implemented 00576 void operator=(const vtkExodusReader&); // Not implemented 00577 00578 void AddDisplacements(vtkUnstructuredGrid* output); 00579 void RemoveBeginningAndTrailingSpaces(char **names, int len); 00580 00581 void FixMetadataTruthTable(int *table, int len); 00582 }; 00583 00584 #endif