aboutsummaryrefslogtreecommitdiff
path: root/src/AMRPlus/AMRgridreaderPlus.C
blob: f2e64374fd07186c61385699308af94b24b18f31 (plain)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
#include "AMRgridreaderPlus.h"
#include <stdio.h>
#include <stdlib.h>
#include <iostream.h>





AMRgridPlus *AMRgridreaderPlus::getGridInfo(AMRgridPlus &g,int index){
  g.dataveclen=1;
  if(file.seek(index)<index)
    return 0; // fail if index past end
  IObase::DataType dt;
  file.readInfo(dt,g.rank,g.dims);
  g.datatype = dt;
  g.nbytes = IObase::nBytes(dt,g.rank,g.dims);
  // find the deepest level (finest time resolution)
  // Attrib Names?
  IObase::DataType atype;
  int length;
  int attrnum=file.readAttributeInfo("level",atype,length);
  if(attrnum>=0){
    int lev; // should be Int level
    file.readAttribute(attrnum,&lev);
    if(lev>g.maxlevel) g.maxlevel=lev;
    g.level=lev;
  }
  attrnum=file.readAttributeInfo("time",atype,length);
  if(attrnum>=0){
    file.readAttribute(attrnum,&(g.time));
  }attrnum=file.readAttributeInfo("time_refinement",atype,length);
  if(attrnum>=0){
    file.readAttribute(attrnum,&(g.timerefinement));
  }
  attrnum=file.readAttributeInfo("timestep",atype,length);
  if(attrnum>=0){
    file.readAttribute(attrnum,&(g.timestep));
  }
  attrnum=file.readAttributeInfo("origin",atype,length);
  if(attrnum>=0)
    file.readAttribute(attrnum,(g.origin));
  attrnum=file.readAttributeInfo("delta",atype,length);
  if(attrnum>=0)
    file.readAttribute(attrnum,(g.delta));
  attrnum=file.readAttributeInfo("persistence",atype,length);
  if(attrnum>=0){
    file.readAttribute(attrnum,&(g.persistence));
    g.maxtime = g.timestep + g.persistence;
  }
  
  attrnum=file.readAttributeInfo("time_refinement",atype,length);
  if(attrnum>=0){
    file.readAttribute(attrnum,&(g.timeref));
  }
  attrnum=file.readAttributeInfo("spatial_refinement",atype,length);
  if(attrnum>=0){
    file.readAttribute(attrnum,&(g.spaceref));
  }
  attrnum=file.readAttributeInfo("grid_placement_refinement",atype,length);
  if(attrnum>=0){
    file.readAttribute(attrnum,&(g.placeref));
  }
  
  attrnum=file.readAttributeInfo("range",atype,length);
  if(attrnum>=0){
    double range[2];
    file.readAttribute(attrnum,&range);
    g.scalarmin=range[0]; g.scalarmax=range[1];
  }
  else {
      g.data=malloc(g.nbytes);
      file.read(g.data);
      double range[2];
      
      switch (dt){
	  case IObase::Float32:{
	    long nscal=g.nbytes/sizeof(float);
	    float *scalars=(float *)g.data;
    
	    range[0]=range[1]=scalars[0];
	    for (int ii=1;ii<nscal;ii++){
		float scal=scalars[ii];
		range[0]=(range[0]<scal)?range[0]:scal;
		range[1]=(range[1]>scal)?range[1]:scal;
	    }
	  } break;
	  case IObase::Float64:{
	    long nscal=g.nbytes/sizeof(double);
	    double *scalars=(double *)g.data;
    
	    range[0]=range[1]=scalars[0];
	    for (int ii=1;ii<nscal;ii++){
		double scal=scalars[ii];
		range[0]=(range[0]<scal)?range[0]:scal;
		range[1]=(range[1]>scal)?range[1]:scal;
	    }
	  } break;
	  default:
	    cout<<"Data type is not float or double."<<endl;
	    range[0]=0;range[1]=1;
	}
	g.scalarmin=range[0];
	g.scalarmax=range[1];
	free(g.data);
    }
  
  g.data=0;
  return &g;
} // done 
  
AMRgridPlus *AMRgridreaderPlus::getGridData(AMRgridPlus &g,int index){
  IObase::DataType atype;
  // if(data) free(data); data=0; // make certain it is empty first
  g.data = malloc(g.nbytes);
  file.seek(index);
  file.readInfo(atype,g.rank,g.dims);
  g.datatype=atype;
  file.read(g.data);
  return &g;
}