ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
ERF_NCWpsFile.H File Reference
#include <sstream>
#include <string>
#include <atomic>
#include "AMReX_FArrayBox.H"
#include "AMReX_IArrayBox.H"
#include "ERF_EpochTime.H"
#include "ERF_NCInterface.H"
Include dependency graph for ERF_NCWpsFile.H:
This graph shows which files directly or indirectly include this file:

Go to the source code of this file.

Classes

struct  NDArray< DataType >
 

Typedefs

using PlaneVector = amrex::Vector< amrex::FArrayBox >
 

Enumerations

enum class  NC_Data_Dims_Type {
  Time_BT_SN_WE , Time_SN_WE , Time_BT , Time ,
  Time_BdyWidth_BT_SN , Time_BdyWidth_BT_WE , Time_BdyWidth_SN , Time_BdyWidth_WE
}
 

Functions

int BuildFABsFromWRFBdyFile (const std::string &fname, amrex::Vector< amrex::Vector< amrex::FArrayBox >> &bdy_data_xlo, amrex::Vector< amrex::Vector< amrex::FArrayBox >> &bdy_data_xhi, amrex::Vector< amrex::Vector< amrex::FArrayBox >> &bdy_data_ylo, amrex::Vector< amrex::Vector< amrex::FArrayBox >> &bdy_data_yhi)
 
template<typename DType >
void ReadTimeSliceFromNetCDFFile (const std::string &fname, const int tidx, amrex::Vector< std::string > names, amrex::Vector< NDArray< DType > > &arrays, amrex::Vector< int > &success)
 
template<typename DType >
void ReadNetCDFFile (const std::string &fname, amrex::Vector< std::string > names, amrex::Vector< NDArray< DType > > &arrays, amrex::Vector< int > &success)
 
template<class FAB , typename DType >
void fill_fab_from_arrays (int iv, amrex::Vector< NDArray< float >> &nc_arrays, const std::string &var_name, NC_Data_Dims_Type &NC_dim_type, FAB &temp)
 
template<class FAB , typename DType >
void BuildFABsFromNetCDFFile (const amrex::Box &domain, const std::string &fname, amrex::Vector< std::string > nc_var_names, amrex::Vector< enum NC_Data_Dims_Type > NC_dim_types, amrex::Vector< FAB * > fab_vars, amrex::Vector< int > &success)
 

Typedef Documentation

◆ PlaneVector

using PlaneVector = amrex::Vector<amrex::FArrayBox>

Enumeration Type Documentation

◆ NC_Data_Dims_Type

enum NC_Data_Dims_Type
strong
Enumerator
Time_BT_SN_WE 
Time_SN_WE 
Time_BT 
Time 
Time_BdyWidth_BT_SN 
Time_BdyWidth_BT_WE 
Time_BdyWidth_SN 
Time_BdyWidth_WE 

Function Documentation

◆ BuildFABsFromNetCDFFile()

template<class FAB , typename DType >
void BuildFABsFromNetCDFFile ( const amrex::Box &  domain,
const std::string &  fname,
amrex::Vector< std::string >  nc_var_names,
amrex::Vector< enum NC_Data_Dims_Type NC_dim_types,
amrex::Vector< FAB * >  fab_vars,
amrex::Vector< int > &  success 
)

Function to read NetCDF variables and fill the corresponding Array4's

Parameters
fnameName of the NetCDF file to be read
nc_var_namesVariable names in the NetCDF file
NC_dim_typesNetCDF data dimension types
fab_varsFab data we are to fill
344 {
345  int ioproc = amrex::ParallelDescriptor::IOProcessorNumber(); // I/O rank
346 
347  amrex::Vector<NDArray<float>> nc_arrays(nc_var_names.size());
348 
349  if (amrex::ParallelDescriptor::IOProcessor())
350  {
351  ReadNetCDFFile(fname, nc_var_names, nc_arrays, success);
352  }
353 
354  amrex::ParallelDescriptor::Bcast(success.dataPtr(), success.size(), ioproc);
355 
356  for (int iv = 0; iv < nc_var_names.size(); iv++)
357  {
358  if (success[iv] == 1) {
359  FAB tmp;
360  if (amrex::ParallelDescriptor::IOProcessor()) {
361  fill_fab_from_arrays<FAB,DType>(iv, nc_arrays, nc_var_names[iv], NC_dim_types[iv], tmp);
362  }
363 
364  int ncomp = tmp.nComp();
365  amrex::Box box = tmp.box();
366 
367  amrex::ParallelDescriptor::Bcast(&box, 1, ioproc);
368  amrex::ParallelDescriptor::Bcast(&ncomp, 1, ioproc);
369 
370  if (!amrex::ParallelDescriptor::IOProcessor()) {
371 #ifdef AMREX_USE_GPU
372  tmp.resize(box,ncomp,amrex::The_Pinned_Arena());
373 #else
374  tmp.resize(box,ncomp);
375 #endif
376  }
377 
378  amrex::ParallelDescriptor::Bcast(tmp.dataPtr(), tmp.size(), ioproc);
379 
380  // Shift box by the domain lower corner
381  amrex::Box fab_bx = tmp.box();
382  amrex::Dim3 dom_lb = lbound(domain);
383  fab_bx += amrex::IntVect(dom_lb.x,dom_lb.y,dom_lb.z);
384  // fab_vars points to data on device
385  fab_vars[iv]->resize(fab_bx,1);
386 #ifdef AMREX_USE_GPU
387  amrex::Gpu::copy(amrex::Gpu::hostToDevice,
388  tmp.dataPtr(), tmp.dataPtr() + tmp.size(),
389  fab_vars[iv]->dataPtr());
390 #else
391  // Provided by BaseFab inheritance through FArrayBox
392  fab_vars[iv]->copy(tmp,tmp.box(),0,fab_bx,0,1);
393 #endif
394  } // success
395  } // iv
396 }
void ReadNetCDFFile(const std::string &fname, amrex::Vector< std::string > names, amrex::Vector< NDArray< DType > > &arrays, amrex::Vector< int > &success)
Definition: ERF_NCWpsFile.H:179
Here is the call graph for this function:

◆ BuildFABsFromWRFBdyFile()

int BuildFABsFromWRFBdyFile ( const std::string &  fname,
amrex::Vector< amrex::Vector< amrex::FArrayBox >> &  bdy_data_xlo,
amrex::Vector< amrex::Vector< amrex::FArrayBox >> &  bdy_data_xhi,
amrex::Vector< amrex::Vector< amrex::FArrayBox >> &  bdy_data_ylo,
amrex::Vector< amrex::Vector< amrex::FArrayBox >> &  bdy_data_yhi 
)

◆ fill_fab_from_arrays()

template<class FAB , typename DType >
void fill_fab_from_arrays ( int  iv,
amrex::Vector< NDArray< float >> &  nc_arrays,
const std::string &  var_name,
NC_Data_Dims_Type NC_dim_type,
FAB &  temp 
)

Helper function for reading data from NetCDF file into a provided FAB.

Parameters
ivIndex for which variable we are going to fill
nc_arraysArrays of data from NetCDF file
var_nameVariable name
NC_dim_typeDimension type for the variable as stored in the NetCDF file
tempFAB where we store the variable data from the NetCDF Arrays
270 {
271  int ns1, ns2, ns3; // bottom_top, south_north, west_east (these can be staggered or unstaggered)
272  if (NC_dim_type == NC_Data_Dims_Type::Time_BT) {
273  ns1 = nc_arrays[iv].get_vshape()[1];
274  ns2 = 1;
275  ns3 = 1;
276  // amrex::Print() << "TYPE BT " << ns1 << std::endl;
277  } else if (NC_dim_type == NC_Data_Dims_Type::Time_SN_WE) {
278  ns1 = 1;
279  ns2 = nc_arrays[iv].get_vshape()[1];
280  ns3 = nc_arrays[iv].get_vshape()[2];
281  // amrex::Print() << "TYPE SN WE " << ns2 << " " << ns3 << std::endl;
282  } else if (NC_dim_type == NC_Data_Dims_Type::Time_BT_SN_WE) {
283  ns1 = nc_arrays[iv].get_vshape()[1];
284  ns2 = nc_arrays[iv].get_vshape()[2];
285  ns3 = nc_arrays[iv].get_vshape()[3];
286  // amrex::Print() << "TYPE BT SN WE " << ns1 << " " << ns2 << " " << ns3 << std::endl;
287  } else {
288  amrex::Abort("Dont know this NC_Data_Dims_Type");
289  }
290 
291  // TODO: The box will only start at (0,0,0) at level 0 -- we need to generalize this
292  amrex::Box my_box(amrex::IntVect(0,0,0), amrex::IntVect(ns3-1,ns2-1,ns1-1));
293 
294  if (var_name == "PH" || var_name == "PHB") {
295  my_box.setType(amrex::IndexType(amrex::IntVect(0,0,1)));
296  }
297  else if (var_name == "U" || var_name == "UU" || var_name == "MAPFAC_U") {
298  my_box.setType(amrex::IndexType(amrex::IntVect(1,0,0)));
299  }
300  else if (var_name == "V" || var_name == "VV" || var_name == "MAPFAC_V") {
301  my_box.setType(amrex::IndexType(amrex::IntVect(0,1,0)));
302  }
303  else if (var_name == "W" || var_name == "WW") {
304  my_box.setType(amrex::IndexType(amrex::IntVect(0,0,1)));
305  }
306 
307  amrex::Arena* Arena_Used = amrex::The_Arena();
308 #ifdef AMREX_USE_GPU
309  // Make sure temp lives on CPU since nc_arrays lives on CPU only
310  Arena_Used = amrex::The_Pinned_Arena();
311 #endif
312  temp.resize(my_box,1, Arena_Used);
313  amrex::Array4<DType> fab_arr = temp.array();
314 
315  int ioff = temp.box().smallEnd()[0];
316  int joff = temp.box().smallEnd()[1];
317 
318  auto num_pts = my_box.numPts();
319 
320  for (int n(0); n < num_pts; ++n) {
321  int k = n / (ns2*ns3);
322  int j = (n - k*(ns2*ns3)) / ns3 + joff;
323  int i = n - k*(ns2*ns3) - (j-joff) * ns3 + ioff;
324  fab_arr(i,j,k,0) = static_cast<DType>(*(nc_arrays[iv].get_data()+n));
325  }
326 }

◆ ReadNetCDFFile()

template<typename DType >
void ReadNetCDFFile ( const std::string &  fname,
amrex::Vector< std::string >  names,
amrex::Vector< NDArray< DType > > &  arrays,
amrex::Vector< int > &  success 
)
181 {
182  AMREX_ASSERT(arrays.size() == names.size());
183 
184  if (amrex::ParallelDescriptor::IOProcessor())
185  {
186  auto ncf = ncutils::NCFile::open(fname, NC_CLOBBER | NC_NETCDF4);
187 
188  /*
189  // get the dimension information
190  int Time = static_cast<int>(ncf.dim("Time").len());
191  int DateStrLen = static_cast<int>(ncf.dim("DateStrLen").len());
192  int west_east = static_cast<int>(ncf.dim("west_east").len());
193  int south_north = static_cast<int>(ncf.dim("south_north").len());
194  int bottom_top = static_cast<int>(ncf.dim("bottom_top").len());
195  int bottom_top_stag = static_cast<int>(ncf.dim("bottom_top_stag").len());
196  int west_east_stag = static_cast<int>(ncf.dim("west_east_stag").len());
197  int south_north_stag = static_cast<int>(ncf.dim("south_north_stag").len());
198  int bdy_width = static_cast<int>(ncf.dim("bdy_width").len());
199  */
200 
201  // amrex::Print() << "Reading the dimensions from the netcdf file " << "\n";
202 
203  for (auto n=0; n<arrays.size(); ++n)
204  {
205  std::string vname_to_write = names[n];
206  std::string vname_to_read = names[n];
207 
208  if (vname_to_read.substr(0,2) == "R_") {
209  vname_to_read = names[n+4]; // This allows us to read "T" instead -- we will over-write this later
210  }
211 
212  success[n] = ncf.has_var(vname_to_read);
213  if (success[n] == 0) {
214  amrex::Print() << " Skipping " << vname_to_read << std::endl;
215  } else {
216  amrex::Print() << " Reading " << vname_to_read << std::endl;
217  }
218 
219  if (success[n] == 1) {
220 
221  std::vector<std::string> dimnames = ncf.var(vname_to_read).dimnames();
222  AMREX_ALWAYS_ASSERT(dimnames[0] == "Time" || dimnames[0] == "time");
223 
224  std::vector<size_t> shape = ncf.var(vname_to_read).shape();
225  arrays[n] = NDArray<DType>(vname_to_read,shape);
226  DType* dataPtr = arrays[n].get_data();
227 
228  std::vector<size_t> start(shape.size(), 0);
229 
230 #if 0
231  auto numPts = arrays[n].ndim();
232  amrex::Print() << "NetCDF Variable name = " << vname_to_read << std::endl;
233  amrex::Print() << "numPts read from NetCDF file/var = " << numPts << std::endl;
234  amrex::Print() << "Dims in var = " << ncf.var(vname_to_read).ndim() << std::endl;
235  amrex::Print() << "Dim names = (";
236  for (auto &dim:dimnames)
237  amrex::Print() << dim << ", " ;
238  amrex::Print() << ")" << std::endl;
239  amrex::Print() << "Dims of the variable = (";
240  for (auto &dim:shape)
241  amrex::Print() << dim << ", " ;
242  amrex::Print() << ")" << std::endl;
243 #endif
244 
245  ncf.var(vname_to_read).get(dataPtr, start, shape);
246 
247  } // has_var
248  }
249  ncf.close();
250  }
251 }
static NCFile open(const std::string &name, const int cmode=NC_NOWRITE)
Definition: ERF_NCInterface.cpp:707
Definition: ERF_NCWpsFile.H:51

Referenced by BuildFABsFromNetCDFFile().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ ReadTimeSliceFromNetCDFFile()

template<typename DType >
void ReadTimeSliceFromNetCDFFile ( const std::string &  fname,
const int  tidx,
amrex::Vector< std::string >  names,
amrex::Vector< NDArray< DType > > &  arrays,
amrex::Vector< int > &  success 
)
137 {
138  amrex::Print() << "Reading time slice " << tidx << " from " << fname << std::endl;
139  AMREX_ASSERT(arrays.size() == names.size());
140 
141  if (amrex::ParallelDescriptor::IOProcessor())
142  {
143  auto ncf = ncutils::NCFile::open(fname, NC_CLOBBER | NC_NETCDF4);
144 
145  int ntimes = ncf.dim("Time").len();
146  AMREX_ALWAYS_ASSERT((tidx >= 0) && (tidx < ntimes));
147 
148  for (auto n=0; n<arrays.size(); ++n)
149  {
150  std::string vname_to_write = names[n];
151  std::string vname_to_read = names[n];
152  if (vname_to_read.substr(0,2) == "R_") {
153  vname_to_read = names[n+4]; // This allows us to read "T" instead -- we will over-write this later
154  }
155 
156  success[n] = ncf.has_var(vname_to_read);
157 
158  if (success[n] == 1)
159  {
160  std::vector<std::string> dimnames = ncf.var(vname_to_read).dimnames();
161  AMREX_ALWAYS_ASSERT(dimnames[0] == "Time");
162 
163  std::vector<size_t> count = ncf.var(vname_to_read).shape();
164  std::vector<size_t> start(count.size(), 0);
165  start[0] = tidx;
166  count[0] = 1;
167 
168  arrays[n] = NDArray<DType>(vname_to_read, count);
169  DType* dataPtr = arrays[n].get_data();
170 
171  ncf.var(vname_to_read).get(dataPtr, start, count);
172  } // has_var
173  }
174  ncf.close();
175  }
176 }
Here is the call graph for this function: