1 #ifndef ERF_NCWPSFILE_H_
2 #define ERF_NCWPSFILE_H_
8 #include "AMReX_FArrayBox.H"
9 #include "AMReX_IArrayBox.H"
48 template <
typename DataType>
51 using DType =
typename std::remove_const<DataType>::type;
54 explicit NDArray (
const std::string vname,
const std::vector<size_t>& vshape)
68 ref_counted.fetch_add(1, std::memory_order_relaxed);
77 ref_counted.fetch_add(1, std::memory_order_relaxed);
83 ref_counted.fetch_sub(1, std::memory_order_acq_rel);
89 ref_counted.fetch_add(1, std::memory_order_relaxed);
106 int isize =
static_cast<int>(
shape.size());
107 for (
auto i=0; i<isize; ++i)
num *=
shape[i];
125 amrex::Vector<amrex::Vector<amrex::FArrayBox>>& bdy_data_xlo,
126 amrex::Vector<amrex::Vector<amrex::FArrayBox>>& bdy_data_xhi,
127 amrex::Vector<amrex::Vector<amrex::FArrayBox>>& bdy_data_ylo,
128 amrex::Vector<amrex::Vector<amrex::FArrayBox>>& bdy_data_yhi);
130 template<
typename DType>
133 amrex::Vector<std::string> names,
135 amrex::Vector<int>& success)
137 amrex::Print() <<
"Reading time slice " << tidx <<
" from " << fname << std::endl;
138 AMREX_ASSERT(arrays.size() == names.size());
140 if (amrex::ParallelDescriptor::IOProcessor())
144 int ntimes = ncf.dim(
"Time").len();
145 AMREX_ALWAYS_ASSERT((tidx >= 0) && (tidx < ntimes));
147 for (
auto n=0; n<arrays.size(); ++n)
149 std::string vname_to_write = names[n];
150 std::string vname_to_read = names[n];
151 if (vname_to_read.substr(0,2) ==
"R_") {
152 vname_to_read = names[n+4];
155 success[n] = ncf.has_var(vname_to_read);
159 std::vector<std::string> dimnames = ncf.var(vname_to_read).dimnames();
160 AMREX_ALWAYS_ASSERT(dimnames[0] ==
"Time");
162 std::vector<size_t> count = ncf.var(vname_to_read).shape();
163 std::vector<size_t> start(count.size(), 0);
168 DType* dataPtr = arrays[n].get_data();
170 ncf.var(vname_to_read).get(dataPtr, start, count);
177 template<
typename DType>
179 amrex::Vector<
NDArray<DType> >& arrays, amrex::Vector<int>& success)
181 AMREX_ASSERT(arrays.size() == names.size());
183 if (amrex::ParallelDescriptor::IOProcessor())
202 for (
auto n=0; n<arrays.size(); ++n)
204 std::string vname_to_write = names[n];
205 std::string vname_to_read = names[n];
207 if (vname_to_read.substr(0,2) ==
"R_") {
208 vname_to_read = names[n+4];
211 success[n] = ncf.has_var(vname_to_read);
212 if (success[n] == 0) {
213 amrex::Print() <<
" Skipping " << vname_to_read << std::endl;
215 amrex::Print() <<
" Reading " << vname_to_read << std::endl;
218 if (success[n] == 1) {
220 std::vector<std::string> dimnames = ncf.var(vname_to_read).dimnames();
221 AMREX_ALWAYS_ASSERT(dimnames[0] ==
"Time" || dimnames[0] ==
"time");
223 std::vector<size_t> shape = ncf.var(vname_to_read).shape();
225 DType* dataPtr = arrays[n].get_data();
227 std::vector<size_t> start(shape.size(), 0);
230 auto numPts = arrays[n].ndim();
231 amrex::Print() <<
"NetCDF Variable name = " << vname_to_read << std::endl;
232 amrex::Print() <<
"numPts read from NetCDF file/var = " << numPts << std::endl;
233 amrex::Print() <<
"Dims in var = " << ncf.var(vname_to_read).ndim() << std::endl;
234 amrex::Print() <<
"Dim names = (";
235 for (
auto &dim:dimnames)
236 amrex::Print() << dim <<
", " ;
237 amrex::Print() <<
")" << std::endl;
238 amrex::Print() <<
"Dims of the variable = (";
239 for (
auto &dim:shape)
240 amrex::Print() << dim <<
", " ;
241 amrex::Print() <<
")" << std::endl;
244 ncf.var(vname_to_read).get(dataPtr, start, shape);
262 template<
class FAB,
typename DType>
266 const std::string& var_name,
272 ns1 = nc_arrays[iv].get_vshape()[1];
278 ns2 = nc_arrays[iv].get_vshape()[1];
279 ns3 = nc_arrays[iv].get_vshape()[2];
282 ns1 = nc_arrays[iv].get_vshape()[1];
283 ns2 = nc_arrays[iv].get_vshape()[2];
284 ns3 = nc_arrays[iv].get_vshape()[3];
287 amrex::Abort(
"Dont know this NC_Data_Dims_Type");
291 amrex::Box my_box(amrex::IntVect(0,0,0), amrex::IntVect(ns3-1,ns2-1,ns1-1));
293 if (var_name ==
"PH" || var_name ==
"PHB") {
294 my_box.setType(amrex::IndexType(amrex::IntVect(0,0,1)));
296 else if (var_name ==
"U" || var_name ==
"UU" || var_name ==
"MAPFAC_U") {
297 my_box.setType(amrex::IndexType(amrex::IntVect(1,0,0)));
299 else if (var_name ==
"V" || var_name ==
"VV" || var_name ==
"MAPFAC_V") {
300 my_box.setType(amrex::IndexType(amrex::IntVect(0,1,0)));
302 else if (var_name ==
"W" || var_name ==
"WW") {
303 my_box.setType(amrex::IndexType(amrex::IntVect(0,0,1)));
306 amrex::Arena* Arena_Used = amrex::The_Arena();
309 Arena_Used = amrex::The_Pinned_Arena();
311 temp.resize(my_box,1, Arena_Used);
312 amrex::Array4<DType> fab_arr = temp.array();
314 int ioff = temp.box().smallEnd()[0];
315 int joff = temp.box().smallEnd()[1];
317 auto num_pts = my_box.numPts();
319 for (
int n(0); n < num_pts; ++n) {
320 int k = n / (ns2*ns3);
321 int j = (n - k*(ns2*ns3)) / ns3 + joff;
322 int i = n - k*(ns2*ns3) - (j-joff) * ns3 + ioff;
323 fab_arr(i,j,k,0) =
static_cast<DType
>(*(nc_arrays[iv].get_data()+n));
335 template<
class FAB,
typename DType>
338 const std::string &fname,
339 amrex::Vector<std::string> nc_var_names,
340 amrex::Vector<enum NC_Data_Dims_Type> NC_dim_types,
341 amrex::Vector<FAB*> fab_vars,
342 amrex::Vector<int>& success)
344 int ioproc = amrex::ParallelDescriptor::IOProcessorNumber();
346 amrex::Vector<NDArray<float>> nc_arrays(nc_var_names.size());
348 if (amrex::ParallelDescriptor::IOProcessor())
353 amrex::ParallelDescriptor::Bcast(success.dataPtr(), success.size(), ioproc);
355 for (
int iv = 0; iv < nc_var_names.size(); iv++)
357 if (success[iv] == 1) {
359 if (amrex::ParallelDescriptor::IOProcessor()) {
360 fill_fab_from_arrays<FAB,DType>(iv, nc_arrays, nc_var_names[iv], NC_dim_types[iv], tmp);
363 int ncomp = tmp.nComp();
364 amrex::Box box = tmp.box();
366 amrex::ParallelDescriptor::Bcast(&box, 1, ioproc);
367 amrex::ParallelDescriptor::Bcast(&ncomp, 1, ioproc);
369 if (!amrex::ParallelDescriptor::IOProcessor()) {
371 tmp.resize(box,ncomp,amrex::The_Pinned_Arena());
373 tmp.resize(box,ncomp);
377 amrex::ParallelDescriptor::Bcast(tmp.dataPtr(), tmp.size(), ioproc);
380 amrex::Box fab_bx = tmp.box();
381 amrex::Dim3 dom_lb = lbound(domain);
382 fab_bx += amrex::IntVect(dom_lb.x,dom_lb.y,dom_lb.z);
384 fab_vars[iv]->resize(fab_bx,1);
386 amrex::Gpu::copy(amrex::Gpu::hostToDevice,
387 tmp.dataPtr(), tmp.dataPtr() + tmp.size(),
388 fab_vars[iv]->dataPtr());
391 fab_vars[iv]->copy(tmp,tmp.box(),0,fab_bx,0,1);
@ num
Definition: ERF_DataStruct.H:22
void ReadNetCDFFile(const std::string &fname, amrex::Vector< std::string > names, amrex::Vector< NDArray< DType > > &arrays, amrex::Vector< int > &success)
Definition: ERF_NCWpsFile.H:178
NC_Data_Dims_Type
Definition: ERF_NCWpsFile.H:29
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)
Definition: ERF_NCWpsFile.H:337
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)
Definition: ERF_NCWpsFile.H:264
amrex::Vector< amrex::FArrayBox > PlaneVector
Definition: ERF_NCWpsFile.H:13
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)
void ReadTimeSliceFromNetCDFFile(const std::string &fname, const int tidx, amrex::Vector< std::string > names, amrex::Vector< NDArray< DType > > &arrays, amrex::Vector< int > &success)
Definition: ERF_NCWpsFile.H:131
static NCFile open(const std::string &name, const int cmode=NC_NOWRITE)
Definition: ERF_NCInterface.cpp:707
Definition: ERF_NCWpsFile.H:50
bool owned
Definition: ERF_NCWpsFile.H:120
std::string name
Definition: ERF_NCWpsFile.H:117
size_t ndim()
Definition: ERF_NCWpsFile.H:104
NDArray()
Definition: ERF_NCWpsFile.H:60
std::string get_vname()
Definition: ERF_NCWpsFile.H:94
NDArray(const NDArray &array)
Definition: ERF_NCWpsFile.H:63
NDArray & operator=(const NDArray &array)
Definition: ERF_NCWpsFile.H:72
void set_vshape(std::vector< size_t > vshape)
Definition: ERF_NCWpsFile.H:112
NDArray(const std::string vname, const std::vector< size_t > &vshape)
Definition: ERF_NCWpsFile.H:54
std::atomic< size_t > ref_counted
Definition: ERF_NCWpsFile.H:119
typename std::remove_const< DataType >::type DType
Definition: ERF_NCWpsFile.H:51
std::vector< size_t > shape
Definition: ERF_NCWpsFile.H:118
DType * data
Definition: ERF_NCWpsFile.H:121
decltype(auto) get_data()
Definition: ERF_NCWpsFile.H:88
std::vector< size_t > get_vshape()
Definition: ERF_NCWpsFile.H:99
~NDArray()
Definition: ERF_NCWpsFile.H:82