Browse Source

Reading and caching multiple variables from the source.

interpolate
Michael Uleysky 1 year ago
parent
commit
df85e7ca51
  1. 57
      actions/actiontsc.h
  2. 73
      include/actiondep.h
  3. 2
      include/layereddata.h
  4. 2
      include/simple2ddata.h
  5. 28
      include/traits.h
  6. 42
      sources/AVISOLOCAL.cpp
  7. 2
      sources/AVISOLOCAL.h
  8. 31
      sources/BINFILE.cpp
  9. 2
      sources/BINFILE.h
  10. 10
      sources/MODISBINLOCAL.cpp
  11. 2
      sources/MODISBINLOCAL.h
  12. 95
      src/layereddata.cpp

57
actions/actiontsc.h

@ -24,10 +24,20 @@ template<class D> MString ActionTSC::DoAction(const CLArgs& args, D& ds)
auto [tindexes, err] = GetTIndexes(ds, args, pars); auto [tindexes, err] = GetTIndexes(ds, args, pars);
if(err.Exist()) return err; if(err.Exist()) return err;
if(!args.contains("var")) return "Variable not specified"; if(!args.contains("var")) return "Variables not specified";
MString vname = args.at("var"); auto vlist = michlib::Split_on_words(args.at("var"), " ,", false);
std::vector<MString> vnames{std::make_move_iterator(std::begin(vlist)), std::make_move_iterator(std::end(vlist))};
{
std::set<MString> used;
for(const auto& vname: vnames)
{
if(used.contains(vname)) return "Duplicate variable " + vname + " in list " + args.at("var");
if(!ds.CheckVar(vname)) return "Variable " + vname + " not exists in this dataset"; if(!ds.CheckVar(vname)) return "Variable " + vname + " not exists in this dataset";
pars.SetParameter("variable", vname); used.insert(vname);
}
}
pars.SetParameter("variables", args.at("var"));
std::unique_ptr<const BaseParameters> sourcepars; std::unique_ptr<const BaseParameters> sourcepars;
if constexpr(ParametersSupported<D>) if constexpr(ParametersSupported<D>)
@ -47,11 +57,15 @@ template<class D> MString ActionTSC::DoAction(const CLArgs& args, D& ds)
} }
auto p = sourcepars.get(); auto p = sourcepars.get();
auto data = Read(ds, vname, p, tindexes); auto data = Read(ds, vnames, p, tindexes);
if(!data) return "Can't read data"; if(data.size() != vnames.size()) return "Can't read data";
if(!data.Unit().Exist()) michlib::errmessage("Unknown measurement unit!");
if(!data.StandartName().Exist()) data.SetStandartName(NCFuncs::Name2StName(vname)); for(size_t i = 0; i < vnames.size(); i++)
if(!data.LongName().Exist()) data.SetStandartName(NCFuncs::Name2LongName(vname)); {
if(!data[i].Unit().Exist()) michlib::errmessage("Unknown measurement unit for variable " + vnames[i] + "!");
if(!data[i].StandartName().Exist()) data[i].SetStandartName(NCFuncs::Name2StName(vnames[i]));
if(!data[i].LongName().Exist()) data[i].SetStandartName(NCFuncs::Name2LongName(vnames[i]));
}
MString name = args.contains("out") ? args.at("out") : "out.bin"; MString name = args.contains("out") ? args.at("out") : "out.bin";
MString outfmt = args.contains("format") ? args.at("format") : (GetExt(name) == "nc" ? "nc" : "bin"); MString outfmt = args.contains("format") ? args.at("format") : (GetExt(name) == "nc" ? "nc" : "bin");
@ -60,16 +74,17 @@ template<class D> MString ActionTSC::DoAction(const CLArgs& args, D& ds)
{ {
BFileW fw; BFileW fw;
fw.Create(name, 3); size_t ic = 0;
fw.SetColumnName(1, "Longitude"); fw.Create(name, 2 + data.size());
fw.SetColumnName(2, "Latitude"); fw.SetColumnName(++ic, "Longitude");
fw.SetColumnName(3, vname + ", " + (data.Unit().Exist() ? data.Unit() : "unknown")); fw.SetColumnName(++ic, "Latitude");
for(size_t i = 0; i < data.size(); i++) fw.SetColumnName(++ic, vnames[i] + ", " + (data[i].Unit().Exist() ? data[i].Unit() : "unknown"));
fw.SetParameters(pars); fw.SetParameters(pars);
for(size_t i = 0; i < data.N(); i++) for(size_t r = 0; r < data[0].N(); r++)
{ {
fw.Write(data.Lon(i)); fw.Write(data[0].Lon(r));
fw.Write(data.Lat(i)); fw.Write(data[0].Lat(r));
fw.Write(data.IsFill(i) ? NAN : data(i)); for(size_t i = 0; i < data.size(); i++) fw.Write(data[i].IsFill(r) ? NAN : data[i](r));
} }
fw.Finalize(); fw.Finalize();
fw.Close(); fw.Close();
@ -84,10 +99,12 @@ template<class D> MString ActionTSC::DoAction(const CLArgs& args, D& ds)
if(args.contains("compress")) compress = args.at("compress").ToInt(); if(args.contains("compress")) compress = args.at("compress").ToInt();
if(!err.Exist()) err = fw.Create(data, name, args.at("_cmdline"), compress); if(!err.Exist()) err = fw.Create(data[0], name, args.at("_cmdline"), compress);
if(!err.Exist()) err = fw.AddVariable(vname, data.StandartName(), data.LongName(), data.Unit(), data.Comment()); for(size_t i = 0; i < data.size(); i++)
if(!err.Exist()) err = fw.WriteGrid(data); if(!err.Exist()) err = fw.AddVariable(vnames[i], data[i].StandartName(), data[i].LongName(), data[i].Unit(), data[i].Comment());
if(!err.Exist()) err = fw.WriteVariable(data, vname); if(!err.Exist()) err = fw.WriteGrid(data[0]);
for(size_t i = 0; i < data.size(); i++)
if(!err.Exist()) err = fw.WriteVariable(data[i], vnames[i]);
if(err.Exist()) return err; if(err.Exist()) return err;
fw.Close(); fw.Close();

73
include/actiondep.h

@ -105,41 +105,50 @@ template<class D> std::pair<TIndex, MString> GetTIndexes(const D& data, const CL
return {tindexes, ""}; return {tindexes, ""};
} }
template<class D> ReadType<D> Read(const D& data, const MString& vname, const BaseParameters* p, const TIndex& tindex) template<class D> std::vector<ReadType<D>> Read(const D& data, const std::vector<MString>& vnames, const BaseParameters* p, const TIndex& tindex)
{ {
using RT = ReadType<D>; using RT = ReadType<D>;
size_t ind; size_t ind;
std::vector<RT> out;
if(tindex.size() == 1) if(tindex.size() == 1)
{ {
ind = tindex[0]; ind = tindex[0];
michlib::message("Time: " + data.Time(ind).ToTString()); michlib::message("Time: " + data.Time(ind).ToTString());
std::map<MString, RT> cache;
for(const auto& vname: vnames)
{
bool res;
if constexpr(ReadPSupported<D>) if constexpr(ReadPSupported<D>)
return data.Read(vname, p, ind); res = data.Read(vname, cache, p, ind);
else if constexpr(ReadSupported<D>) else if constexpr(ReadSupported<D>)
return data.Read(vname, ind); res = data.Read(vname, cache, ind);
if(!res) return out;
}
for(size_t i = 0; i < vnames.size(); i++) out.emplace_back(std::move(cache[vnames[i]]));
} }
else else
{ {
Averager<RT> out; std::vector<Averager<RT>> aver(vnames.size());
bool ok = true;
for(size_t i = 0; i < tindex.size(); i++) for(size_t i = 0; i < tindex.size(); i++)
{ {
if(!ok) break;
ind = tindex[i]; ind = tindex[i];
michlib::message("Time: " + data.Time(ind).ToTString()); michlib::message("Time: " + data.Time(ind).ToTString());
RT dat; std::map<MString, RT> cache;
for(const auto& vname: vnames)
{
bool res;
if constexpr(ReadPSupported<D>) if constexpr(ReadPSupported<D>)
dat = data.Read(vname, p, ind); res = data.Read(vname, cache, p, ind);
else if constexpr(ReadSupported<D>) else if constexpr(ReadSupported<D>)
dat = data.Read(vname, ind); res = data.Read(vname, cache, ind);
if(dat) if(!res) return out;
out.Add(std::move(dat)); }
else for(size_t i = 0; i < vnames.size(); i++) aver[i].Add(std::move(cache[vnames[i]]));
ok = false;
} }
if(ok) return std::move(out.Div()); for(size_t i = 0; i < vnames.size(); i++) out.emplace_back(std::move(aver[i].Div()));
return out;
} }
return RT(); return out;
} }
template<class D> UVData<ReadType<D>> ReadUV(const D& data, const BaseParameters* p, size_t ind) template<class D> UVData<ReadType<D>> ReadUV(const D& data, const BaseParameters* p, size_t ind)
@ -147,18 +156,15 @@ template<class D> UVData<ReadType<D>> ReadUV(const D& data, const BaseParameters
using RT = ReadType<D>; using RT = ReadType<D>;
using UV = UVData<RT>; using UV = UVData<RT>;
michlib::message("Time: " + data.Time(ind).ToTString()); michlib::message("Time: " + data.Time(ind).ToTString());
RT u, v; std::map<MString, RT> cache;
bool res = false;
if constexpr(ReadPSupported<D>) if constexpr(ReadPSupported<D>)
{ res = data.Read("u", cache, p, ind) && data.Read("v", cache, p, ind);
u = data.Read("u", p, ind);
v = data.Read("v", p, ind);
}
else if constexpr(ReadSupported<D>) else if constexpr(ReadSupported<D>)
{ res = data.Read("u", cache, ind) && data.Read("v", cache, ind);
u = data.Read("u", ind); if(!res) return UV();
v = data.Read("v", ind);
} return UV(cache.at("u"), cache.at("v"));
return UV(u, v);
} }
template<class D> UVData<ReadType<D>> ReadUV(const D& data, const BaseParameters* p, const TIndex& tindex) template<class D> UVData<ReadType<D>> ReadUV(const D& data, const BaseParameters* p, const TIndex& tindex)
@ -178,18 +184,15 @@ template<class D> UVData<ReadType<D>> ReadUV(const D& data, const BaseParameters
if(!ok) break; if(!ok) break;
ind = tindex[i]; ind = tindex[i];
michlib::message("Time: " + data.Time(ind).ToTString()); michlib::message("Time: " + data.Time(ind).ToTString());
RT u, v; std::map<MString, RT> cache;
bool res = false;
if constexpr(ReadPSupported<D>) if constexpr(ReadPSupported<D>)
{ res = data.Read("u", cache, p, ind) && data.Read("v", cache, p, ind);
u = data.Read("u", p, ind);
v = data.Read("v", p, ind);
}
else if constexpr(ReadSupported<D>) else if constexpr(ReadSupported<D>)
{ res = data.Read("u", cache, ind) && data.Read("v", cache, ind);
u = data.Read("u", ind); if(!res) return UV();
v = data.Read("v", ind);
} UV dat(cache.at("u"), cache.at("v"));
UV dat(u, v);
if(dat) if(dat)
out.Add(std::move(dat)); out.Add(std::move(dat));
else else

2
include/layereddata.h

@ -135,7 +135,7 @@ class LayeredData: public NCFuncs
std::pair<const BaseParameters*, MString> Parameters(michlib_internal::ParameterListEx& pars, const CLArgs& args, const struct Region& reg) const; std::pair<const BaseParameters*, MString> Parameters(michlib_internal::ParameterListEx& pars, const CLArgs& args, const struct Region& reg) const;
Data Read(const MString& vname, const BaseParameters* ip, size_t i) const; bool Read(const MString& vname, std::map<MString, Data>& cache, const BaseParameters* ip, size_t i) const;
bool isOk() const { return nc.size() > 0; } bool isOk() const { return nc.size() > 0; }

2
include/simple2ddata.h

@ -18,6 +18,8 @@ class Simple2DData: public BaseData
{ {
} }
Simple2DData CopyGrid() const { return {nx, ny, x0, y0, xstep, ystep}; }
const real& V(size_t ix, size_t iy) const { return V(iy * nx + ix); } const real& V(size_t ix, size_t iy) const { return V(iy * nx + ix); }
real& V(size_t ix, size_t iy) { return V(iy * nx + ix); } real& V(size_t ix, size_t iy) { return V(iy * nx + ix); }

28
include/traits.h

@ -31,16 +31,16 @@ template<class T>
concept ParametersSupported = ParametersRequiredRegion<T> || ParametersNotRequiredRegion<T>; concept ParametersSupported = ParametersRequiredRegion<T> || ParametersNotRequiredRegion<T>;
template<class T> template<class T>
concept ReadPSupported = requires(T t, const MString& vname, const BaseParameters* ip, size_t i) { concept ReadPSupported = requires(T t, const MString& vname, std::map<MString, typename T::Data>& cache, const BaseParameters* ip, size_t i) {
{ {
t.Read(vname, ip, i)(0) t.Read(vname, cache, ip, i)
} -> std::convertible_to<real>; } -> std::same_as<bool>;
}; };
template<class T> template<class T>
concept ReadSupported = requires(T t, const MString& vname, size_t i) { concept ReadSupported = requires(T t, const MString& vname, std::map<MString, typename T::Data>& cache, size_t i) {
{ {
t.Read(vname, i)(0) t.Read(vname, cache, i)
} -> std::convertible_to<real>; } -> std::same_as<bool>;
}; };
template<class T> template<class T>
@ -70,20 +70,14 @@ concept HaveXYStep = requires(T t) {
template<typename...> using void_ = void; template<typename...> using void_ = void;
template<class D, bool RP, bool R> struct GetReadType_; template<class D, bool R> struct GetReadType_;
template<class D, bool R> struct GetReadType_<D, true, R>
{
using p = BaseParameters*;
using type = decltype(D().Read(MString(), p(), 0));
};
template<class D> struct GetReadType_<D, false, true> template<class D> struct GetReadType_<D, true>
{ {
using type = decltype(D().Read(MString(), 0)); using type = typename D::Data;
}; };
template<class D> struct GetReadType_<D, false, false> template<class D> struct GetReadType_<D, false>
{ {
using type = std::decay_t<D>; using type = std::decay_t<D>;
}; };
@ -125,7 +119,7 @@ template<class Act, class S> consteval bool IsDisabled()
} }
} }
template<class D> using ReadType = internal::GetReadType_<D, ReadPSupported<D>, ReadSupported<D>>::type; template<class D> using ReadType = internal::GetReadType_<D, ReadPSupported<D> || ReadSupported<D>>::type;
template<class T> template<class T>
concept ReadIsGrid = requires { concept ReadIsGrid = requires {

42
sources/AVISOLOCAL.cpp

@ -86,9 +86,10 @@ std::pair<const BaseParameters*, MString> AVISOLOCALData::Parameters(michlib_int
return {ppar.release(), ""}; return {ppar.release(), ""};
} }
AVISOLOCALData::Data AVISOLOCALData::Read(const MString& vname, const BaseParameters* ip, size_t i) const bool AVISOLOCALData::Read(const MString& vname, std::map<MString, AVISOLOCALData::Data>& cache, const BaseParameters* ip, size_t i) const
{ {
if(!isOk()) return Data(); if(cache.contains(vname)) return true;
if(!isOk()) return false;
auto p = dynamic_cast<const struct Parameters*>(ip); auto p = dynamic_cast<const struct Parameters*>(ip);
NCFileA nc; NCFileA nc;
@ -96,7 +97,7 @@ AVISOLOCALData::Data AVISOLOCALData::Read(const MString& vname, const BaseParame
bool isfloat = false, isint2 = false, isint4 = false; bool isfloat = false, isint2 = false, isint4 = false;
nc.Reset(datapath + "/uv-" + times[i].ToString() + ".nc"); nc.Reset(datapath + "/uv-" + times[i].ToString() + ".nc");
if(!nc) return Data(); if(!nc) return false;
{ {
auto head = nc.Header(); auto head = nc.Header();
for(const auto& v: head.Variables()) for(const auto& v: head.Variables())
@ -119,27 +120,36 @@ AVISOLOCALData::Data AVISOLOCALData::Read(const MString& vname, const BaseParame
if(vname == "U" || vname == "U2") if(vname == "U" || vname == "U2")
{ {
bool square = vname == "U2"; bool square = vname == "U2";
auto u = Read("u", ip, i); if(!(Read("u", cache, ip, i) && Read("v", cache, ip, i))) return false;
auto v = Read("v", ip, i); cache[vname] = cache.at("u").CopyGrid();
if(!(u && v)) return Data(); auto& U = cache.at(vname);
if(square && u.Unit().Exist()) u.SetUnit("(" + u.Unit() + ")2"); const auto& u = cache.at("u");
for(size_t ind = 0; ind < u.N(); ind++) const auto& v = cache.at("v");
if(u.Unit().Exist()) U.SetUnit(square ? ("(" + u.Unit() + ")2") : u.Unit());
for(size_t ind = 0; ind < U.N(); ind++)
{ {
if(u.IsFill(ind) || v.IsFill(ind)) if(u.IsFill(ind) || v.IsFill(ind))
u.V(ind) = u.Fillval(); U.V(ind) = U.Fillval();
else else
u.V(ind) = square ? (u(ind) * u(ind) + v(ind) * v(ind)) : michlib::Hypot(u(ind), v(ind)); U.V(ind) = square ? (u(ind) * u(ind) + v(ind) * v(ind)) : michlib::Hypot(u(ind), v(ind));
} }
return u; return true;
} }
return Data(); return false;
} }
// Direct read // Direct read
if(isint2) return ReadVarRaw<int2>(nc, name, p); Data data;
if(isint4) return ReadVarRaw<int4>(nc, name, p); if(isint2) data = ReadVarRaw<int2>(nc, name, p);
if(isfloat) return ReadVarRaw<float>(nc, name, p); if(isint4) data = ReadVarRaw<int4>(nc, name, p);
return Data(); if(isfloat) data = ReadVarRaw<float>(nc, name, p);
if(data)
{
cache[vname] = std::move(data);
return true;
}
return false;
} }
template<class DataType> AVISOLOCALData::Data AVISOLOCALData::ReadVarRaw(const NCFileA& nc, const MString& name, const struct AVISOLOCALData::Parameters* p) const template<class DataType> AVISOLOCALData::Data AVISOLOCALData::ReadVarRaw(const NCFileA& nc, const MString& name, const struct AVISOLOCALData::Parameters* p) const

2
sources/AVISOLOCAL.h

@ -60,7 +60,7 @@ class AVISOLOCALData: public NCFuncs
return NCFuncs::CheckVar(vname, [&nc = std::as_const(nc)](const MString& vn) { return HaveVar(nc, vn); }); return NCFuncs::CheckVar(vname, [&nc = std::as_const(nc)](const MString& vn) { return HaveVar(nc, vn); });
} }
Data Read(const MString& vname, const BaseParameters* ip, size_t i) const; bool Read(const MString& vname, std::map<MString, Data>& cache, const BaseParameters* ip, size_t i) const;
template<class DataType> Data ReadVarRaw(const NCFileA& nc, const MString& name, const struct Parameters* p) const; template<class DataType> Data ReadVarRaw(const NCFileA& nc, const MString& name, const struct Parameters* p) const;
}; };

31
sources/BINFILE.cpp

@ -53,9 +53,10 @@ MString BINFILEData::Open(const CLArgs& args)
return ""; return "";
} }
BINFILEData::Data BINFILEData::Read(const MString& vname, size_t i) const bool BINFILEData::Read(const MString& vname, std::map<MString, BINFILEData::Data>& cache, size_t i) const
{ {
if(!isOk()) return Data(); if(cache.contains(vname)) return true;
if(!isOk()) return false;
// Only rectangular grids are supported // Only rectangular grids are supported
real xs = (data->Lon(data->Nx() - 1, data->Ny() - 1) - data->Lon(0, 0)) / (data->Nx() - 1); real xs = (data->Lon(data->Nx() - 1, data->Ny() - 1) - data->Lon(0, 0)) / (data->Nx() - 1);
@ -66,18 +67,20 @@ BINFILEData::Data BINFILEData::Read(const MString& vname, size_t i) const
if(vname == "U" || vname == "U2") if(vname == "U" || vname == "U2")
{ {
bool square = vname == "U2"; bool square = vname == "U2";
auto u = Read("u", i); if(!(Read("u", cache, i) && Read("v", cache, i))) return false;
auto v = Read("v", i); cache[vname] = cache.at("u").CopyGrid();
if(!(u && v)) return Data(); auto& U = cache.at(vname);
if(square && u.Unit().Exist()) u.SetUnit("(" + u.Unit() + ")2"); const auto& u = cache.at("u");
for(size_t ind = 0; ind < u.N(); ind++) const auto& v = cache.at("v");
if(u.Unit().Exist()) U.SetUnit(square ? ("(" + u.Unit() + ")2") : u.Unit());
for(size_t ind = 0; ind < U.N(); ind++)
{ {
if(u.IsFill(ind) || v.IsFill(ind)) if(u.IsFill(ind) || v.IsFill(ind))
u.V(ind) = u.Fillval(); U.V(ind) = U.Fillval();
else else
u.V(ind) = square ? (u(ind) * u(ind) + v(ind) * v(ind)) : michlib::Hypot(u(ind), v(ind)); U.V(ind) = square ? (u(ind) * u(ind) + v(ind) * v(ind)) : michlib::Hypot(u(ind), v(ind));
} }
return u; return true;
} }
if(vname == "u" || vname == "v") if(vname == "u" || vname == "v")
@ -91,7 +94,13 @@ BINFILEData::Data BINFILEData::Read(const MString& vname, size_t i) const
else else
out.V(ix, iy) = out.Fillval(); out.V(ix, iy) = out.Fillval();
} }
if(out)
{
cache[vname] = std::move(out);
return true;
}
} }
return out; return false;
} }

2
sources/BINFILE.h

@ -44,5 +44,5 @@ class BINFILEData
bool CheckVar(const MString& vname) const { return vname == "u" || vname == "v" || vname == "U" || vname == "U2"; } bool CheckVar(const MString& vname) const { return vname == "u" || vname == "v" || vname == "U" || vname == "U2"; }
Data Read(const MString& vname, size_t i) const; bool Read(const MString& vname, std::map<MString, Data>& cache, size_t i) const;
}; };

10
sources/MODISBINLOCAL.cpp

@ -83,8 +83,9 @@ std::pair<const BaseParameters*, MString> MODISBINLOCALData::Parameters(michlib_
return {ppar.release(), ""}; return {ppar.release(), ""};
} }
MODISBINLOCALData::Data MODISBINLOCALData::Read(const MString& var, const BaseParameters* ip, size_t tind) const bool MODISBINLOCALData::Read(const MString& vname, std::map<MString, MODISBINLOCALData::Data>& cache, const BaseParameters* ip, size_t tind) const
{ {
if(cache.contains(vname)) return true;
auto p = pointer_cast<const struct Parameters*>(ip); auto p = pointer_cast<const struct Parameters*>(ip);
Averager<Data> out; Averager<Data> out;
@ -94,7 +95,12 @@ MODISBINLOCALData::Data MODISBINLOCALData::Read(const MString& var, const BasePa
Data dat = ReadFile(suf, p, tind); Data dat = ReadFile(suf, p, tind);
if(dat) out.Add(std::move(dat)); if(dat) out.Add(std::move(dat));
} }
return std::move(out.Div()); if(out)
{
cache[vname] = std::move(out.Div());
return true;
}
return false;
} }
MODISBINLOCALData::Data MODISBINLOCALData::ReadFile(const MString& suf, const struct MODISBINLOCALData::Parameters* p, size_t tind) const MODISBINLOCALData::Data MODISBINLOCALData::ReadFile(const MString& suf, const struct MODISBINLOCALData::Parameters* p, size_t tind) const

2
sources/MODISBINLOCAL.h

@ -40,7 +40,7 @@ class MODISBINLOCALData
MString Info() const; MString Info() const;
MString Open(const CLArgs& args); MString Open(const CLArgs& args);
Data Read(const MString& var, const BaseParameters* ip, size_t tind) const; bool Read(const MString& vname, std::map<MString, Data>& cache, const BaseParameters* ip, size_t tind) const;
Data ReadFile(const MString& suf, const struct Parameters* p, size_t tind) const; Data ReadFile(const MString& suf, const struct Parameters* p, size_t tind) const;
std::pair<const BaseParameters*, MString> Parameters(michlib_internal::ParameterListEx& pars, const CLArgs& args, const struct Region& reg) const; std::pair<const BaseParameters*, MString> Parameters(michlib_internal::ParameterListEx& pars, const CLArgs& args, const struct Region& reg) const;

95
src/layereddata.cpp

@ -178,10 +178,10 @@ std::pair<const BaseParameters*, MString> LayeredData::Parameters(michlib_intern
return {ppar.release(), ""}; return {ppar.release(), ""};
} }
LayeredData::Data LayeredData::Read(const MString& vname, const BaseParameters* ip, size_t i) const bool LayeredData::Read(const MString& vname, std::map<MString, LayeredData::Data>& cache, const BaseParameters* ip, size_t i) const
{ {
if(!isOk()) return Data(); if(cache.contains(vname)) return true;
bool nodepth = false; if(!isOk()) return false;
auto p = dynamic_cast<const struct Parameters*>(ip); auto p = dynamic_cast<const struct Parameters*>(ip);
auto [name, id, tid] = VarNameLoc(vname, times[i]); auto [name, id, tid] = VarNameLoc(vname, times[i]);
@ -190,86 +190,101 @@ LayeredData::Data LayeredData::Read(const MString& vname, const BaseParameters*
// ptemp from temp and sal // ptemp from temp and sal
if(vname == "ptemp") if(vname == "ptemp")
{ {
auto temp = Read("temp", ip, i); if(!(Read("temp", cache, ip, i) && Read("sal", cache, ip, i))) return false;
auto sal = Read("sal", ip, i); cache[vname] = cache.at("temp").CopyGrid();
if(!(temp && sal)) return Data(); auto& ptemp = cache.at(vname);
temp.SetUnit("degrees_C"); const auto& temp = cache.at("temp");
for(size_t ind = 0; ind < temp.N(); ind++) const auto& sal = cache.at("sal");
ptemp.SetUnit("degrees_C");
for(size_t ind = 0; ind < ptemp.N(); ind++)
{ {
if(temp.IsFill(ind) || sal.IsFill(ind)) if(temp.IsFill(ind) || sal.IsFill(ind))
temp.V(ind) = temp.Fillval(); ptemp.V(ind) = ptemp.Fillval();
else else
temp.V(ind) = Temp2PTemp(temp.V(ind), sal.V(ind), Depth(p->layer), temp.Lon(ind), temp.Lat(ind)); ptemp.V(ind) = Temp2PTemp(temp.V(ind), sal.V(ind), Depth(p->layer), ptemp.Lon(ind), ptemp.Lat(ind));
} }
return temp; return true;
} }
// temp from ptemp and sal // temp from ptemp and sal
if(vname == "temp") if(vname == "temp")
{ {
auto temp = Read("ptemp", ip, i); if(!(Read("ptemp", cache, ip, i) && Read("sal", cache, ip, i))) return false;
auto sal = Read("sal", ip, i); cache[vname] = cache.at("ptemp").CopyGrid();
if(!(temp && sal)) return Data(); auto& temp = cache.at(vname);
const auto& ptemp = cache.at("ptemp");
const auto& sal = cache.at("sal");
temp.SetUnit("degrees_C"); temp.SetUnit("degrees_C");
for(size_t ind = 0; ind < temp.N(); ind++) for(size_t ind = 0; ind < temp.N(); ind++)
{ {
if(temp.IsFill(ind) || sal.IsFill(ind)) if(ptemp.IsFill(ind) || sal.IsFill(ind))
temp.V(ind) = temp.Fillval(); temp.V(ind) = temp.Fillval();
else else
temp.V(ind) = PTemp2Temp(temp.V(ind), sal.V(ind), Depth(p->layer), temp.Lon(ind), temp.Lat(ind)); temp.V(ind) = PTemp2Temp(ptemp.V(ind), sal.V(ind), Depth(p->layer), temp.Lon(ind), temp.Lat(ind));
} }
return temp; return true;
} }
// pdens from temp and sal // pdens from temp and sal
if(vname == "pdens") if(vname == "pdens")
{ {
bool tempispot = HaveVar("ptemp"); bool tempispot = HaveVar("ptemp");
auto temp = Read(tempispot ? "ptemp" : "temp", ip, i); if(!(Read(tempispot ? "ptemp" : "temp", cache, ip, i) && Read("sal", cache, ip, i))) return false;
auto sal = Read("sal", ip, i); cache[vname] = cache.at("sal").CopyGrid();
if(!(temp && sal)) return Data(); auto& pdens = cache.at(vname);
temp.SetUnit("kg m-3"); const auto& temp = cache.at(tempispot ? "ptemp" : "temp");
for(size_t ind = 0; ind < temp.N(); ind++) const auto& sal = cache.at("sal");
pdens.SetUnit("kg m-3");
for(size_t ind = 0; ind < pdens.N(); ind++)
{ {
if(temp.IsFill(ind) || sal.IsFill(ind)) if(temp.IsFill(ind) || sal.IsFill(ind))
temp.V(ind) = temp.Fillval(); pdens.V(ind) = pdens.Fillval();
else else
temp.V(ind) = tempispot ? PTemp2PDens(temp.V(ind), sal.V(ind), Depth(p->layer), temp.Lon(ind), temp.Lat(ind)) pdens.V(ind) = tempispot ? PTemp2PDens(temp.V(ind), sal.V(ind), Depth(p->layer), pdens.Lon(ind), pdens.Lat(ind))
: Temp2PDens(temp.V(ind), sal.V(ind), Depth(p->layer), temp.Lon(ind), temp.Lat(ind)); : Temp2PDens(temp.V(ind), sal.V(ind), Depth(p->layer), pdens.Lon(ind), pdens.Lat(ind));
} }
return temp; return true;
} }
// U and U2 from u and v // U and U2 from u and v
if(vname == "U" || vname == "U2") if(vname == "U" || vname == "U2")
{ {
bool square = vname == "U2"; bool square = vname == "U2";
auto u = Read("u", ip, i); if(!(Read("u", cache, ip, i) && Read("v", cache, ip, i))) return false;
auto v = Read("v", ip, i); cache[vname] = cache.at("u").CopyGrid();
if(!(u && v)) return Data(); auto& U = cache.at(vname);
if(square && u.Unit().Exist()) u.SetUnit("(" + u.Unit() + ")2"); const auto& u = cache.at("u");
for(size_t ind = 0; ind < u.N(); ind++) const auto& v = cache.at("v");
if(u.Unit().Exist()) U.SetUnit(square ? ("(" + u.Unit() + ")2") : u.Unit());
for(size_t ind = 0; ind < U.N(); ind++)
{ {
if(u.IsFill(ind) || v.IsFill(ind)) if(u.IsFill(ind) || v.IsFill(ind))
u.V(ind) = u.Fillval(); U.V(ind) = U.Fillval();
else else
u.V(ind) = square ? (u(ind) * u(ind) + v(ind) * v(ind)) : michlib::Hypot(u(ind), v(ind)); U.V(ind) = square ? (u(ind) * u(ind) + v(ind) * v(ind)) : michlib::Hypot(u(ind), v(ind));
} }
return u; return true;
} }
return Data(); return false;
} }
// Direct read // Direct read
bool nodepth = false;
Data data;
auto head = nc[id]->Header(); auto head = nc[id]->Header();
for(const auto& v: head.Variables()) for(const auto& v: head.Variables())
if(v.Name() == name) if(v.Name() == name)
{ {
if(v.Dimensions().size() == 3) nodepth = true; if(v.Dimensions().size() == 3) nodepth = true;
if(v.Type().Id() == NC_SHORT) return ReadVarRaw<int2>(nc[id], name, tid, nodepth, p); if(v.Type().Id() == NC_SHORT) data = ReadVarRaw<int2>(nc[id], name, tid, nodepth, p);
if(v.Type().Id() == NC_INT) return ReadVarRaw<int>(nc[id], name, tid, nodepth, p); if(v.Type().Id() == NC_INT) data = ReadVarRaw<int>(nc[id], name, tid, nodepth, p);
if(v.Type().Id() == NC_FLOAT) return ReadVarRaw<float>(nc[id], name, tid, nodepth, p); if(v.Type().Id() == NC_FLOAT) data = ReadVarRaw<float>(nc[id], name, tid, nodepth, p);
if(v.Type().Id() == NC_DOUBLE) return ReadVarRaw<double>(nc[id], name, tid, nodepth, p); if(v.Type().Id() == NC_DOUBLE) data = ReadVarRaw<double>(nc[id], name, tid, nodepth, p);
if(data)
{
cache[vname] = std::move(data);
return true;
}
} }
return Data(); return false;
} }
template<class DataType> LayeredData::Data LayeredData::ReadVarRaw(const NC& f, const MString& name, size_t i, bool nodepth, const struct LayeredData::Parameters* p) const template<class DataType> LayeredData::Data LayeredData::ReadVarRaw(const NC& f, const MString& name, size_t i, bool nodepth, const struct LayeredData::Parameters* p) const

Loading…
Cancel
Save