Compare commits
67 Commits
interpolat
...
master
57 changed files with 7146 additions and 821 deletions
@ -0,0 +1,130 @@
|
||||
#define MICHLIB_NOSOURCE |
||||
#include "actiongrad.h" |
||||
|
||||
MString GradMethods::NCFileW::Create(const MString& name, const MString& history, const std::vector<MString>& vnames, const std::vector<MString>& lnames, |
||||
const std::vector<GradMethods::DataType>& lons, const std::vector<GradMethods::DataType>& lats, int compress) |
||||
{ |
||||
const float node_offset = 0.0; |
||||
|
||||
auto nx = lons.size(); |
||||
auto ny = lats.size(); |
||||
|
||||
const auto fill = static_cast<GradMethods::DataType>(std::numeric_limits<GradMethods::Matrix::MDataType>::max()); |
||||
|
||||
// Creating
|
||||
Open(name); |
||||
if(!*this) return "Can't create netcdf file " + name + ": " + ErrMessage(); |
||||
|
||||
AddAtt("history", history); |
||||
AddAtt("node_offset", node_offset); |
||||
|
||||
AddDim("longitude", nx); |
||||
AddDim("latitude", ny); |
||||
AddVar("longitude", NC_FLOAT, "longitude"); |
||||
AddVar("latitude", NC_FLOAT, "latitude"); |
||||
SetComp("longitude", compress); |
||||
SetComp("latitude", compress); |
||||
AddAtt("longitude", "standard_name", "longitude"); |
||||
AddAtt("longitude", "long_name", "Longitude"); |
||||
AddAtt("latitude", "standard_name", "latitude"); |
||||
AddAtt("latitude", "long_name", "Latitude"); |
||||
|
||||
// Variables
|
||||
for(size_t i = 0; i < vnames.size(); i++) |
||||
{ |
||||
AddVar(vnames[i], NC_FLOAT, "latitude", "longitude"); |
||||
SetComp(vnames[i], compress); |
||||
if(lnames[i].Exist()) AddAtt(vnames[i], "long_name", lnames[i]); |
||||
AddAtt(vnames[i], "_FillValue", fill); |
||||
} |
||||
|
||||
// End definitions
|
||||
EndDef(); |
||||
|
||||
// Writing lon, lat
|
||||
WriteVar("longitude", lons.data()); |
||||
WriteVar("latitude", lats.data()); |
||||
|
||||
if(!*this) return "Can't set grid in the netcdf file " + name + ": " + ErrMessage(); |
||||
|
||||
return ""; |
||||
} |
||||
|
||||
MString GradMethods::NCFileW::WriteVariable(const MString& name, const GradMethods::Matrix& data) |
||||
{ |
||||
WriteVar(name, data.Data().data()); |
||||
if(!*this) return "Can't write variable " + name + ": " + ErrMessage(); |
||||
|
||||
return ""; |
||||
} |
||||
|
||||
GradMethods::Matrix::Matrix(const std::vector<GradMethods::DataType>& in, size_t nx_, size_t ny_, struct GradMethods::MinMax minmax): nx(nx_), ny(ny_), data(nx_ * ny_) |
||||
{ |
||||
if(minmax.automin || minmax.automax) |
||||
{ |
||||
DataType min = in[0]; |
||||
DataType max = in[0]; |
||||
for(size_t i = 1; i < in.size(); i++) |
||||
if(in[i] != minmax.fill) |
||||
{ |
||||
min = std::min(min, in[i]); |
||||
max = std::max(max, in[i]); |
||||
} |
||||
if(minmax.automin) minmax.min = min; |
||||
if(minmax.automax) minmax.max = max; |
||||
} |
||||
|
||||
if(minmax.log) |
||||
{ |
||||
minmax.min = michlib_internal::RealType<sizeof(DataType)>::Log(minmax.min); |
||||
minmax.max = michlib_internal::RealType<sizeof(DataType)>::Log(minmax.max); |
||||
} |
||||
|
||||
DataType a = (std::numeric_limits<MDataType>::max() - 1) / (minmax.max - minmax.min); |
||||
for(size_t i = 1; i < in.size(); i++) |
||||
{ |
||||
DataType v = minmax.log ? michlib_internal::RealType<sizeof(DataType)>::Log(in[i]) : in[i]; |
||||
if(in[i] == minmax.fill) |
||||
data[i] = std::numeric_limits<MDataType>::max(); |
||||
else if(v <= minmax.min) |
||||
data[i] = 0; |
||||
else if(v >= minmax.max) |
||||
data[i] = std::numeric_limits<MDataType>::max() - 1; |
||||
else |
||||
data[i] = static_cast<MDataType>(michlib_internal::RealType<sizeof(DataType)>::Round(a * (v - minmax.min))); |
||||
} |
||||
} |
||||
|
||||
void GradMethods::Matrix::Grad() |
||||
{ |
||||
std::vector<MDataType> out(data.size()); |
||||
const auto bad = std::numeric_limits<MDataType>::max(); |
||||
|
||||
for(size_t iy = 0; iy < ny; iy++) |
||||
for(size_t ix = 0; ix < nx; ix++) |
||||
{ |
||||
if(iy < 1 || ix < 1 || iy > ny - 2 || ix > nx - 2) |
||||
out[iy * nx + ix] = bad; |
||||
else if(V(ix - 1, iy - 1) == bad || V(ix, iy - 1) == bad || V(ix + 1, iy - 1) == bad || V(ix - 1, iy) == bad || V(ix, iy) == bad || V(ix + 1, iy) == bad || |
||||
V(ix - 1, iy + 1) == bad || V(ix, iy + 1) == bad || V(ix + 1, iy + 1) == bad) |
||||
out[iy * nx + ix] = bad; |
||||
else |
||||
{ |
||||
using IT = michlib::int4; |
||||
// Possible but unlikely overflow
|
||||
const IT m1 = -1; |
||||
const IT m2 = -2; |
||||
const IT p1 = 1; |
||||
const IT p2 = 2; |
||||
|
||||
IT gx = m1 * V(ix - 1, iy + 1) + p1 * V(ix + 1, iy + 1) + m2 * V(ix - 1, iy) + p2 * V(ix + 1, iy) + m1 * V(ix - 1, iy - 1) + p1 * V(ix + 1, iy - 1); |
||||
IT gy = m1 * V(ix - 1, iy - 1) + p1 * V(ix - 1, iy + 1) + m2 * V(ix, iy - 1) + p2 * V(ix, iy + 1) + m1 * V(ix + 1, iy - 1) + p1 * V(ix + 1, iy + 1); |
||||
|
||||
auto sq = static_cast<IT>(michlib::Round(michlib::Hypot(gx, gy))); |
||||
if(sq >= bad) sq = bad - 1; |
||||
out[iy * nx + ix] = static_cast<MDataType>(sq); |
||||
} |
||||
} |
||||
|
||||
data = std::move(out); |
||||
} |
||||
@ -0,0 +1,133 @@
|
||||
#pragma once |
||||
#include "actiondep.h" |
||||
#include "ncfilew.h" |
||||
#include "ncfuncs.h" |
||||
|
||||
class GradMethods |
||||
{ |
||||
public: |
||||
using DataType = float; |
||||
|
||||
static DataType ToNum(const MString& str) { return michlib_internal::RealType<sizeof(DataType)>::String2Real(str.Buf()); } |
||||
|
||||
struct MinMax |
||||
{ |
||||
bool automin, automax, log; |
||||
DataType min, max; |
||||
DataType fill; |
||||
}; |
||||
|
||||
class Matrix; |
||||
class NCFileW; |
||||
}; |
||||
|
||||
class GradMethods::Matrix |
||||
{ |
||||
public: |
||||
using MDataType = michlib::uint2; |
||||
|
||||
private: |
||||
size_t nx, ny; |
||||
std::vector<MDataType> data; |
||||
|
||||
public: |
||||
Matrix(const std::vector<DataType>& in, size_t nx_, size_t ny_, struct MinMax minmax); |
||||
|
||||
void Grad(); |
||||
|
||||
auto Nx() const { return nx; } |
||||
auto Ny() const { return ny; } |
||||
|
||||
const auto& V(size_t ix, size_t iy) const { return data[iy * nx + ix]; } |
||||
|
||||
auto& V(size_t ix, size_t iy) { return data[iy * nx + ix]; } |
||||
|
||||
const auto& Data() const { return data; } |
||||
}; |
||||
|
||||
class GradMethods::NCFileW: public NCFileWBase |
||||
{ |
||||
public: |
||||
MString Create(const MString& name, const MString& history, const std::vector<MString>& vnames, const std::vector<MString>& lnames, const std::vector<GradMethods::DataType>& lons, |
||||
const std::vector<GradMethods::DataType>& lats, int compress); |
||||
MString WriteVariable(const MString& name, const GradMethods::Matrix& data); |
||||
}; |
||||
|
||||
template<class T> |
||||
concept GradSupported = requires(T t, const MString& vname) { |
||||
{ |
||||
t.ReadVar(vname) |
||||
} -> std::same_as<std::vector<GradMethods::DataType>>; |
||||
}; |
||||
|
||||
ADD_ACTION(GRAD, grad, GradSupported<Source>, GradMethods); |
||||
|
||||
template<class D> MString ActionGRAD::DoAction(const CLArgs& args, D& ds) |
||||
{ |
||||
auto resop = ds.Open(args); |
||||
if(resop.Exist()) return "Can't open source: " + resop; |
||||
|
||||
MString name = args.contains("out") ? args.at("out") : "out.nc"; |
||||
|
||||
MString min = args.contains("min") ? args.at("min") : "auto"; |
||||
MString max = args.contains("max") ? args.at("max") : "auto"; |
||||
|
||||
int compress = args.contains("compress") ? args.at("compress").ToInt() : 3; |
||||
|
||||
std::vector<Matrix> data; |
||||
|
||||
std::vector<MString> lnames; |
||||
|
||||
// Read data
|
||||
for(size_t i = 0; i < ds.NVar(); i++) |
||||
{ |
||||
const MString& name = ds.VarNames()[i]; |
||||
const MString& lname = ds.LongNames()[i]; |
||||
bool hmin = args.contains(name + "_min"); |
||||
bool hmax = args.contains(name + "_max"); |
||||
struct MinMax minmax; |
||||
|
||||
minmax.log = args.contains(name + "_log"); |
||||
|
||||
if(hmin) |
||||
{ |
||||
MString vmin = args.at(name + "_min"); |
||||
minmax.automin = (vmin == "auto"); |
||||
minmax.min = ToNum(vmin); |
||||
} |
||||
else |
||||
{ |
||||
minmax.automin = (min == "auto"); |
||||
minmax.min = ToNum(min); |
||||
} |
||||
|
||||
if(hmax) |
||||
{ |
||||
MString vmax = args.at(name + "_max"); |
||||
minmax.automax = (vmax == "auto"); |
||||
minmax.max = ToNum(vmax); |
||||
} |
||||
else |
||||
{ |
||||
minmax.automax = (max == "auto"); |
||||
minmax.max = ToNum(max); |
||||
} |
||||
|
||||
minmax.fill = ds.FillVal(name); |
||||
|
||||
data.emplace_back(ds.ReadVar(name), ds.Nx(), ds.Ny(), minmax); |
||||
lnames.emplace_back(lname + ", gradient"); |
||||
} |
||||
|
||||
NCFileW fw; |
||||
fw.Create(name, (ds.History().Exist() ? (ds.History() + "; ") : "") + args.at("_cmdline"), ds.VarNames(), lnames, ds.ReadLons(), ds.ReadLats(), compress); |
||||
|
||||
for(size_t i = 0; i < ds.NVar(); i++) |
||||
{ |
||||
const MString& name = ds.VarNames()[i]; |
||||
data[i].Grad(); |
||||
fw.WriteVariable(name, data[i]); |
||||
} |
||||
|
||||
return ""; |
||||
}; |
||||
@ -0,0 +1,210 @@
|
||||
#pragma once |
||||
#include "BFileR.h" |
||||
#include "BFileW.h" |
||||
#include "actiondep.h" |
||||
#include "interpolation.h" |
||||
#include "ncfilew.h" |
||||
#include "ncfuncs.h" |
||||
#include <memory> |
||||
|
||||
template<class T> |
||||
concept ReadIsInterpolable = requires { |
||||
{ |
||||
std::declval<ReadType<T>>().GridPos(0.0, 0.0) |
||||
} -> std::convertible_to<struct GridPoint>; |
||||
}; |
||||
|
||||
ADD_ACTION(INT, int, (ReadPSupported<Source> || ReadSupported<Source>)&&ReadIsInterpolable<Source>); |
||||
|
||||
template<class D> MString ActionINT::DoAction(const CLArgs& args, D& ds) |
||||
{ |
||||
struct TPoint |
||||
{ |
||||
struct Point2D p; |
||||
MDateTime t; |
||||
real lon0, lat0, dtime; |
||||
}; |
||||
|
||||
MDateTime refdate("1980-01-01"); |
||||
|
||||
if(!args.contains("input")) return "No input file given"; |
||||
michlib::BFileR in; |
||||
if(in.Open(args.at("input")) != ERR_NOERR) return "Can't open input file " + args.at("input"); |
||||
std::vector<TPoint> points; |
||||
struct Region reg; |
||||
|
||||
{ |
||||
michlib::uint loncol = in.Columns(), latcol = in.Columns(), timecol = in.Columns(); |
||||
michlib::uint lon0col = in.Columns(), lat0col = in.Columns(), dtimecol = in.Columns(); |
||||
for(uint i = 0; i < in.Columns(); i++) |
||||
{ |
||||
if(in.ColumnName(i + 1) == "lon") loncol = i; |
||||
if(in.ColumnName(i + 1) == "lat") latcol = i; |
||||
if(in.ColumnName(i + 1) == "rtime") timecol = i; |
||||
if(in.ColumnName(i + 1) == "lon0") lon0col = i; |
||||
if(in.ColumnName(i + 1) == "lat0") lat0col = i; |
||||
if(in.ColumnName(i + 1) == "time") dtimecol = i; |
||||
} |
||||
if(loncol >= in.Columns()) return "Lon column not found in input file " + args.at("input"); |
||||
if(latcol >= in.Columns()) return "Lat column not found in input file " + args.at("input"); |
||||
if(timecol >= in.Columns()) return "Time column not found in input file " + args.at("input"); |
||||
if(lon0col >= in.Columns()) return "Lon0 column not found in input file " + args.at("input"); |
||||
if(lat0col >= in.Columns()) return "Lat0 column not found in input file " + args.at("input"); |
||||
if(dtimecol >= in.Columns()) return "Age column not found in input file " + args.at("input"); |
||||
|
||||
points.resize(in.Rows()); |
||||
reg.lonb = reg.lone = in[loncol][0]; |
||||
reg.latb = reg.late = in[latcol][0]; |
||||
for(uint i = 0; i < in.Rows(); i++) |
||||
{ |
||||
points[i].p.x = in[loncol][i]; |
||||
points[i].p.y = in[latcol][i]; |
||||
points[i].t = refdate + static_cast<time_t>(michlib::Round(in[timecol][i] * 86400)); |
||||
points[i].lon0 = in[lon0col][i]; |
||||
points[i].lat0 = in[lat0col][i]; |
||||
points[i].dtime = in[dtimecol][i]; |
||||
|
||||
reg.lonb = std::min(reg.lonb, points[i].p.x); |
||||
reg.lone = std::max(reg.lone, points[i].p.x); |
||||
reg.latb = std::min(reg.latb, points[i].p.y); |
||||
reg.late = std::max(reg.late, points[i].p.y); |
||||
} |
||||
std::ranges::sort(points, {}, &TPoint::t); |
||||
} |
||||
in.Close(); |
||||
|
||||
auto resop = ds.Open(args); |
||||
if(resop.Exist()) return "Can't open source: " + resop; |
||||
|
||||
michlib_internal::ParameterListEx pars; |
||||
pars.UsePrefix(""); |
||||
pars.SetParameter("source", args.at("source")); |
||||
pars.SetParameter("history", args.at("_cmdline")); |
||||
|
||||
MString varstring; |
||||
if(args.contains("var")) |
||||
varstring = args.at("var"); |
||||
else |
||||
{ |
||||
if(args.contains("vars")) |
||||
varstring = args.at("vars"); |
||||
else |
||||
{ |
||||
if constexpr(HasDefVars<D>) |
||||
varstring = ds.DefaultVars(); |
||||
else |
||||
return "Variables not specified"; |
||||
} |
||||
} |
||||
|
||||
auto vlist = michlib::Split_on_words(varstring, " ,", 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 " + varstring; |
||||
if(ds.CheckVar(vname) == VarPresence::NONE) return "Variable " + vname + " not exists in this dataset"; |
||||
used.insert(vname); |
||||
} |
||||
} |
||||
|
||||
pars.SetParameter("variables", varstring); |
||||
|
||||
//int compress = 3;
|
||||
//if(args.contains("compress")) compress = args.at("compress").ToInt();
|
||||
|
||||
std::unique_ptr<const BaseParameters> sourcepars; |
||||
if constexpr(ParametersSupported<D>) |
||||
{ |
||||
if constexpr(ParametersRequiredRegion<D>) |
||||
{ |
||||
auto [p, err] = ds.Parameters(pars, args, reg); |
||||
if(err.Exist()) return err; |
||||
sourcepars.reset(p); |
||||
} |
||||
else |
||||
{ |
||||
auto [p, err] = ds.Parameters(pars, args); |
||||
if(err.Exist()) return err; |
||||
sourcepars.reset(p); |
||||
} |
||||
} |
||||
auto p = sourcepars.get(); |
||||
|
||||
MString name = args.contains("out") ? args.at("out") : "out.bin"; |
||||
//MString outfmt = args.contains("format") ? args.at("format") : (GetExt(name) == "nc" ? "nc" : "bin");
|
||||
|
||||
size_t loind = ds.NTimes(), hiind = ds.NTimes(), curind = 0; |
||||
michlib::BFileW fw; |
||||
const MDateTime first = ds.Time(0), last = ds.Time(ds.NTimes() - 1); |
||||
|
||||
std::vector<LinearInterpolator<ReadType<D>>> lo, hi; |
||||
|
||||
{ |
||||
size_t ic = 0; |
||||
fw.Create(name, 6 + vnames.size()); |
||||
fw.SetColumnName(++ic, "lon0"); |
||||
fw.SetColumnName(++ic, "lat0"); |
||||
fw.SetColumnName(++ic, "lon"); |
||||
fw.SetColumnName(++ic, "lat"); |
||||
fw.SetColumnName(++ic, "time"); |
||||
fw.SetColumnName(++ic, "rtime"); |
||||
for(size_t i = 0; i < vnames.size(); i++) fw.SetColumnName(++ic, vnames[i]); |
||||
fw.SetParameters(pars); |
||||
} |
||||
|
||||
for(size_t i = 0; i < points.size(); i++) |
||||
{ |
||||
if(points[i].t < first || points[i].t > last) continue; |
||||
while(!(points[i].t >= ds.Time(curind) && points[i].t <= ds.Time(curind + 1))) curind++; |
||||
if(curind != loind && curind != hiind) |
||||
{ |
||||
loind = curind; |
||||
hiind = curind + 1; |
||||
{ |
||||
auto temp = Read(ds, vnames, p, loind); |
||||
if(temp.size() != vnames.size()) return "Can't read data"; |
||||
//lo.resize(temp.size());
|
||||
for(size_t j = 0; j < temp.size(); j++) lo.emplace_back(std::move(temp[j])); |
||||
} |
||||
{ |
||||
auto temp = Read(ds, vnames, p, hiind); |
||||
if(temp.size() != vnames.size()) return "Can't read data"; |
||||
//hi.resize(temp.size());
|
||||
for(size_t j = 0; j < temp.size(); j++) hi.emplace_back(std::move(temp[j])); |
||||
} |
||||
} |
||||
else if(curind == hiind) |
||||
{ |
||||
loind = curind; |
||||
hiind = curind + 1; |
||||
lo = std::move(hi); |
||||
{ |
||||
auto temp = Read(ds, vnames, p, hiind); |
||||
if(temp.size() != vnames.size()) return "Can't read data"; |
||||
//hi.resize(temp.size());
|
||||
for(size_t j = 0; j < temp.size(); j++) hi.emplace_back(std::move(temp[j])); |
||||
} |
||||
} |
||||
auto step = (ds.Time(hiind) - ds.Time(loind)).S(); |
||||
auto delta = (points[i].t - ds.Time(loind)).S(); |
||||
real trel = delta / step; |
||||
|
||||
fw.Write(points[i].lon0); |
||||
fw.Write(points[i].lat0); |
||||
fw.Write(points[i].p.x); |
||||
fw.Write(points[i].p.y); |
||||
fw.Write(points[i].dtime); |
||||
fw.Write(DeltaDates(refdate, points[i].t)); |
||||
for(size_t iv = 0; iv < lo.size(); iv++) |
||||
{ |
||||
real vlo = lo[iv](points[i].p); |
||||
real vhi = hi[iv](points[i].p); |
||||
fw.Write(vlo + (vhi - vlo) * trel); |
||||
} |
||||
} |
||||
fw.Finalize(); |
||||
fw.Close(); |
||||
return ""; |
||||
}; |
||||
@ -1,10 +0,0 @@
|
||||
#define MICHLIB_NOSOURCE |
||||
#include "actioninterpolate2d.h" |
||||
|
||||
InterpolateMethods::Mode InterpolateMethods::GetMode(const CLArgs& args) |
||||
{ |
||||
MString mode = args.contains("mode") ? args.at("mode") : "vylet_start"; |
||||
if(mode == "vylet_start") return Mode::VYLETSTART; |
||||
if(mode == "vylet_end") return Mode::VYLETEND; |
||||
return Mode::INVALID; |
||||
} |
||||
@ -1,181 +0,0 @@
|
||||
#pragma once |
||||
#include "BFileR.h" |
||||
#include "actiondep.h" |
||||
|
||||
//using michlib::message;
|
||||
|
||||
class InterpolateMethods |
||||
{ |
||||
protected: |
||||
enum class Mode |
||||
{ |
||||
INVALID, |
||||
VYLETSTART, |
||||
VYLETEND |
||||
}; |
||||
struct DataPoint |
||||
{ |
||||
real lon, lat; |
||||
MDateTime t; |
||||
}; |
||||
|
||||
class VyletFile |
||||
{ |
||||
michlib::BFileR fr; |
||||
bool invtime; |
||||
bool usestart; |
||||
MDateTime beg, end, t0; |
||||
michlib::CompiledParser xy2lon, xy2lat; |
||||
real x, y; |
||||
|
||||
public: |
||||
MString Open(const MString& name, Mode mode); |
||||
}; |
||||
|
||||
static Mode GetMode(const CLArgs& args); |
||||
}; |
||||
|
||||
ADD_ACTION(Interpolate, interpolate, CanInterpolate2D<Source>, InterpolateMethods); |
||||
|
||||
template<class D> MString ActionInterpolate::DoAction(const CLArgs& args, D& ds) |
||||
{ |
||||
Mode mode = GetMode(args); |
||||
if(mode == Mode::INVALID) return "Unknown mode"; |
||||
|
||||
if(!args.contains("in")) return "No vylet data file specified"; |
||||
MString in = args.at("in"); |
||||
|
||||
if(!args.contains("out")) return "No output file specified"; |
||||
MString out = args.at("out"); |
||||
|
||||
michlib::BFileR fr; |
||||
if(fr.Open(in) != ERR_NOERR) return "Can't open file " + in; |
||||
|
||||
std::vector<struct DataPoint> datapoints; |
||||
{ |
||||
bool invtime; |
||||
{ |
||||
fr.UsePrefix(""); |
||||
MString method = fr.ParameterSValue("Method", ""); |
||||
if(method == "Bicubic" || method == "BicubicL") |
||||
invtime = false; |
||||
else if(method == "BicubicI" || method == "BicubicIL") |
||||
invtime = true; |
||||
else |
||||
return "Unknown method in the file " + in; |
||||
} |
||||
|
||||
MDateTime beg, end; |
||||
{ |
||||
fr.UsePrefix("Datafile_Info"); |
||||
MString s; |
||||
s = fr.ParameterSValue("BeginDate", ""); |
||||
beg.FromString(s); |
||||
s = fr.ParameterSValue("EndDate", ""); |
||||
end.FromString(s); |
||||
} |
||||
|
||||
MDateTime t0; |
||||
{ |
||||
fr.UsePrefix(""); |
||||
auto tbeg = static_cast<time_t>(fr.ParameterRValue("tbeg", 0.0) * 86400); |
||||
if(invtime) |
||||
{ |
||||
t0 = end; |
||||
t0.AddSeconds(-tbeg); |
||||
} |
||||
else |
||||
{ |
||||
t0 = beg; |
||||
t0.AddSeconds(tbeg); |
||||
} |
||||
} |
||||
|
||||
michlib::CompiledParser xy2lon, xy2lat; |
||||
real x, y; |
||||
michlib::ParserVars pv; |
||||
|
||||
pv["x"] = &x; |
||||
pv["y"] = &y; |
||||
|
||||
fr.UsePrefix("Datafile_Info"); |
||||
if(!ArifmeticCompiler(fr.ParameterSValue("xy2lon", ""), xy2lon, &pv)) return "Can't find xy2lon in the file " + in; |
||||
if(!ArifmeticCompiler(fr.ParameterSValue("xy2lat", ""), xy2lat, &pv)) return "Can't find xy2lat in the file " + in; |
||||
|
||||
auto resop = ds.Open(args); |
||||
if(resop.Exist()) return "Can't open source: " + resop; |
||||
|
||||
for(size_t i = 0; i < fr.Rows(); i++) |
||||
{ |
||||
if(mode == Mode::VYLETSTART) |
||||
{ |
||||
x = fr[0][i]; |
||||
y = fr[1][i]; |
||||
real lon, lat; |
||||
xy2lon.Run(lon); |
||||
xy2lat.Run(lat); |
||||
datapoints.emplace_back(lon, lat, t0); |
||||
} |
||||
if(mode == Mode::VYLETEND) |
||||
{ |
||||
x = fr[2][i]; |
||||
y = fr[3][i]; |
||||
real lon, lat; |
||||
xy2lon.Run(lon); |
||||
xy2lat.Run(lat); |
||||
time_t t = static_cast<time_t>(fr[4][i] * 86400); |
||||
MDateTime time = invtime ? end : beg; |
||||
time.AddSeconds((invtime ? -1 : 1) * t); |
||||
datapoints.emplace_back(lon, lat, time); |
||||
} |
||||
} |
||||
} |
||||
|
||||
if(datapoints.size() == 0) return "No data readed from file " + in; |
||||
|
||||
struct Region reg; |
||||
reg.lonb = datapoints[0].lon; |
||||
reg.lone = datapoints[0].lon; |
||||
reg.latb = datapoints[0].lat; |
||||
reg.late = datapoints[0].lat; |
||||
|
||||
for(const auto& p: datapoints) |
||||
{ |
||||
if(p.lon < reg.lonb) reg.lonb = p.lon; |
||||
if(p.lat < reg.latb) reg.latb = p.lat; |
||||
if(p.lon > reg.lone) reg.lone = p.lon; |
||||
if(p.lat > reg.late) reg.late = p.lat; |
||||
} |
||||
|
||||
return ""; |
||||
/*
|
||||
|
||||
if(!args.contains("var")) return "Variable not specified"; |
||||
MString vname = args.at("var"); |
||||
if(!ds.CheckVar(vname)) return "Variable " + vname + " not exists in this dataset"; |
||||
pars.SetParameter("variable", vname); |
||||
|
||||
std::unique_ptr<const BaseParameters> sourcepars; |
||||
if constexpr(ParametersSupported<D>) |
||||
{ |
||||
if constexpr(ParametersRequiredRegion<D>) |
||||
{ |
||||
auto [p, err] = ds.Parameters(pars, args, reg); |
||||
if(err.Exist()) return err; |
||||
sourcepars.reset(p); |
||||
} |
||||
else |
||||
{ |
||||
auto [p, err] = ds.Parameters(pars, args); |
||||
if(err.Exist()) return err; |
||||
sourcepars.reset(p); |
||||
} |
||||
} |
||||
auto p = sourcepars.get(); |
||||
|
||||
auto data = Read(ds, vname, p, tindexes); |
||||
if(!data) return "Can't read data"; |
||||
|
||||
return ""; |
||||
*/ |
||||
}; |
||||
@ -0,0 +1,22 @@
|
||||
#pragma once |
||||
#include "actiondep.h" |
||||
#include "merrors.h" |
||||
|
||||
using michlib::message; |
||||
using michlib::Error; |
||||
|
||||
template<class T> |
||||
concept MirrorSupported = requires(T t, const CLArgs& args) { |
||||
{ |
||||
t.Mirror(args) |
||||
} -> std::convertible_to<Error>; |
||||
}; |
||||
|
||||
ADD_ACTION(Mirror, mirror, MirrorSupported<Source>); |
||||
|
||||
template<class D> MString ActionMirror::DoAction(const CLArgs& args, D& data) |
||||
{ |
||||
auto res = data.Mirror(args); |
||||
if(!res) return "Mirroring failed"; |
||||
return ""; |
||||
}; |
||||
@ -0,0 +1,447 @@
|
||||
#pragma once |
||||
#include "GPL.h" |
||||
#include "mirrorfuncs.h" |
||||
#include <functional> |
||||
#include <libpq-fe.h> |
||||
#include <optional> |
||||
#include <sqlite3.h> |
||||
#include <time.h> |
||||
#include <variant> |
||||
|
||||
using michlib::GPL; |
||||
using michlib::int1; |
||||
using michlib::int4; |
||||
using michlib::int_cast; |
||||
using michlib::MString; |
||||
using michlib::pointer_cast; |
||||
|
||||
class SQLiteConnection |
||||
{ |
||||
public: |
||||
using DBType = sqlite3*; |
||||
using FuncType = std::function<void(DBType)>; |
||||
|
||||
private: |
||||
static DBType db; |
||||
static size_t count; |
||||
static std::vector<FuncType> destructs; |
||||
|
||||
public: |
||||
SQLiteConnection() |
||||
{ |
||||
count++; |
||||
if(db == nullptr) |
||||
{ |
||||
MString oldprefix = GPL.UsePrefix("SQLITE"); |
||||
MString name = GPL.ParameterSValue("db", ""); |
||||
GPL.UsePrefix(oldprefix); |
||||
|
||||
auto ret = sqlite3_open(name.Buf(), &db); |
||||
if(ret != SQLITE_OK) |
||||
{ |
||||
sqlite3_close(db); |
||||
db = nullptr; |
||||
} |
||||
} |
||||
} |
||||
|
||||
SQLiteConnection([[maybe_unused]] const SQLiteConnection& sq): SQLiteConnection() {} |
||||
SQLiteConnection(SQLiteConnection&&) = default; |
||||
|
||||
SQLiteConnection& operator=([[maybe_unused]] const SQLiteConnection& sq) |
||||
{ |
||||
*this = {}; |
||||
return *this; |
||||
} |
||||
SQLiteConnection& operator=(SQLiteConnection&&) = default; |
||||
|
||||
~SQLiteConnection() |
||||
{ |
||||
if(count == 0) michlib::errmessage("Destructor of SQLiteConnection called on count==0"); |
||||
if(count > 1) |
||||
count--; |
||||
else |
||||
{ |
||||
count = 0; |
||||
if(db != nullptr) |
||||
{ |
||||
for(const auto& f: destructs) f(db); |
||||
sqlite3_close(db); |
||||
} |
||||
db = nullptr; |
||||
} |
||||
} |
||||
|
||||
static void AddDestructor(FuncType&& f) { destructs.emplace_back(std::move(f)); } |
||||
|
||||
operator DBType() const { return db; } |
||||
|
||||
static DBType GetDB() { return db; } |
||||
|
||||
explicit operator bool() const { return db != nullptr; } |
||||
}; |
||||
|
||||
class PostgreSQLConnection |
||||
{ |
||||
public: |
||||
using DBType = PGconn*; |
||||
using FuncType = std::function<void(DBType)>; |
||||
|
||||
private: |
||||
static DBType conn; |
||||
static size_t count; |
||||
static std::vector<FuncType> destructs; |
||||
|
||||
public: |
||||
PostgreSQLConnection() |
||||
{ |
||||
count++; |
||||
if(conn == nullptr) |
||||
{ |
||||
MString oldprefix = GPL.UsePrefix("POSTGRES"); |
||||
MString name = GPL.ParameterSValue("connection", ""); |
||||
GPL.UsePrefix(oldprefix); |
||||
|
||||
conn = PQconnectdb(name.Buf()); |
||||
if(PQstatus(conn) != CONNECTION_OK) |
||||
{ |
||||
michlib::errmessage(PQerrorMessage(conn)); |
||||
PQfinish(conn); |
||||
conn = nullptr; |
||||
} |
||||
} |
||||
} |
||||
|
||||
PostgreSQLConnection([[maybe_unused]] const PostgreSQLConnection& pq): PostgreSQLConnection() {} |
||||
PostgreSQLConnection(PostgreSQLConnection&&) = default; |
||||
|
||||
PostgreSQLConnection& operator=([[maybe_unused]] const PostgreSQLConnection& pq) |
||||
{ |
||||
*this = {}; |
||||
return *this; |
||||
} |
||||
PostgreSQLConnection& operator=(PostgreSQLConnection&&) = default; |
||||
|
||||
~PostgreSQLConnection() |
||||
{ |
||||
if(count == 0) michlib::errmessage("Destructor of PostgreSQLConnection called on count==0"); |
||||
if(count > 1) |
||||
count--; |
||||
else |
||||
{ |
||||
count = 0; |
||||
if(conn != nullptr) |
||||
{ |
||||
for(const auto& f: destructs) f(conn); |
||||
PQfinish(conn); |
||||
} |
||||
conn = nullptr; |
||||
} |
||||
} |
||||
|
||||
static void AddDestructor(FuncType&& f) { destructs.emplace_back(std::move(f)); } |
||||
|
||||
operator DBType() const { return conn; } |
||||
|
||||
static DBType GetDB() { return conn; } |
||||
|
||||
explicit operator bool() const { return conn != nullptr; } |
||||
}; |
||||
|
||||
class GenericCache |
||||
{ |
||||
public: |
||||
virtual bool Put(const MString& key, const MString& value, size_t ttl) const = 0; |
||||
virtual std::pair<MString, bool> Get(const MString& key) const = 0; |
||||
virtual ~GenericCache() {} |
||||
}; |
||||
|
||||
class FakeCache: public GenericCache |
||||
{ |
||||
public: |
||||
virtual bool Put([[maybe_unused]] const MString& key, [[maybe_unused]] const MString& value, [[maybe_unused]] size_t ttl) const override { return false; } |
||||
virtual std::pair<MString, bool> Get([[maybe_unused]] const MString& key) const override { return {"", false}; } |
||||
virtual ~FakeCache() override {} |
||||
}; |
||||
|
||||
class SQLiteCache: public GenericCache |
||||
{ |
||||
static bool regdest; |
||||
|
||||
SQLiteConnection db; |
||||
|
||||
public: |
||||
bool Init() |
||||
{ |
||||
if(!db) return false; |
||||
|
||||
if(!regdest) |
||||
{ |
||||
// Create table
|
||||
sqlite3_stmt* sqst; |
||||
int i; |
||||
|
||||
i = sqlite3_prepare_v2(db, |
||||
"CREATE TABLE IF NOT EXISTS `cache`('key' TEXT PRIMARY KEY ON CONFLICT REPLACE NOT NULL ON CONFLICT FAIL, 'value' BLOB NOT NULL ON CONFLICT FAIL, " |
||||
"'exptime' INTEGER NOT NULL ON CONFLICT FAIL) WITHOUT ROWID, STRICT;", |
||||
-1, &sqst, 0); |
||||
i = sqlite3_step(sqst); |
||||
if(i != SQLITE_DONE) |
||||
{ |
||||
sqlite3_finalize(sqst); |
||||
return false; |
||||
} |
||||
sqlite3_finalize(sqst); |
||||
sqlite3_busy_timeout(db, 1000); |
||||
|
||||
db.AddDestructor( |
||||
[](SQLiteConnection::DBType db) |
||||
{ |
||||
sqlite3_stmt* sqst = nullptr; |
||||
int i = SQLITE_OK; |
||||
|
||||
if(i == SQLITE_OK) i = sqlite3_prepare_v2(db, "DELETE from `cache` WHERE exptime<?1;", -1, &sqst, 0); |
||||
if(i == SQLITE_OK) i = sqlite3_bind_int64(sqst, 1, time(nullptr)); |
||||
if(i == SQLITE_OK) i = sqlite3_step(sqst); |
||||
sqlite3_finalize(sqst); |
||||
}); |
||||
regdest = true; |
||||
} |
||||
return true; |
||||
} |
||||
|
||||
virtual bool Put(const MString& key, const MString& value, size_t ttl) const override |
||||
{ |
||||
if(!*this) return false; |
||||
sqlite3_stmt* sqst = nullptr; |
||||
int i = SQLITE_OK; |
||||
|
||||
if(i == SQLITE_OK) i = sqlite3_prepare_v2(db, "INSERT OR REPLACE into `cache` VALUES(?1,?2,?3);", -1, &sqst, 0); |
||||
if(i == SQLITE_OK) i = sqlite3_bind_text(sqst, 1, key.Buf(), -1, SQLITE_STATIC); |
||||
if(i == SQLITE_OK) i = sqlite3_bind_blob64(sqst, 2, value.Buf(), value.Len(), SQLITE_STATIC); |
||||
if(i == SQLITE_OK) i = sqlite3_bind_int64(sqst, 3, time(nullptr) + ttl); |
||||
if(i == SQLITE_OK) i = sqlite3_step(sqst); |
||||
sqlite3_finalize(sqst); |
||||
|
||||
return i == SQLITE_OK; |
||||
} |
||||
|
||||
virtual std::pair<MString, bool> Get(const MString& key) const override |
||||
{ |
||||
if(!*this) return {"", false}; |
||||
|
||||
sqlite3_stmt* sqst = nullptr; |
||||
int i = SQLITE_OK; |
||||
|
||||
if(i == SQLITE_OK) i = sqlite3_prepare_v2(db, "SELECT value from `cache` WHERE key=?1 AND exptime>?2;", -1, &sqst, 0); |
||||
if(i == SQLITE_OK) i = sqlite3_bind_text(sqst, 1, key.Buf(), -1, SQLITE_STATIC); |
||||
if(i == SQLITE_OK) i = sqlite3_bind_int64(sqst, 2, time(nullptr)); |
||||
if(i == SQLITE_OK) i = sqlite3_step(sqst); |
||||
if(i == SQLITE_ROW) |
||||
{ |
||||
auto p = sqlite3_column_blob(sqst, 0); |
||||
auto sz = sqlite3_column_bytes(sqst, 0); |
||||
if(p != nullptr) |
||||
{ |
||||
MString out(p, sz); |
||||
sqlite3_finalize(sqst); |
||||
return {std::move(out), true}; |
||||
} |
||||
} |
||||
sqlite3_finalize(sqst); |
||||
return {"", false}; |
||||
} |
||||
|
||||
virtual ~SQLiteCache() override = default; |
||||
|
||||
explicit operator bool() const { return db != nullptr; } |
||||
}; |
||||
|
||||
class PostgreSQLHelpers |
||||
{ |
||||
static constexpr time_t postgresepoch = 946648800; |
||||
|
||||
class PGResultRAIIDT |
||||
{ |
||||
public: |
||||
// TODO: make static
|
||||
void operator()(PGresult* res) { PQclear(res); } |
||||
}; |
||||
|
||||
protected: |
||||
PostgreSQLConnection conn; |
||||
|
||||
// Convert Postgres binary representation of timestamp to Unix epoch seconds. Microseconds ignored
|
||||
static time_t raw2epoch(time_t raw) { return Invert(raw) / 1000000 + postgresepoch; } |
||||
|
||||
// Convert Unix epoch time to Postres binary representation
|
||||
static time_t epoch2raw(time_t epoch) { return Invert((epoch - postgresepoch) * 1000000); } |
||||
|
||||
class PGresultRAII: public std::unique_ptr<PGresult, PGResultRAIIDT> |
||||
{ |
||||
public: |
||||
PGresultRAII() = default; |
||||
PGresultRAII(PGresult* res): std::unique_ptr<PGresult, PGResultRAIIDT>(res) {} |
||||
|
||||
operator PGresult*() const { return get(); } |
||||
}; |
||||
|
||||
bool CheckCon() const |
||||
{ |
||||
if(!*this) return false; |
||||
if(PQstatus(conn) == CONNECTION_OK) return true; |
||||
PQreset(conn); |
||||
return PQstatus(conn) == CONNECTION_OK; |
||||
} |
||||
|
||||
template<class D> static D Invert(D d) |
||||
{ |
||||
using michlib::int1; |
||||
D out; |
||||
Invert(&d, &out, sizeof(D)); |
||||
return out; |
||||
} |
||||
|
||||
static void Invert(const void* src, void* dst, size_t sz) |
||||
{ |
||||
if(sz == 0) return; |
||||
const int1* pin = pointer_cast<const int1*>(src); |
||||
int1* pout = pointer_cast<int1*>(dst); |
||||
for(size_t i = 0; i < sz; i++) pout[sz - i - 1] = pin[i]; |
||||
} |
||||
|
||||
public: |
||||
explicit operator bool() const { return conn != nullptr; } |
||||
}; |
||||
|
||||
class PostgreSQLCache: public GenericCache, public PostgreSQLHelpers |
||||
{ |
||||
static bool regdest; |
||||
|
||||
public: |
||||
bool Init() |
||||
{ |
||||
if(!conn) return false; |
||||
|
||||
if(!regdest) |
||||
{ |
||||
PGresultRAII res; |
||||
|
||||
// Create table
|
||||
res = PQexec(conn, "SET client_min_messages=WARNING;"); |
||||
res = PQexec(conn, "CREATE TABLE IF NOT EXISTS cache(key TEXT PRIMARY KEY NOT NULL, value BYTEA, exptime TIMESTAMP(0) NOT NULL);"); |
||||
if(PQresultStatus(res) != PGRES_COMMAND_OK) |
||||
{ |
||||
michlib::errmessage(PQresStatus(PQresultStatus(res))); |
||||
michlib::errmessage(PQerrorMessage(conn)); |
||||
} |
||||
|
||||
res = PQexec(conn, "SET client_min_messages=NOTICE;"); |
||||
|
||||
conn.AddDestructor([](PostgreSQLConnection::DBType conn) { PGresultRAII res = PQexec(conn, "DELETE FROM cache WHERE exptime<localtimestamp;"); }); |
||||
|
||||
regdest = true; |
||||
} |
||||
return true; |
||||
} |
||||
|
||||
virtual bool Put(const MString& key, const MString& value, size_t ttl) const override |
||||
{ |
||||
if(!CheckCon()) return false; |
||||
|
||||
auto interval = michlib::int_cast<michlib::int8>(ttl); |
||||
michlib::int8 rinterval = Invert(interval); |
||||
const char* params[] = {key.Buf(), value.Buf(), pointer_cast<const char*>(&rinterval)}; |
||||
int plens[] = {int_cast<int>(key.Len()), int_cast<int>(value.Len()), sizeof(rinterval)}; |
||||
int pfor[] = {0, 1, 1}; |
||||
|
||||
PGresultRAII res = PQexecParams(conn, |
||||
"INSERT INTO cache(key,value,exptime) VALUES($1,$2,localtimestamp + ($3::bigint ||' seconds')::interval)" |
||||
"ON CONFLICT(key) DO UPDATE SET value=EXCLUDED.value, exptime=EXCLUDED.exptime;", |
||||
3, nullptr, params, plens, pfor, 1); |
||||
if(PQresultStatus(res) != PGRES_COMMAND_OK) |
||||
{ |
||||
michlib::errmessage(PQresStatus(PQresultStatus(res))); |
||||
michlib::errmessage(PQerrorMessage(conn)); |
||||
return false; |
||||
} |
||||
return true; |
||||
} |
||||
|
||||
virtual std::pair<MString, bool> Get(const MString& key) const override |
||||
{ |
||||
if(!CheckCon()) return {"", false}; |
||||
|
||||
const char* params[] = {key.Buf()}; |
||||
int plens[] = {int_cast<int>(key.Len())}; |
||||
int pfor[] = {0}; |
||||
|
||||
PGresultRAII res = PQexecParams(conn, "SELECT value from cache WHERE key=$1::text AND exptime>localtimestamp;", 1, nullptr, params, plens, pfor, 1); |
||||
if(PQresultStatus(res) != PGRES_TUPLES_OK) |
||||
{ |
||||
michlib::errmessage(PQresStatus(PQresultStatus(res))); |
||||
michlib::errmessage(PQerrorMessage(conn)); |
||||
return {"", false}; |
||||
} |
||||
else if(PQntuples(res) == 0) |
||||
return {"", false}; |
||||
|
||||
MString val(PQgetvalue(res, 0, 0), PQgetlength(res, 0, 0)); |
||||
return {std::move(val), true}; |
||||
} |
||||
|
||||
virtual ~PostgreSQLCache() override = default; |
||||
}; |
||||
|
||||
inline GenericCache* CreateCache(const MString& cachedesc) |
||||
{ |
||||
auto i = cachedesc.GetPos(':'); |
||||
auto name = i == 0 ? cachedesc : cachedesc.SubStr(1, i - 1); |
||||
auto par = i == 0 ? "" : cachedesc.SubStr(i + 1, cachedesc.Len() - i); |
||||
|
||||
if(name == "no") return new FakeCache; |
||||
|
||||
if(name == "sqlite") |
||||
{ |
||||
auto ret = new SQLiteCache; |
||||
ret->Init(); |
||||
if(*ret) return ret; |
||||
delete ret; |
||||
} |
||||
|
||||
if(name == "postgre" || name == "postgres" || name == "postgresql") |
||||
{ |
||||
auto ret = new PostgreSQLCache; |
||||
ret->Init(); |
||||
if(*ret) return ret; |
||||
delete ret; |
||||
} |
||||
|
||||
return nullptr; |
||||
} |
||||
|
||||
class FileInfoCache: public PostgreSQLHelpers |
||||
{ |
||||
public: |
||||
using DataType = std::optional<MString>; |
||||
using CallbackType = std::function<DataType(const MString&)>; |
||||
|
||||
private: |
||||
static bool regdest; |
||||
CallbackType readfunc; |
||||
MString dir; |
||||
int4 dirid; |
||||
|
||||
FileInfoCache() = delete; |
||||
|
||||
CallbackType::result_type GetData(const MString& fname) const { return readfunc(dir + "/" + fname); } |
||||
|
||||
void GetDirId(); |
||||
|
||||
public: |
||||
FileInfoCache(CallbackType&& readfunc_, const MString& dir_); |
||||
|
||||
Error UpdateCache(bool force = false) const; |
||||
|
||||
DataType GetInfo(const MString& name) const; |
||||
}; |
||||
@ -0,0 +1,54 @@
|
||||
#pragma once |
||||
#include "cache.h" |
||||
#include "curlfuncs.h" |
||||
#include "merrors.h" |
||||
#include <json/json.h> |
||||
|
||||
using michlib::Error; |
||||
using michlib::RetVal; |
||||
|
||||
class CopernicusCatalog |
||||
{ |
||||
static const MString caturl; |
||||
|
||||
std::unique_ptr<GenericCache> cache; |
||||
CURLRAII chandle; |
||||
Json::Value catalog; |
||||
|
||||
// Download catalog
|
||||
Error GetCatalog(); |
||||
|
||||
// Asset url from dataset
|
||||
RetVal<MString> AssetURL(const MString& prod, const MString& dataset, const MString& asset) const; |
||||
|
||||
public: |
||||
CopernicusCatalog(); |
||||
|
||||
// Download JSON from url
|
||||
RetVal<Json::Value> GetJSON(const MString& url) const; |
||||
|
||||
// List of products
|
||||
RetVal<std::vector<MString>> ProductList() const; |
||||
|
||||
// List of datasets in product
|
||||
RetVal<std::vector<MString>> DatasetList(const MString& prod) const; |
||||
|
||||
// URL of product
|
||||
RetVal<MString> ProductURL(const MString& prod) const; |
||||
|
||||
// URL of dataset
|
||||
RetVal<MString> DatasetURL(const MString& prod, const MString& dataset) const; |
||||
|
||||
// URL of native data (files) in dataset
|
||||
RetVal<MString> DatasetNativeURL(const MString& prod, const MString& dataset) const { return AssetURL(prod, dataset, "native"); } |
||||
|
||||
// URL of timechuncked data (files) in dataset
|
||||
RetVal<MString> DatasetTimeURL(const MString& prod, const MString& dataset) const { return AssetURL(prod, dataset, "timeChunked"); } |
||||
|
||||
// URL of geochuncked data (files) in dataset
|
||||
RetVal<MString> DatasetGeoURL(const MString& prod, const MString& dataset) const { return AssetURL(prod, dataset, "geoChunked"); } |
||||
|
||||
bool Valid() const { return catalog.isObject(); } |
||||
|
||||
explicit operator bool() const { return Valid(); } |
||||
}; |
||||
@ -0,0 +1,36 @@
|
||||
#pragma once |
||||
#include "MString.h" |
||||
#include <curl/curl.h> |
||||
#include <memory> |
||||
|
||||
using michlib::MString; |
||||
|
||||
class CURLRAIIDT |
||||
{ |
||||
public: |
||||
// TODO: make static
|
||||
void operator()(CURL* c) { curl_easy_cleanup(c); } |
||||
}; |
||||
|
||||
class CURLRAII: public std::unique_ptr<CURL, CURLRAIIDT> |
||||
{ |
||||
char err[CURL_ERROR_SIZE]; |
||||
|
||||
public: |
||||
CURLRAII() |
||||
{ |
||||
reset(curl_easy_init()); |
||||
curl_easy_setopt(*this, CURLOPT_ERRORBUFFER, err); |
||||
} |
||||
operator CURL*() const { return get(); } |
||||
const char* Err() const { return err; } |
||||
}; |
||||
|
||||
// Curl writeback function, write to MString
|
||||
size_t Write2String(char* ptr, size_t size, size_t n, void* data); |
||||
|
||||
// Curl writeback function, write to file descriptor
|
||||
size_t Write2File(char* ptr, size_t size, size_t n, void* data); |
||||
|
||||
// Get content of url to MString
|
||||
std::pair<MString, CURLcode> GetUrl(const CURLRAII& chandle, const MString& url); |
||||
@ -0,0 +1,51 @@
|
||||
#pragma once |
||||
#include "basedata.h" |
||||
#include "comdefs.h" |
||||
|
||||
using michlib::real; |
||||
|
||||
template<class D> class LinearInterpolator: public D |
||||
{ |
||||
public: |
||||
LinearInterpolator(D&& d): D(std::move(d)) {} |
||||
|
||||
real operator()(const struct Point2D& in) const |
||||
{ |
||||
const auto gp = D::GridPos(in.x, in.y); |
||||
if(!gp.Valid()) return NAN; |
||||
|
||||
if(gp.x == 0.0 && gp.y == 0.0) return D::V(gp.ix, gp.iy); |
||||
|
||||
real v00; |
||||
real v10; |
||||
real v01; |
||||
real v11; |
||||
bool isfill = false; |
||||
size_t fx, fy; |
||||
|
||||
// Count fills
|
||||
for(size_t ix = 0; ix <= 1; ix++) |
||||
for(size_t iy = 0; iy <= 1; iy++) |
||||
{ |
||||
if(isfill && D::IsFill(gp.ix + ix, gp.iy + iy)) return NAN; |
||||
if(D::IsFill(gp.ix + ix, gp.iy + iy)) |
||||
{ |
||||
fx = ix; |
||||
fy = iy; |
||||
isfill = true; |
||||
} |
||||
} |
||||
|
||||
v00 = D::V(gp.ix, gp.iy); |
||||
v10 = D::V(gp.ix + 1, gp.iy); |
||||
v01 = D::V(gp.ix, gp.iy + 1); |
||||
v11 = D::V(gp.ix + 1, gp.iy + 1); |
||||
|
||||
if(isfill && fx == 0 && fy == 0) v00 = (v10 + v01 + v11) / 3.0; |
||||
if(isfill && fx == 1 && fy == 0) v10 = (v00 + v01 + v11) / 3.0; |
||||
if(isfill && fx == 0 && fy == 1) v01 = (v10 + v00 + v11) / 3.0; |
||||
if(isfill && fx == 1 && fy == 1) v11 = (v10 + v01 + v00) / 3.0; |
||||
|
||||
return v00 + gp.y * (v01 - v00) + gp.x * ((v10 - v00) + gp.y * (v00 - v01 + v11 - v10)); |
||||
}; |
||||
}; |
||||
@ -0,0 +1,207 @@
|
||||
#pragma once |
||||
#include "gsw.h" |
||||
#include "ncfuncs.h" |
||||
#include "nczarr.h" |
||||
#include "simple2ddata.h" |
||||
#include <memory> |
||||
|
||||
using michlib::Ceil; |
||||
using michlib::DetGeoDomain; |
||||
using michlib::Floor; |
||||
using michlib::GPL; |
||||
using michlib::int2; |
||||
|
||||
class LayeredDataZ: public NCFuncs |
||||
{ |
||||
public: |
||||
using Data = Simple2DData; |
||||
|
||||
private: |
||||
class NC: public NCZarr |
||||
{ |
||||
std::vector<MDateTime> times; |
||||
|
||||
public: |
||||
Error ReadTimes(const MString& tname) |
||||
{ |
||||
static const MString pref = "LayeredDataZ::NC::ReadTimes"; |
||||
if(!*this) return Error(pref, "Dataset not open"); |
||||
std::vector<double> time; |
||||
|
||||
{ |
||||
auto ret = Read(tname, time); |
||||
if(!ret) return ret.Add(pref, "Can't read time"); |
||||
} |
||||
|
||||
MDateTime refdate; |
||||
time_t step = 0; |
||||
{ |
||||
auto units = AttString(tname, "units"); |
||||
if(!units.Exist()) return Error(pref, "Can't read refdate"); |
||||
auto [rd, st, suc] = Refdate(units); |
||||
if(!suc) return Error(pref, "Can't parse " + units + " to refdate"); |
||||
if(st == 0) return Error(pref, "Can't get timestep from string " + units); |
||||
refdate = rd; |
||||
step = st; |
||||
} |
||||
|
||||
times.resize(time.size()); |
||||
for(size_t i = 0; i < time.size(); i++) times[i] = refdate + static_cast<time_t>(time[i] * step); |
||||
return Error(); |
||||
} |
||||
|
||||
MDateTime Begin() const { return times.front(); } |
||||
MDateTime End() const { return times.back(); } |
||||
|
||||
const std::vector<MDateTime>& Times() const { return times; } |
||||
|
||||
size_t Index(MDateTime tm) const |
||||
{ |
||||
if(tm < Begin() || tm > End()) return 0; |
||||
size_t b = 0, e = times.size() - 1; |
||||
if(tm == times[b]) return b + 1; |
||||
if(tm == times[e]) return e + 1; |
||||
while(e - b > 1) |
||||
{ |
||||
size_t c = (e + b) / 2; |
||||
if(tm == times[c]) return c + 1; |
||||
if(tm > times[c]) |
||||
b = c; |
||||
else |
||||
e = c; |
||||
} |
||||
return 0; |
||||
} |
||||
}; |
||||
std::vector<NC> nc; |
||||
std::vector<real> depths; |
||||
bool depthinv; |
||||
std::vector<MDateTime> times; |
||||
struct DimNames dname; |
||||
real lonb, latb, lone, late; |
||||
real lonstep, latstep; |
||||
MString title; |
||||
|
||||
class EnvVar |
||||
{ |
||||
MString name, oldvalue; |
||||
bool activated, saved; |
||||
|
||||
public: |
||||
EnvVar(): activated(false) {} |
||||
~EnvVar() { Deactivate(); } |
||||
|
||||
void Activate(const MString& var, const MString& val) |
||||
{ |
||||
if(activated) Deactivate(); |
||||
name = var; |
||||
char* curval = getenv(name.Buf()); |
||||
if(nullptr == curval) |
||||
saved = false; |
||||
else |
||||
{ |
||||
oldvalue = curval; |
||||
saved = true; |
||||
} |
||||
setenv(name.Buf(), val.Buf(), 1); |
||||
} |
||||
void Deactivate() |
||||
{ |
||||
if(!activated) return; |
||||
if(saved) |
||||
setenv(name.Buf(), oldvalue.Buf(), 1); |
||||
else |
||||
unsetenv(name.Buf()); |
||||
activated = false; |
||||
} |
||||
}; |
||||
|
||||
EnvVar proxy; |
||||
|
||||
protected: |
||||
struct Parameters: public BaseParameters |
||||
{ |
||||
size_t xb, yb, xe, ye, layer; |
||||
virtual ~Parameters() override = default; |
||||
}; |
||||
|
||||
// TODO: RetVal
|
||||
MString Open(const MString& dataset); |
||||
|
||||
void SetTitle(const MString& newtitle) { title = newtitle; } |
||||
|
||||
public: |
||||
MString Info() const; |
||||
|
||||
std::pair<const BaseParameters*, MString> Parameters(michlib_internal::ParameterListEx& pars, const CLArgs& args, const struct Region& reg) 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; } |
||||
|
||||
explicit operator bool() const { return nc.size() > 0; } |
||||
|
||||
real Depth(size_t l) const { return isOk() ? depths[l] : -1000.0; } |
||||
|
||||
real Depth(const BaseParameters* ip) const { return Depth(dynamic_cast<const struct Parameters*>(ip)->layer); } |
||||
|
||||
real Lon(size_t ix) const { return isOk() ? (lonb + ix * lonstep) : -1000.0; } |
||||
|
||||
real Lat(size_t iy) const { return isOk() ? (latb + iy * latstep) : -1000.0; } |
||||
|
||||
size_t NDepths() const { return depths.size(); } |
||||
|
||||
size_t NTimes() const { return times.size(); } |
||||
|
||||
MDateTime Time(size_t i) const |
||||
{ |
||||
if(!isOk() || i >= times.size()) return MDateTime(); |
||||
return times[i]; |
||||
} |
||||
|
||||
time_t Timestep() const { return isOk() ? (times[1] - times[0]).Seconds() : 0; } |
||||
|
||||
MString Title() const { return title; } |
||||
|
||||
MString Dump(const struct Parameters* ppar) const |
||||
{ |
||||
// clang-format off
|
||||
return |
||||
"Current settings:\n" + MString() + |
||||
" Longitudes: from " + Lon(ppar->xb) + " (" + ppar->xb + ") to "+ Lon(ppar->xe) + " (" + ppar->xe + ")\n" + |
||||
" Latitudes: from " + Lat(ppar->yb) + " (" + ppar->yb + ") to "+ Lat(ppar->ye) + " (" + ppar->ye + ")\n" + |
||||
" Depth: layer " + ppar->layer + ", depth " + Depth(ppar->layer) + " m\n"; |
||||
// clang-format on
|
||||
} |
||||
|
||||
VarPresence CheckVar(const MString& vname) const |
||||
{ |
||||
return NCFuncs::CheckVar(vname, [this](const MString& vn) { return HaveVar(vn); }); |
||||
} |
||||
|
||||
private: |
||||
Data ReadVarRaw(const NC& f, const MString& name, size_t i, bool nodepth, const struct Parameters* p) const; |
||||
|
||||
bool HaveVar(const MString& vname) const |
||||
{ |
||||
for(size_t i = 0; i < nc.size(); i++) |
||||
if(NCFuncs::HaveVar(nc[i], vname)) return true; |
||||
return false; |
||||
} |
||||
|
||||
std::tuple<MString, size_t, size_t> VarNameLoc(const MString vname, MDateTime tm) const |
||||
{ |
||||
for(size_t i = 0; i < nc.size(); i++) |
||||
{ |
||||
auto tind = nc[i].Index(tm); |
||||
if(tind == 0) continue; |
||||
for(const auto& v: nc[i].VarNames()) |
||||
{ |
||||
auto stname = nc[i].AttString(v, "standard_name"); |
||||
if(!stname.Exist()) continue; |
||||
if(StName2Name(stname) == vname) return {v, i, tind - 1}; |
||||
} |
||||
} |
||||
return {"", 0, 0}; |
||||
} |
||||
}; |
||||
@ -0,0 +1,64 @@
|
||||
#pragma once |
||||
#include "curlfuncs.h" |
||||
#include "mdatetime.h" |
||||
#include <dirent.h> |
||||
#include <fcntl.h> |
||||
#include <sys/stat.h> |
||||
#include <sys/types.h> |
||||
#include <vector> |
||||
|
||||
using michlib::Error; |
||||
using michlib::MDateTime; |
||||
using michlib::RetVal; |
||||
|
||||
class DIRRAIIDT |
||||
{ |
||||
public: |
||||
// TODO: make static
|
||||
void operator()(DIR* d) { closedir(d); } |
||||
}; |
||||
|
||||
class DIRRAII: public std::unique_ptr<DIR, DIRRAIIDT> |
||||
{ |
||||
public: |
||||
operator DIR*() const { return get(); } |
||||
}; |
||||
|
||||
struct FileInfo |
||||
{ |
||||
MString url; |
||||
MString name; |
||||
MDateTime mtime; |
||||
size_t size; |
||||
}; |
||||
|
||||
// Remove last element from path
|
||||
inline MString DirName(const MString& name) |
||||
{ |
||||
auto p = name.GetPos('/', false); |
||||
if(p == 0) return name; |
||||
return name.SubStr(1, p - 1); |
||||
} |
||||
|
||||
// Get last element from path
|
||||
inline MString FileName(const MString& name) |
||||
{ |
||||
auto p = name.GetPos('/', false); |
||||
if(p == 0) return name; |
||||
return name.SubStr(p + 1, name.Len() - p); |
||||
} |
||||
|
||||
// Check and, if necessary, create the path to the file
|
||||
bool MakePath(const MString& dname); |
||||
|
||||
// Get local file list
|
||||
RetVal<std::vector<struct FileInfo>> ReadLocalFileList(const MString& dir, const bool nofollow = true, const MString& path = ""); |
||||
|
||||
// Download file to the local mirror
|
||||
Error DownloadFile(const CURLRAII& chandle, const struct FileInfo& rinfo, const MString& root); |
||||
|
||||
// Remove file from the local mirror
|
||||
Error RemoveFile(const struct FileInfo& linfo); |
||||
|
||||
// Updare file in the local mirror
|
||||
Error UpdateFile(const CURLRAII& chandle, const struct FileInfo& rinfo, const struct FileInfo& linfo, const MString& root); |
||||
@ -0,0 +1,87 @@
|
||||
#pragma once |
||||
#include "nczarrcommon.h" |
||||
#include <netcdf.h> |
||||
#include <variant> |
||||
|
||||
class NCSimpleTypes: public NcZarrTypes |
||||
{ |
||||
protected: |
||||
template<class VType> class ReadedData |
||||
{ |
||||
using Vec = std::vector<size_t>; |
||||
|
||||
private: |
||||
std::unique_ptr<VType[]> data; |
||||
|
||||
public: |
||||
ReadedData() = default; |
||||
|
||||
ReadedData(std::unique_ptr<VType[]>&& d): data(std::move(d)) {} |
||||
|
||||
VType operator()(size_t lini) const { return data[lini]; } |
||||
}; |
||||
}; |
||||
|
||||
class NCSimpleFunctions: public NCSimpleTypes |
||||
{ |
||||
int ncid; |
||||
|
||||
RetVal<std::vector<Attribute>> ReadAtts(int vid) const; |
||||
|
||||
protected: |
||||
NCSimpleFunctions(): ncid(0) {} |
||||
|
||||
template<class VType> RetVal<ReadedData<VType>> Read(const MString& var, const size_t* start, const size_t* count) const |
||||
{ |
||||
static const MString pref = "NCSimpleFunctions::Read"; |
||||
|
||||
size_t ind = FindInd(var, vars); |
||||
const size_t N = vars[ind].NDim(); |
||||
|
||||
std::unique_ptr<VType[]> cdata; |
||||
|
||||
size_t dsize = 1; |
||||
for(size_t i = 0; i < N; i++) dsize *= count[i]; |
||||
cdata.reset(new VType[dsize]); |
||||
|
||||
int vid; |
||||
int res = nc_inq_varid(ncid, var.Buf(), &vid); |
||||
if(res != NC_NOERR) return Error(pref, MString("nc_inq_varid error: ") + nc_strerror(res)); |
||||
|
||||
if constexpr(std::is_same_v<VType, float>) |
||||
res = nc_get_vara_float(ncid, vid, start, count, cdata.get()); |
||||
else if constexpr(std::is_same_v<VType, double>) |
||||
res = nc_get_vara_double(ncid, vid, start, count, cdata.get()); |
||||
else if constexpr(std::is_same_v<VType, int>) |
||||
res = nc_get_vara_int(ncid, vid, start, count, cdata.get()); |
||||
else if constexpr(std::is_same_v<VType, long>) |
||||
res = nc_get_vara_long(ncid, vid, start, count, cdata.get()); |
||||
else if constexpr(std::is_same_v<VType, long long>) |
||||
res = nc_get_vara_longlong(ncid, vid, start, count, cdata.get()); |
||||
else if constexpr(std::is_same_v<VType, short>) |
||||
res = nc_get_vara_short(ncid, vid, start, count, cdata.get()); |
||||
else if constexpr(std::is_same_v<VType, signed char>) |
||||
res = nc_get_vara_schar(ncid, vid, start, count, cdata.get()); |
||||
else if constexpr(std::is_same_v<VType, unsigned int>) |
||||
res = nc_get_vara_uint(ncid, vid, start, count, cdata.get()); |
||||
else if constexpr(std::is_same_v<VType, unsigned long long>) |
||||
res = nc_get_vara_ulonglong(ncid, vid, start, count, cdata.get()); |
||||
else if constexpr(std::is_same_v<VType, unsigned short>) |
||||
res = nc_get_vara_ushort(ncid, vid, start, count, cdata.get()); |
||||
else if constexpr(std::is_same_v<VType, unsigned char>) |
||||
res = nc_get_vara_ubyte(ncid, vid, start, count, cdata.get()); |
||||
else |
||||
return Error(pref, "Unsupported variable type"); |
||||
|
||||
if(res != NC_NOERR) return Error(pref, MString("nc_get_vara error: ") + nc_strerror(res)); |
||||
|
||||
return ReadedData<VType>(std::move(cdata)); |
||||
} |
||||
|
||||
public: |
||||
~NCSimpleFunctions() { nc_close(ncid); } |
||||
|
||||
Error Open(const MString& filename); |
||||
}; |
||||
|
||||
using NCSimple = NcZarrRead<NCSimpleFunctions>; |
||||
@ -0,0 +1,166 @@
|
||||
#pragma once |
||||
#include "nczarrmulti.h" |
||||
#include <variant> |
||||
|
||||
using NCZarrBase = std::variant<Zarr, NCSimple, ZarrMulti, NCMulti>; |
||||
|
||||
#define DEFFUNC(NAME) \ |
||||
auto NAME() const \
|
||||
{ \
|
||||
return V([](const auto& arg) { return arg.NAME(); }); \
|
||||
} |
||||
#define DEFFUNC1(NAME) \ |
||||
auto NAME(const MString& v) const \
|
||||
{ \
|
||||
return V([&v = std::as_const(v)](const auto& arg) { return arg.NAME(v); }); \
|
||||
} |
||||
#define DEFFUNC2(NAME) \ |
||||
auto NAME(const MString& v1, const MString& v2) const \
|
||||
{ \
|
||||
return V([&v1 = std::as_const(v1), &v2 = std::as_const(v2)](const auto& arg) { return arg.NAME(v1, v2); }); \
|
||||
} |
||||
|
||||
class NCZarr: private NCZarrBase, private DimReqDef |
||||
{ |
||||
template<class Visitor> auto V(Visitor&& visitor) const { return std::visit(std::forward<Visitor>(visitor), *static_cast<const NCZarrBase*>(this)); } |
||||
|
||||
public: |
||||
explicit operator bool() const |
||||
{ |
||||
return V([](const auto& arg) { return static_cast<bool>(arg); }); |
||||
} |
||||
|
||||
DEFFUNC1(NDim) |
||||
|
||||
DEFFUNC(NAtt) |
||||
DEFFUNC1(NAtt) |
||||
|
||||
DEFFUNC(AttNames) |
||||
DEFFUNC1(AttNames) |
||||
|
||||
DEFFUNC(VarNames) |
||||
DEFFUNC1(VarT) |
||||
DEFFUNC1(VarFill) |
||||
|
||||
DEFFUNC1(DimNames) |
||||
|
||||
DEFFUNC2(DimSize) |
||||
|
||||
DEFFUNC2(AttT) |
||||
DEFFUNC2(AttInt) |
||||
DEFFUNC2(AttUInt) |
||||
DEFFUNC2(AttReal) |
||||
DEFFUNC2(AttString) |
||||
DEFFUNC2(AttBool) |
||||
|
||||
DEFFUNC1(AttT) |
||||
DEFFUNC1(AttInt) |
||||
DEFFUNC1(AttUInt) |
||||
DEFFUNC1(AttReal) |
||||
DEFFUNC1(AttString) |
||||
DEFFUNC1(AttBool) |
||||
|
||||
DEFFUNC2(HasDim) |
||||
DEFFUNC1(HasVar) |
||||
DEFFUNC1(HasAtt) |
||||
DEFFUNC2(HasAtt) |
||||
|
||||
template<class Data, class Transform> Error Read(const MString& vname, Data& data, Transform transform, const char* request) const |
||||
{ |
||||
return V([&vname = std::as_const(vname), &data = data, &transform = transform, request = request](const auto& arg) { return arg.Read(vname, data, transform, request); }); |
||||
} |
||||
|
||||
template<class Data, class Transform> Error Read(const MString& vname, Data& data, Transform transform, DimReq&& req1) const |
||||
{ |
||||
return V([&vname = std::as_const(vname), &data = data, &transform = transform, &req1 = req1](const auto& arg) { return arg.Read(vname, data, transform, std::move(req1)); }); |
||||
} |
||||
|
||||
template<class Data, class Transform> Error Read(const MString& vname, Data& data, Transform transform, DimReq&& req1, DimReq&& req2) const |
||||
{ |
||||
return V([&vname = std::as_const(vname), &data = data, &transform = transform, &req1 = req1, &req2 = req2](const auto& arg) |
||||
{ return arg.Read(vname, data, transform, std::move(req1), std::move(req2)); }); |
||||
} |
||||
|
||||
template<class Data, class Transform> Error Read(const MString& vname, Data& data, Transform transform, DimReq&& req1, DimReq&& req2, DimReq&& req3) const |
||||
{ |
||||
return V([&vname = std::as_const(vname), &data = data, &transform = transform, &req1 = req1, &req2 = req2, &req3 = req3](const auto& arg) |
||||
{ return arg.Read(vname, data, transform, std::move(req1), std::move(req2), std::move(req3)); }); |
||||
} |
||||
|
||||
template<class Data, class Transform> Error Read(const MString& vname, Data& data, Transform transform, DimReq&& req1, DimReq&& req2, DimReq&& req3, DimReq&& req4) const |
||||
{ |
||||
return V([&vname = std::as_const(vname), &data = data, &transform = transform, &req1 = req1, &req2 = req2, &req3 = req3, &req4 = req4](const auto& arg) |
||||
{ return arg.Read(vname, data, transform, std::move(req1), std::move(req2), std::move(req3), std::move(req4)); }); |
||||
} |
||||
|
||||
template<class Data, class Transform> Error Read(const MString& vname, Data& data, Transform transform) const |
||||
{ |
||||
return V([&vname = std::as_const(vname), &data = data, &transform = transform](const auto& arg) { return arg.Read(vname, data, transform); }); |
||||
} |
||||
|
||||
template<class Data, class Transform> Error Read(const MString& vname, Data& data, Transform transform, const std::vector<DimReqDef::DimReq>& reqs) const |
||||
{ |
||||
return V([&vname = std::as_const(vname), &data = data, &transform = transform, &reqs = std::as_const(reqs)](const auto& arg) { return arg.Read(vname, data, transform, reqs); }); |
||||
} |
||||
|
||||
template<class Data, class Transform> Error Read(const MString& vname, Data& data, Transform transform, const MString& request) const |
||||
{ |
||||
return V([&vname = std::as_const(vname), &data = data, &transform = transform, &request = std::as_const(request)](const auto& arg) |
||||
{ return arg.Read(vname, data, transform, request); }); |
||||
} |
||||
|
||||
template<class Type> Error Read(const MString& vname, std::vector<Type>& out) const |
||||
{ |
||||
return V([&vname = std::as_const(vname), &out = out](const auto& arg) { return arg.Read(vname, out); }); |
||||
} |
||||
|
||||
Error OpenNC(const MString& filename) |
||||
{ |
||||
auto pv = static_cast<NCZarrBase*>(this); |
||||
|
||||
*pv = NCZarrBase(NCSimple()); |
||||
auto p = std::get_if<NCSimple>(pv); |
||||
if(p == nullptr) return Error("NCZarr::OpenNC", "Impossible error!"); |
||||
return p->Open(filename); |
||||
} |
||||
|
||||
Error OpenZarr(const MString& product, const MString& dataset, bool time = true) |
||||
{ |
||||
auto pv = static_cast<NCZarrBase*>(this); |
||||
|
||||
*pv = NCZarrBase(Zarr()); |
||||
auto p = std::get_if<Zarr>(pv); |
||||
if(p == nullptr) return Error("NCZarr::OpenNC", "Impossible error!"); |
||||
return p->Open(product, dataset, time); |
||||
} |
||||
|
||||
Error OpenMultiNC(const std::vector<MString>& names) |
||||
{ |
||||
if(names.size() == 0) return Error("NCZarr::OpenMultiNC", "no names"); |
||||
if(names.size() == 1) return OpenNC(names[0]); |
||||
|
||||
auto pv = static_cast<NCZarrBase*>(this); |
||||
|
||||
*pv = NCZarrBase(NCMulti()); |
||||
auto p = std::get_if<NCMulti>(pv); |
||||
if(p == nullptr) return Error("NCZarr::OpenMultiNC", "Impossible error!"); |
||||
return p->Open(names); |
||||
} |
||||
|
||||
Error OpenMultiZarr(const MString& product, const std::vector<MString>& dsets, bool time = true) |
||||
{ |
||||
if(dsets.size() == 0) return Error("NCZarr::OpenMultiZarr", "no datasets"); |
||||
if(dsets.size() == 1) return OpenZarr(product, dsets[0], time); |
||||
|
||||
auto pv = static_cast<NCZarrBase*>(this); |
||||
|
||||
*pv = NCZarrBase(ZarrMulti()); |
||||
auto p = std::get_if<ZarrMulti>(pv); |
||||
if(p == nullptr) return Error("NCZarr::OpenMultiZarr", "Impossible error!"); |
||||
return p->Open(product, dsets, time); |
||||
} |
||||
}; |
||||
|
||||
#undef DEFFUNC |
||||
#undef DEFFUNC1 |
||||
#undef DEFFUNC2 |
||||
@ -0,0 +1,826 @@
|
||||
#pragma once |
||||
#include "merrors.h" |
||||
#include <utility> |
||||
#include <variant> |
||||
|
||||
using michlib::Error; |
||||
using michlib::int1; |
||||
using michlib::int2; |
||||
using michlib::int4; |
||||
using michlib::int8; |
||||
using michlib::int_cast; |
||||
using michlib::MString; |
||||
using michlib::RetVal; |
||||
using michlib::uint1; |
||||
using michlib::uint8; |
||||
|
||||
class NcZarrTypes |
||||
{ |
||||
protected: |
||||
using AttVT = std::variant<std::monostate, int8, uint8, double, MString, bool>; |
||||
|
||||
class ArrCounter |
||||
{ |
||||
using VT = std::vector<size_t>; |
||||
const VT count; |
||||
VT ind; |
||||
bool end; |
||||
|
||||
public: |
||||
static size_t Index(const VT& i, const VT& c) |
||||
{ |
||||
size_t out = 0; |
||||
size_t mul = 1; |
||||
for(size_t ii = i.size(); ii != 0; ii--) |
||||
{ |
||||
out += mul * i[ii - 1]; |
||||
mul *= c[ii - 1]; |
||||
} |
||||
return out; |
||||
} |
||||
|
||||
static VT Index(size_t lind, const VT& c) |
||||
{ |
||||
VT out(c.size()); |
||||
|
||||
size_t j = lind; |
||||
|
||||
for(auto i = c.size(); i > 0; i--) |
||||
{ |
||||
out[i - 1] = j % c[i - 1]; |
||||
j = j / c[i - 1]; |
||||
} |
||||
|
||||
return out; |
||||
} |
||||
|
||||
ArrCounter() = delete; |
||||
ArrCounter(const VT& cnt): count(cnt), ind(cnt.size(), 0), end(false) {} |
||||
|
||||
size_t operator[](size_t i) const { return ind[i]; } |
||||
|
||||
ArrCounter& operator++() |
||||
{ |
||||
size_t curind = count.size(); |
||||
while(curind != 0) |
||||
{ |
||||
ind[curind - 1]++; |
||||
if(ind[curind - 1] >= count[curind - 1]) |
||||
{ |
||||
ind[curind - 1] = 0; |
||||
curind--; |
||||
} |
||||
else |
||||
return *this; |
||||
} |
||||
ind = count; |
||||
end = true; |
||||
return *this; |
||||
} |
||||
|
||||
explicit operator bool() const { return !end; } |
||||
|
||||
size_t Index() const { return Index(ind, count); } |
||||
size_t Index(const VT& i) const { return Index(i, count); } |
||||
VT Index(size_t lind) const { return Index(lind, count); } |
||||
|
||||
size_t Count(size_t i) const { return count[i]; } |
||||
|
||||
const VT& VIndex() const { return ind; } |
||||
VT VIndex(const VT& start) const |
||||
{ |
||||
VT out(ind.size()); |
||||
for(size_t i = 0; i < ind.size(); i++) out[i] = ind[i] + start[i]; |
||||
return out; |
||||
} |
||||
|
||||
const auto& Count() const { return count; } |
||||
|
||||
size_t N() const |
||||
{ |
||||
size_t out = 1; |
||||
for(size_t i = 0; i < count.size(); i++) out *= count[i]; |
||||
return out; |
||||
} |
||||
}; |
||||
|
||||
public: |
||||
enum class AttType |
||||
{ |
||||
UNDEF, |
||||
INT, |
||||
UINT, |
||||
REAL, |
||||
STRING, |
||||
BOOL |
||||
}; |
||||
enum class VarType |
||||
{ |
||||
UNDEF, |
||||
FLOAT, |
||||
DOUBLE, |
||||
INT1, |
||||
INT2, |
||||
INT4, |
||||
INT8, |
||||
UINT1 |
||||
}; |
||||
|
||||
protected: |
||||
template<VarType VT, class Dummy = void> struct VarType2Type; |
||||
|
||||
template<class Dummy> struct VarType2Type<VarType::FLOAT, Dummy> |
||||
{ |
||||
using type = float; |
||||
}; |
||||
template<class Dummy> struct VarType2Type<VarType::DOUBLE, Dummy> |
||||
{ |
||||
using type = double; |
||||
}; |
||||
template<class Dummy> struct VarType2Type<VarType::INT1, Dummy> |
||||
{ |
||||
using type = int1; |
||||
}; |
||||
template<class Dummy> struct VarType2Type<VarType::INT2, Dummy> |
||||
{ |
||||
using type = int2; |
||||
}; |
||||
template<class Dummy> struct VarType2Type<VarType::INT4, Dummy> |
||||
{ |
||||
using type = int4; |
||||
}; |
||||
template<class Dummy> struct VarType2Type<VarType::INT8, Dummy> |
||||
{ |
||||
using type = int8; |
||||
}; |
||||
template<class Dummy> struct VarType2Type<VarType::UINT1, Dummy> |
||||
{ |
||||
using type = uint1; |
||||
}; |
||||
|
||||
template<VarType VT> using Type = VarType2Type<VT>::type; |
||||
|
||||
static constexpr size_t SizeOf(VarType vt) |
||||
{ |
||||
switch(vt) |
||||
{ |
||||
case(VarType::UNDEF): return 0; |
||||
case(VarType::FLOAT): return sizeof(Type<VarType::FLOAT>); |
||||
case(VarType::DOUBLE): return sizeof(Type<VarType::DOUBLE>); |
||||
case(VarType::INT1): return sizeof(Type<VarType::INT1>); |
||||
case(VarType::INT2): return sizeof(Type<VarType::INT2>); |
||||
case(VarType::INT4): return sizeof(Type<VarType::INT4>); |
||||
case(VarType::INT8): return sizeof(Type<VarType::INT8>); |
||||
case(VarType::UINT1): return sizeof(Type<VarType::UINT1>); |
||||
} |
||||
return 0; |
||||
} |
||||
|
||||
template<class T> static size_t FindInd(const MString& name, const std::vector<T>& arr) |
||||
{ |
||||
for(size_t i = 0; i < arr.size(); i++) |
||||
if(arr[i].Name() == name) return i; |
||||
return arr.size(); |
||||
} |
||||
|
||||
class Attribute: public AttVT |
||||
{ |
||||
MString name; |
||||
|
||||
public: |
||||
Attribute(const MString& n, AttVT&& v): AttVT(std::move(v)), name(n) {} |
||||
Attribute(const std::string& n, AttVT&& v): AttVT(std::move(v)), name(n.c_str(), n.size()) {} |
||||
|
||||
const MString& Name() const { return name; } |
||||
|
||||
AttType Type() const |
||||
{ |
||||
if(std::holds_alternative<int8>(*this)) |
||||
return AttType::INT; |
||||
else if(std::holds_alternative<uint8>(*this)) |
||||
return AttType::UINT; |
||||
else if(std::holds_alternative<double>(*this)) |
||||
return AttType::REAL; |
||||
else if(std::holds_alternative<MString>(*this)) |
||||
return AttType::STRING; |
||||
else if(std::holds_alternative<bool>(*this)) |
||||
return AttType::BOOL; |
||||
|
||||
return AttType::UNDEF; |
||||
} |
||||
|
||||
int8 I() const |
||||
{ |
||||
if(std::holds_alternative<int8>(*this)) |
||||
return std::get<int8>(*this); |
||||
else if(std::holds_alternative<uint8>(*this)) |
||||
return int_cast<int8>(std::get<uint8>(*this)); |
||||
else if(std::holds_alternative<double>(*this)) |
||||
return static_cast<int8>(std::get<double>(*this)); |
||||
else if(std::holds_alternative<MString>(*this)) |
||||
return std::get<MString>(*this).ToInteger<int8>(); |
||||
else if(std::holds_alternative<bool>(*this)) |
||||
return std::get<bool>(*this) ? 1 : 0; |
||||
return 0; |
||||
} |
||||
|
||||
uint8 U() const |
||||
{ |
||||
if(std::holds_alternative<int8>(*this)) |
||||
return int_cast<uint8>(std::get<int8>(*this)); |
||||
else if(std::holds_alternative<uint8>(*this)) |
||||
return std::get<uint8>(*this); |
||||
else if(std::holds_alternative<double>(*this)) |
||||
return static_cast<uint8>(std::get<double>(*this)); |
||||
else if(std::holds_alternative<MString>(*this)) |
||||
return std::get<MString>(*this).ToInteger<uint8>(); |
||||
else if(std::holds_alternative<bool>(*this)) |
||||
return std::get<bool>(*this) ? 1 : 0; |
||||
return 0; |
||||
} |
||||
|
||||
double D() const |
||||
{ |
||||
if(std::holds_alternative<int8>(*this)) |
||||
return std::get<int8>(*this); |
||||
else if(std::holds_alternative<uint8>(*this)) |
||||
return std::get<uint8>(*this); |
||||
else if(std::holds_alternative<double>(*this)) |
||||
return std::get<double>(*this); |
||||
else if(std::holds_alternative<MString>(*this)) |
||||
return michlib_internal::RealType<sizeof(double)>::String2Real(std::get<MString>(*this).Buf()); |
||||
else if(std::holds_alternative<bool>(*this)) |
||||
return std::get<bool>(*this) ? 1 : 0; |
||||
return 0; |
||||
} |
||||
|
||||
MString S() const |
||||
{ |
||||
if(std::holds_alternative<int8>(*this)) |
||||
return MString().FromInt(std::get<int8>(*this)); |
||||
else if(std::holds_alternative<uint8>(*this)) |
||||
return MString().FromUInt(std::get<uint8>(*this)); |
||||
else if(std::holds_alternative<double>(*this)) |
||||
return MString().FromReal(std::get<double>(*this)); |
||||
else if(std::holds_alternative<MString>(*this)) |
||||
return std::get<MString>(*this); |
||||
else if(std::holds_alternative<bool>(*this)) |
||||
return MString().FromBool(std::get<bool>(*this)); |
||||
return ""; |
||||
} |
||||
|
||||
bool B() const |
||||
{ |
||||
if(std::holds_alternative<int8>(*this)) |
||||
return std::get<int8>(*this) != 0; |
||||
else if(std::holds_alternative<uint8>(*this)) |
||||
return std::get<uint8>(*this) != 0; |
||||
else if(std::holds_alternative<double>(*this)) |
||||
return std::get<double>(*this) != 0.0; |
||||
else if(std::holds_alternative<MString>(*this)) |
||||
return std::get<MString>(*this).ToBool(); |
||||
else if(std::holds_alternative<bool>(*this)) |
||||
return std::get<bool>(*this); |
||||
return false; |
||||
} |
||||
}; |
||||
|
||||
class Dimension |
||||
{ |
||||
MString name; |
||||
size_t size; |
||||
|
||||
public: |
||||
Dimension(const MString& str, size_t num): name(str), size(num) {} |
||||
const MString& Name() const { return name; } |
||||
|
||||
size_t Size() const { return size; } |
||||
}; |
||||
|
||||
class Variable |
||||
{ |
||||
public: |
||||
using FillType = std::variant<std::monostate, int8, uint8, double>; |
||||
|
||||
private: |
||||
MString name; |
||||
VarType type = VarType::UNDEF; |
||||
std::vector<size_t> dims; |
||||
std::vector<Attribute> atts; |
||||
FillType fill; |
||||
|
||||
public: |
||||
Variable(const MString& name_, VarType type_, std::vector<size_t>&& dims_, std::vector<Attribute>&& atts_, FillType fill_ = 0): |
||||
name(name_), type(type_), dims(std::move(dims_)), atts(std::move(atts_)), fill(fill_) |
||||
{ |
||||
} |
||||
|
||||
explicit operator bool() const { return type != VarType::UNDEF; } |
||||
|
||||
const auto& Dims() const { return dims; } |
||||
|
||||
size_t NDim() const { return dims.size(); } |
||||
|
||||
size_t NAtt() const { return atts.size(); } |
||||
|
||||
auto AttNames() const |
||||
{ |
||||
std::vector<MString> out; |
||||
std::transform(atts.cbegin(), atts.cend(), std::back_inserter(out), [](const Attribute& a) { return a.Name(); }); |
||||
return out; |
||||
} |
||||
|
||||
AttType AttT(const MString& name) const |
||||
{ |
||||
size_t ind = FindInd(name, atts); |
||||
return ind < atts.size() ? atts[ind].Type() : AttType::UNDEF; |
||||
} |
||||
|
||||
int8 AttInt(const MString& name) const |
||||
{ |
||||
size_t ind = FindInd(name, atts); |
||||
return ind < atts.size() ? atts[ind].I() : 0; |
||||
} |
||||
|
||||
uint8 AttUInt(const MString& name) const |
||||
{ |
||||
size_t ind = FindInd(name, atts); |
||||
return ind < atts.size() ? atts[ind].U() : 0; |
||||
} |
||||
|
||||
double AttReal(const MString& name) const |
||||
{ |
||||
size_t ind = FindInd(name, atts); |
||||
return ind < atts.size() ? atts[ind].D() : 0.0; |
||||
} |
||||
|
||||
MString AttString(const MString& name) const |
||||
{ |
||||
size_t ind = FindInd(name, atts); |
||||
return ind < atts.size() ? atts[ind].S() : MString(); |
||||
} |
||||
|
||||
bool AttBool(const MString& name) const |
||||
{ |
||||
size_t ind = FindInd(name, atts); |
||||
return ind < atts.size() ? atts[ind].B() : false; |
||||
} |
||||
|
||||
const MString& Name() const { return name; } |
||||
|
||||
auto Type() const { return type; } |
||||
|
||||
const auto& Fill() const { return fill; } |
||||
}; |
||||
|
||||
protected: |
||||
std::vector<Attribute> gats; |
||||
std::vector<Dimension> dims; |
||||
std::vector<Variable> vars; |
||||
|
||||
public: |
||||
explicit operator bool() const { return !vars.empty(); } |
||||
|
||||
size_t NDim(const MString& var) const |
||||
{ |
||||
size_t ind = FindInd(var, vars); |
||||
return ind < vars.size() ? vars[ind].NDim() : 0; |
||||
} |
||||
|
||||
size_t NAtt() const { return gats.size(); } |
||||
|
||||
auto AttNames() const |
||||
{ |
||||
std::vector<MString> out; |
||||
std::transform(gats.cbegin(), gats.cend(), std::back_inserter(out), [](const Attribute& a) { return a.Name(); }); |
||||
return out; |
||||
} |
||||
|
||||
size_t NAtt(const MString& var) const |
||||
{ |
||||
if(!var.Exist()) return NAtt(); |
||||
size_t ind = FindInd(var, vars); |
||||
return ind < vars.size() ? vars[ind].NAtt() : 0; |
||||
} |
||||
|
||||
auto AttNames(const MString& var) const |
||||
{ |
||||
if(!var.Exist()) return AttNames(); |
||||
size_t ind = FindInd(var, vars); |
||||
return ind < vars.size() ? vars[ind].AttNames() : decltype(AttNames())(); |
||||
} |
||||
|
||||
auto VarNames() const |
||||
{ |
||||
std::vector<MString> out; |
||||
std::transform(vars.cbegin(), vars.cend(), std::back_inserter(out), [](const Variable& v) { return v.Name(); }); |
||||
return out; |
||||
} |
||||
|
||||
VarType VarT(const MString& var) const |
||||
{ |
||||
size_t ind = FindInd(var, vars); |
||||
return ind < vars.size() ? vars[ind].Type() : VarType::UNDEF; |
||||
} |
||||
|
||||
auto VarFill(const MString& var) const |
||||
{ |
||||
size_t ind = FindInd(var, vars); |
||||
return ind < vars.size() ? vars[ind].Fill() : Variable::FillType(); |
||||
} |
||||
|
||||
auto DimNames(const MString& var) const |
||||
{ |
||||
size_t ind = FindInd(var, vars); |
||||
|
||||
std::vector<MString> out; |
||||
if(ind >= vars.size()) return out; |
||||
|
||||
auto vdims = vars[ind].Dims(); |
||||
std::transform(vdims.cbegin(), vdims.cend(), std::back_inserter(out), [&dims = std::as_const(dims)](const size_t& i) { return dims[i].Name(); }); |
||||
return out; |
||||
} |
||||
|
||||
size_t DimSize(const MString& var, const MString& dim) const |
||||
{ |
||||
if(!HasDim(var, dim)) return 0; |
||||
|
||||
size_t ind = FindInd(dim, dims); |
||||
return ind < dims.size() ? dims[ind].Size() : 0; |
||||
} |
||||
|
||||
AttType AttT(const MString& var, const MString& name) const |
||||
{ |
||||
if(!var.Exist()) |
||||
{ |
||||
size_t ind = FindInd(name, gats); |
||||
return ind < gats.size() ? gats[ind].Type() : AttType::UNDEF; |
||||
} |
||||
|
||||
size_t ind = FindInd(var, vars); |
||||
return ind < vars.size() ? vars[ind].AttT(name) : AttType::UNDEF; |
||||
} |
||||
|
||||
int8 AttInt(const MString& var, const MString& name) const |
||||
{ |
||||
if(!var.Exist()) |
||||
{ |
||||
size_t ind = FindInd(name, gats); |
||||
return ind < gats.size() ? gats[ind].I() : 0; |
||||
} |
||||
|
||||
size_t ind = FindInd(var, vars); |
||||
return ind < vars.size() ? vars[ind].AttInt(name) : 0; |
||||
} |
||||
|
||||
uint8 AttUInt(const MString& var, const MString& name) const |
||||
{ |
||||
if(!var.Exist()) |
||||
{ |
||||
size_t ind = FindInd(name, gats); |
||||
return ind < gats.size() ? gats[ind].U() : 0; |
||||
} |
||||
|
||||
size_t ind = FindInd(var, vars); |
||||
return ind < vars.size() ? vars[ind].AttUInt(name) : 0; |
||||
} |
||||
|
||||
double AttReal(const MString& var, const MString& name) const |
||||
{ |
||||
if(!var.Exist()) |
||||
{ |
||||
size_t ind = FindInd(name, gats); |
||||
return ind < gats.size() ? gats[ind].D() : 0.0; |
||||
} |
||||
|
||||
size_t ind = FindInd(var, vars); |
||||
return ind < vars.size() ? vars[ind].AttReal(name) : 0.0; |
||||
} |
||||
|
||||
MString AttString(const MString& var, const MString& name) const |
||||
{ |
||||
if(!var.Exist()) |
||||
{ |
||||
size_t ind = FindInd(name, gats); |
||||
return ind < gats.size() ? gats[ind].S() : MString(); |
||||
} |
||||
|
||||
size_t ind = FindInd(var, vars); |
||||
return ind < vars.size() ? vars[ind].AttString(name) : MString(); |
||||
} |
||||
|
||||
bool AttBool(const MString& var, const MString& name) const |
||||
{ |
||||
if(!var.Exist()) |
||||
{ |
||||
size_t ind = FindInd(name, gats); |
||||
return ind < gats.size() ? gats[ind].B() : false; |
||||
} |
||||
|
||||
size_t ind = FindInd(var, vars); |
||||
return ind < vars.size() ? vars[ind].AttBool(name) : false; |
||||
} |
||||
|
||||
auto AttT(const MString& name) const { return AttT("", name); } |
||||
auto AttInt(const MString& name) const { return AttInt("", name); } |
||||
auto AttUInt(const MString& name) const { return AttUInt("", name); } |
||||
auto AttReal(const MString& name) const { return AttReal("", name); } |
||||
auto AttString(const MString& name) const { return AttString("", name); } |
||||
auto AttBool(const MString& name) const { return AttBool("", name); } |
||||
|
||||
bool HasDim(const MString& var, const MString& name) const |
||||
{ |
||||
size_t vind = FindInd(var, vars); |
||||
if(vind >= vars.size()) return false; |
||||
|
||||
for(const auto dind: vars[vind].Dims()) |
||||
if(dims[dind].Name() == name) return true; |
||||
return false; |
||||
} |
||||
|
||||
bool HasVar(const MString& name) const { return FindInd(name, vars) < vars.size(); } |
||||
bool HasAtt(const MString& vname, const MString& aname) const { return AttT(vname, aname) != AttType::UNDEF; } |
||||
bool HasAtt(const MString& aname) const { return AttT(aname) != AttType::UNDEF; } |
||||
}; |
||||
|
||||
class DimReqDef |
||||
{ |
||||
protected: |
||||
struct DimReq |
||||
{ |
||||
static const auto fill = std::numeric_limits<size_t>::max(); |
||||
MString name; |
||||
size_t beg, count; |
||||
|
||||
DimReq(): name(MString()), beg(fill), count(fill) {} |
||||
DimReq(const char* n): name(n), beg(fill), count(fill) {} |
||||
DimReq(const MString& n): name(n), beg(fill), count(fill) {} |
||||
DimReq(MString&& n): name(std::move(n)), beg(fill), count(fill) {} |
||||
DimReq(const char* n, size_t s): name(n), beg(s), count(fill) {} |
||||
DimReq(const MString& n, size_t s): name(n), beg(s), count(fill) {} |
||||
DimReq(MString&& n, size_t s): name(std::move(n)), beg(s), count(fill) {} |
||||
DimReq(const char* n, size_t s, size_t c): name(n), beg(s), count(c) {} |
||||
DimReq(const MString& n, size_t s, size_t c): name(n), beg(s), count(c) {} |
||||
DimReq(MString&& n, size_t s, size_t c): name(std::move(n)), beg(s), count(c) {} |
||||
|
||||
const MString& Name() const { return name; } |
||||
}; |
||||
}; |
||||
|
||||
template<class C> class NcZarrRead: public C, public DimReqDef |
||||
{ |
||||
template<class Data> static constexpr size_t Dimensionity() |
||||
{ |
||||
if constexpr(requires(Data& d) { d(0, 0, 0, 0); }) return 4; |
||||
if constexpr(requires(Data& d) { d(0, 0, 0); }) return 3; |
||||
if constexpr(requires(Data& d) { d(0, 0); }) return 2; |
||||
if constexpr(requires(Data& d) { d(0); }) return 1; |
||||
return 0; |
||||
} |
||||
|
||||
template<class Data, size_t D, class Dummy = void> struct DataTypeExtractorS; |
||||
|
||||
template<class Data, class Dummy> struct DataTypeExtractorS<Data, 1, Dummy> |
||||
{ |
||||
using type = std::decay_t<decltype(std::declval<Data>()(0))>; |
||||
}; |
||||
template<class Data, class Dummy> struct DataTypeExtractorS<Data, 2, Dummy> |
||||
{ |
||||
using type = std::decay_t<decltype(std::declval<Data>()(0, 0))>; |
||||
}; |
||||
template<class Data, class Dummy> struct DataTypeExtractorS<Data, 3, Dummy> |
||||
{ |
||||
using type = std::decay_t<decltype(std::declval<Data>()(0, 0, 0))>; |
||||
}; |
||||
template<class Data, class Dummy> struct DataTypeExtractorS<Data, 4, Dummy> |
||||
{ |
||||
using type = std::decay_t<decltype(std::declval<Data>()(0, 0, 0, 0))>; |
||||
}; |
||||
|
||||
template<class Data> using DataTypeExtractor = DataTypeExtractorS<Data, Dimensionity<Data>()>::type; |
||||
|
||||
template<class VType, class Data, class Transform> |
||||
Error Read(const MString& vname, const std::vector<size_t>& transindex, Data& data, Transform transform, std::vector<DimReq> reqs) const |
||||
{ |
||||
const size_t indim = reqs.size(); |
||||
constexpr size_t outdim = Dimensionity<Data>(); |
||||
|
||||
std::vector<size_t> start; |
||||
std::vector<size_t> count; |
||||
|
||||
start.resize(indim); |
||||
count.resize(indim); |
||||
for(size_t i = 0; i < indim; i++) |
||||
{ |
||||
start[i] = reqs[i].beg; |
||||
count[i] = reqs[i].count; |
||||
} |
||||
|
||||
using DataType = DataTypeExtractor<Data>; |
||||
DataType fillout; |
||||
bool havefill = C::VarFill(vname).index() > 0; |
||||
VType fillin = std::visit( |
||||
[](auto v) |
||||
{ |
||||
if constexpr(std::is_convertible_v<decltype(v), VType>) |
||||
return static_cast<VType>(v); |
||||
else |
||||
return std::numeric_limits<VType>::max(); |
||||
}, |
||||
C::VarFill(vname)); |
||||
|
||||
if constexpr(requires(Data& d) { // Data have own fillvalue
|
||||
{ d.Fillval() } -> std::convertible_to<DataType>; |
||||
}) |
||||
fillout = data.Fillval(); |
||||
else // Data does'nt have own fillvalue, using variable fillvalue
|
||||
fillout = static_cast<DataType>(fillin); |
||||
|
||||
auto ret = C::template Read<VType>(vname, start.data(), count.data()); |
||||
if(!ret) return ret; |
||||
const auto& rawdata = ret.Value(); |
||||
|
||||
std::vector<size_t> mul(indim, 1); |
||||
for(size_t i = indim - 1; i > 0; i--) mul[i - 1] = mul[i] * count[i]; |
||||
|
||||
size_t inind = 0; |
||||
for(typename C::ArrCounter i(count); i; ++i) |
||||
{ |
||||
// TODO: Remove this testing block
|
||||
size_t cind = 0; |
||||
for(size_t j = 0; j < indim; j++) cind += i[j] * mul[j]; |
||||
if(cind != inind) return {"NcZarrRead::Read", "Internal error"}; |
||||
if(i.Index() != inind) return {"NcZarrRead::Read", "Internal error"}; |
||||
if(inind != i.Index(i.Index(inind, count), count)) return {"NcZarrRead::Read", "Internal error"}; |
||||
|
||||
DataType out; |
||||
const VType& in = rawdata(inind); |
||||
if(havefill && in == fillin) |
||||
out = fillout; |
||||
else |
||||
out = transform(in); |
||||
|
||||
if constexpr(outdim == 1) |
||||
data(i[transindex[0]]) = out; |
||||
else if constexpr(outdim == 2) |
||||
data(i[transindex[0]], i[transindex[1]]) = out; |
||||
else if constexpr(outdim == 3) |
||||
data(i[transindex[0]], i[transindex[1]], i[transindex[2]]) = out; |
||||
else if constexpr(outdim == 4) |
||||
data(i[transindex[0]], i[transindex[1]], i[transindex[2]], i[transindex[3]]) = out; |
||||
|
||||
inind++; |
||||
} |
||||
|
||||
return Error(); |
||||
} |
||||
|
||||
public: |
||||
// Request is string
|
||||
template<class Data, class Transform> Error Read(const MString& vname, Data& data, Transform transform, const char* request) const |
||||
{ |
||||
return Read(vname, data, transform, MString(request)); |
||||
} |
||||
|
||||
// Request by one dimension
|
||||
template<class Data, class Transform> Error Read(const MString& vname, Data& data, Transform transform, DimReq&& req1) const |
||||
{ |
||||
return Read(vname, data, transform, std::vector<DimReq>{std::move(req1)}); |
||||
} |
||||
// Request by two dimension
|
||||
template<class Data, class Transform> Error Read(const MString& vname, Data& data, Transform transform, DimReq&& req1, DimReq&& req2) const |
||||
{ |
||||
return Read(vname, data, transform, std::vector<DimReq>{std::move(req1), std::move(req2)}); |
||||
} |
||||
// Request by three dimension
|
||||
template<class Data, class Transform> Error Read(const MString& vname, Data& data, Transform transform, DimReq&& req1, DimReq&& req2, DimReq&& req3) const |
||||
{ |
||||
return Read(vname, data, transform, std::vector<DimReq>{std::move(req1), std::move(req2), std::move(req3)}); |
||||
} |
||||
// Request by four dimension
|
||||
template<class Data, class Transform> Error Read(const MString& vname, Data& data, Transform transform, DimReq&& req1, DimReq&& req2, DimReq&& req3, DimReq&& req4) const |
||||
{ |
||||
return Read(vname, data, transform, std::vector<DimReq>{std::move(req1), std::move(req2), std::move(req3), std::move(req4)}); |
||||
} |
||||
|
||||
// Request full variable
|
||||
template<class Data, class Transform> Error Read(const MString& vname, Data& data, Transform transform) const |
||||
{ |
||||
static const MString pref = "NcZarrRead::Read"; |
||||
if(!C::HasVar(vname)) return {pref, "Variable " + vname + " not found"}; |
||||
|
||||
std::vector<struct DimReq> pdims; |
||||
|
||||
const auto vdims = C::DimNames(vname); |
||||
std::transform(vdims.cbegin(), vdims.cend(), std::back_inserter(pdims), |
||||
[this, &vname = std::as_const(vname)](const MString& n) -> struct DimReq { return {n, 0, C::DimSize(vname, n)}; }); |
||||
|
||||
return Read(vname, data, transform, pdims); |
||||
} |
||||
|
||||
// Base function for all Read's
|
||||
template<class Data, class Transform> Error Read(const MString& vname, Data& data, Transform transform, const std::vector<DimReq>& reqs) const |
||||
{ |
||||
static const MString pref = "NcZarrRead::Read"; |
||||
|
||||
if(!C::HasVar(vname)) return {pref, "Variable " + vname + " not found"}; |
||||
|
||||
std::vector<struct DimReq> pdims; |
||||
{ |
||||
const auto vdims = C::DimNames(vname); |
||||
std::transform(vdims.cbegin(), vdims.cend(), std::back_inserter(pdims), [](const MString& n) -> struct DimReq { return {n, 0, 1}; }); |
||||
} |
||||
|
||||
std::vector<size_t> transindex; |
||||
|
||||
// Parse request
|
||||
if(reqs.size() == 0) return {pref, "Empty request"}; |
||||
for(const auto& req: reqs) |
||||
{ |
||||
size_t ind = C::FindInd(req.name, pdims); |
||||
if(ind >= pdims.size()) return {pref, "Variable " + vname + " has no dimension " + req.name}; |
||||
|
||||
for(size_t i = 0; i < transindex.size(); i++) |
||||
if(transindex[i] == ind) return {pref, "Parameters for dimension " + req.name + " already defined"}; |
||||
transindex.push_back(ind); |
||||
|
||||
size_t dlen = C::DimSize(vname, pdims[ind].name); |
||||
if(req.beg == req.fill && req.count == req.fill) // Only name, so, we request full length
|
||||
{ |
||||
pdims[ind].beg = 0; |
||||
pdims[ind].count = dlen; |
||||
} |
||||
else if(req.count == req.fill) // Name and first index
|
||||
{ |
||||
pdims[ind].beg = req.beg; |
||||
pdims[ind].count = 1; |
||||
} |
||||
else // Name, first index, count
|
||||
{ |
||||
pdims[ind].beg = req.beg; |
||||
pdims[ind].count = req.count; |
||||
} |
||||
// Sanity checks
|
||||
if(pdims[ind].count <= 0) return {pref, "Error parsing request: count must be greter then zero"}; |
||||
if(pdims[ind].beg >= dlen) return {pref, MString("Error parsing request: start index ") + pdims[ind].beg + " must be lesser then " + pdims[ind].name + " size " + dlen}; |
||||
if(pdims[ind].beg + pdims[ind].count > dlen) |
||||
return {pref, MString("Error parsing request: start index ") + pdims[ind].beg + " with count " + pdims[ind].count + " exceeds " + pdims[ind].name + " size " + dlen}; |
||||
|
||||
// Ignore hyperplanes in requests for calculation of data dimensionality
|
||||
if(pdims[transindex.back()].count == 1) transindex.pop_back(); |
||||
} |
||||
|
||||
if(transindex.size() != Dimensionity<Data>()) |
||||
return {pref, MString("Output data dimensions (") + Dimensionity<Data>() + ") not corresponding request dimensions (" + transindex.size() + ")"}; |
||||
switch(C::VarT(vname)) |
||||
{ |
||||
case(C::VarType::UNDEF): return {pref, "No variable with name " + vname + " (impossible)"}; |
||||
case(C::VarType::FLOAT): return Read<typename C::template Type<C::VarType::FLOAT>>(vname, transindex, data, transform, pdims); |
||||
case(C::VarType::DOUBLE): return Read<typename C::template Type<C::VarType::DOUBLE>>(vname, transindex, data, transform, pdims); |
||||
case(C::VarType::INT1): return Read<typename C::template Type<C::VarType::INT1>>(vname, transindex, data, transform, pdims); |
||||
case(C::VarType::INT2): return Read<typename C::template Type<C::VarType::INT2>>(vname, transindex, data, transform, pdims); |
||||
case(C::VarType::INT4): return Read<typename C::template Type<C::VarType::INT4>>(vname, transindex, data, transform, pdims); |
||||
case(C::VarType::INT8): return Read<typename C::template Type<C::VarType::INT8>>(vname, transindex, data, transform, pdims); |
||||
case(C::VarType::UINT1): return Read<typename C::template Type<C::VarType::INT1>>(vname, transindex, data, transform, pdims); |
||||
} |
||||
|
||||
return {pref, "Internal error (impossible)"}; |
||||
} |
||||
|
||||
// Request by string argument
|
||||
template<class Data, class Transform> Error Read(const MString& vname, Data& data, Transform transform, const MString& request) const |
||||
{ |
||||
static const MString pref = "NcZarrRead::Read"; |
||||
|
||||
std::vector<struct DimReq> pdims; |
||||
|
||||
// Parse request
|
||||
const auto dimdesc = request.Split(";, \t"); |
||||
if(dimdesc.size() == 0) return {pref, "Empty request"}; |
||||
for(const auto& dd: dimdesc) |
||||
{ |
||||
const auto dimpar = dd.Split(":", true); |
||||
|
||||
if(dimpar.size() == 1) // Only name, so, we request full length
|
||||
pdims.emplace_back(dimpar[0]); |
||||
else if(dimpar.size() == 2) // Name and first index
|
||||
pdims.emplace_back(dimpar[0], dimpar[1].ToInteger<size_t>()); |
||||
else if(dimpar.size() == 3) // Name, first index, count
|
||||
pdims.emplace_back(dimpar[0], dimpar[1].ToInteger<size_t>(), dimpar[2].ToInteger<size_t>()); |
||||
else |
||||
return {pref, "Can't parse expression " + dd}; |
||||
} |
||||
|
||||
return Read(vname, data, transform, pdims); |
||||
} |
||||
|
||||
// Request full one-dimensional variable
|
||||
template<class Type> Error Read(const MString& vname, std::vector<Type>& out) const |
||||
{ |
||||
const auto& dnames = C::DimNames(vname); |
||||
if(dnames.size() > 0) out.resize(C::DimSize(vname, dnames[0])); |
||||
auto data = [&vec = out](size_t i) -> Type& { return vec[i]; }; |
||||
|
||||
return Read(vname, data, std::identity()); |
||||
} |
||||
}; |
||||
@ -0,0 +1,206 @@
|
||||
#pragma once |
||||
#include "ncsimple.h" |
||||
#include "zarr.h" |
||||
#include <set> |
||||
|
||||
#define DEFFUNC(func, def) \ |
||||
auto func(const MString& var) const \
|
||||
{ \
|
||||
const auto ind = FindVar(var); \
|
||||
return ind < sources.size() ? sources[ind].func(var) : (def); \
|
||||
} |
||||
#define DEFFUNCA(func, def) \ |
||||
auto func(const MString& att) const \
|
||||
{ \
|
||||
const auto ind = FindAtt(att); \
|
||||
return ind < sources.size() ? sources[ind].func(att) : (def); \
|
||||
} |
||||
#define DEFFUNC2(func, def) \ |
||||
auto func(const MString& var, const MString& name) const \
|
||||
{ \
|
||||
const auto ind = FindVar(var); \
|
||||
return ind < sources.size() ? sources[ind].func(var, name) : (def); \
|
||||
} |
||||
#define DEFFUNCNAMES(func) \ |
||||
auto func() const \
|
||||
{ \
|
||||
std::set<MString> names; \
|
||||
for(const auto& s: sources) \
|
||||
{ \
|
||||
auto snames = s.func(); \
|
||||
for(const auto& name: snames) names.insert(name); \
|
||||
} \
|
||||
std::vector<MString> out; \
|
||||
for(const auto& name: names) out.push_back(name); \
|
||||
return out; \
|
||||
} |
||||
|
||||
template<class T> class NCZarrMultiCommon: private DimReqDef |
||||
{ |
||||
protected: |
||||
std::vector<T> sources; |
||||
|
||||
private: |
||||
size_t FindVar(const MString& var) const |
||||
{ |
||||
for(size_t i = 0; i < sources.size(); i++) |
||||
if(sources[i].HasVar(var)) return i; |
||||
return sources.size(); |
||||
} |
||||
|
||||
size_t FindAtt(const MString& att) const |
||||
{ |
||||
for(size_t i = 0; i < sources.size(); i++) |
||||
if(sources[i].HasAtt(att)) return i; |
||||
return sources.size(); |
||||
} |
||||
|
||||
public: |
||||
explicit operator bool() const |
||||
{ |
||||
if(sources.empty()) return false; |
||||
for(const auto& s: sources) |
||||
if(!s) return false; |
||||
return true; |
||||
} |
||||
|
||||
size_t NDim() const { return DimNames().size(); } |
||||
DEFFUNC(NDim, 0) |
||||
|
||||
size_t NAtt() const { return AttNames().size(); } |
||||
DEFFUNC(NAtt, 0) |
||||
|
||||
DEFFUNCNAMES(AttNames) |
||||
DEFFUNC(AttNames, decltype(AttNames())()) |
||||
|
||||
DEFFUNCNAMES(VarNames) |
||||
DEFFUNC(VarT, T::VarType::UNDEF) |
||||
DEFFUNC(VarFill, decltype(T().VarFill(MString()))()) |
||||
|
||||
DEFFUNC(DimNames, decltype(T().DimNames(MString()))()) |
||||
|
||||
DEFFUNC2(DimSize, 0) |
||||
|
||||
DEFFUNC2(AttT, T::AttType::UNDEF) |
||||
DEFFUNC2(AttInt, 0) |
||||
DEFFUNC2(AttUInt, 0) |
||||
DEFFUNC2(AttReal, 0.0) |
||||
DEFFUNC2(AttString, MString()) |
||||
DEFFUNC2(AttBool, false) |
||||
|
||||
DEFFUNCA(AttT, T::AttType::UNDEF) |
||||
DEFFUNCA(AttInt, 0) |
||||
DEFFUNCA(AttUInt, 0) |
||||
DEFFUNCA(AttReal, 0.0) |
||||
DEFFUNCA(AttString, MString()) |
||||
DEFFUNCA(AttBool, false) |
||||
|
||||
DEFFUNC2(HasDim, false) |
||||
DEFFUNC(HasVar, false) |
||||
DEFFUNCA(HasAtt, false) |
||||
DEFFUNC2(HasAtt, false) |
||||
|
||||
// Request is string
|
||||
template<class Data, class Transform> Error Read(const MString& vname, Data& data, Transform transform, const char* request) const |
||||
{ |
||||
return Read(vname, data, transform, MString(request)); |
||||
} |
||||
|
||||
// Request by one dimension
|
||||
template<class Data, class Transform> Error Read(const MString& vname, Data& data, Transform transform, DimReq&& req1) const |
||||
{ |
||||
return Read(vname, data, transform, std::vector<DimReq>{std::move(req1)}); |
||||
} |
||||
// Request by two dimension
|
||||
template<class Data, class Transform> Error Read(const MString& vname, Data& data, Transform transform, DimReq&& req1, DimReq&& req2) const |
||||
{ |
||||
return Read(vname, data, transform, std::vector<DimReq>{std::move(req1), std::move(req2)}); |
||||
} |
||||
// Request by three dimension
|
||||
template<class Data, class Transform> Error Read(const MString& vname, Data& data, Transform transform, DimReq&& req1, DimReq&& req2, DimReq&& req3) const |
||||
{ |
||||
return Read(vname, data, transform, std::vector<DimReq>{std::move(req1), std::move(req2), std::move(req3)}); |
||||
} |
||||
// Request by four dimension
|
||||
template<class Data, class Transform> Error Read(const MString& vname, Data& data, Transform transform, DimReq&& req1, DimReq&& req2, DimReq&& req3, DimReq&& req4) const |
||||
{ |
||||
return Read(vname, data, transform, std::vector<DimReq>{std::move(req1), std::move(req2), std::move(req3), std::move(req4)}); |
||||
} |
||||
|
||||
// Request full variable
|
||||
template<class Data, class Transform> Error Read(const MString& vname, Data& data, Transform transform) const |
||||
{ |
||||
auto ind = FindVar(vname); |
||||
if(ind < sources.size()) return sources[ind].Read(vname, data, transform); |
||||
return Error("NCZarrMultiCommon::Read", "Variable " + vname + " not found"); |
||||
} |
||||
|
||||
// Base function for all Read's
|
||||
template<class Data, class Transform> Error Read(const MString& vname, Data& data, Transform transform, const std::vector<DimReqDef::DimReq>& reqs) const |
||||
{ |
||||
auto ind = FindVar(vname); |
||||
if(ind < sources.size()) return sources[ind].Read(vname, data, transform, reqs); |
||||
return Error("NCZarrMultiCommon::Read", "Variable " + vname + " not found"); |
||||
} |
||||
|
||||
// Request by string argument
|
||||
template<class Data, class Transform> Error Read(const MString& vname, Data& data, Transform transform, const MString& request) const |
||||
{ |
||||
auto ind = FindVar(vname); |
||||
if(ind < sources.size()) return sources[ind].Read(vname, data, transform, request); |
||||
return Error("NCZarrMultiCommon::Read", "Variable " + vname + " not found"); |
||||
} |
||||
|
||||
// Request full one-dimensional variable
|
||||
template<class Type> Error Read(const MString& vname, std::vector<Type>& out) const |
||||
{ |
||||
auto ind = FindVar(vname); |
||||
if(ind < sources.size()) return sources[ind].Read(vname, out); |
||||
return Error("NCZarrMultiCommon::Read", "Variable " + vname + " not found"); |
||||
} |
||||
}; |
||||
|
||||
#undef DEFFUNC |
||||
#undef DEFFUNCA |
||||
#undef DEFFUNC2 |
||||
#undef DEFFUNCNAME |
||||
|
||||
class NCMulti: public NCZarrMultiCommon<NCSimple> |
||||
{ |
||||
public: |
||||
Error Open(const std::vector<MString>& names) |
||||
{ |
||||
sources.clear(); |
||||
decltype(sources) newsources; |
||||
static const MString pref = "NcMulti::Open"; |
||||
if(names.size() == 0) return Error(pref, "empty file list"); |
||||
for(const auto& name: names) |
||||
{ |
||||
newsources.emplace_back(); |
||||
auto ret = newsources.back().Open(name); |
||||
if(!ret) return ret.Add(pref, "Can't open file " + name); |
||||
} |
||||
sources = std::move(newsources); |
||||
return Error(); |
||||
} |
||||
}; |
||||
|
||||
class ZarrMulti: public NCZarrMultiCommon<Zarr> |
||||
{ |
||||
public: |
||||
Error Open(const MString& product, const std::vector<MString>& dsets, bool time = true) |
||||
{ |
||||
sources.clear(); |
||||
decltype(sources) newsources; |
||||
static const MString pref = "ZarrMulti::Open"; |
||||
if(dsets.size() == 0) return Error(pref, "empty datasets list"); |
||||
for(const auto& dset: dsets) |
||||
{ |
||||
newsources.emplace_back(); |
||||
auto ret = newsources.back().Open(product, dset, time); |
||||
if(!ret) return ret.Add(pref, "Can't open dataset " + dset); |
||||
} |
||||
sources = std::move(newsources); |
||||
return Error(); |
||||
} |
||||
}; |
||||
@ -0,0 +1,178 @@
|
||||
#pragma once |
||||
#include "GPL.h" |
||||
#include "cache.h" |
||||
#include "curlfuncs.h" |
||||
#include "nczarrcommon.h" |
||||
#include <json/json.h> |
||||
#include <variant> |
||||
|
||||
class ZarrTypes: public NcZarrTypes |
||||
{ |
||||
protected: |
||||
template<class VType> class ReadedData |
||||
{ |
||||
//public:
|
||||
using Vec = std::vector<size_t>; |
||||
|
||||
private: |
||||
Vec start, chunkstart; |
||||
ArrCounter mainind, chunkind, inchunkind; |
||||
std::vector<std::unique_ptr<VType[]>> data; |
||||
|
||||
public: |
||||
ReadedData(): mainind(Vec()), chunkind(Vec()), inchunkind(Vec()) {} |
||||
|
||||
ReadedData(size_t N, const size_t* start, const size_t* count, const size_t* csize, std::vector<std::unique_ptr<VType[]>>&& d): |
||||
start(start, start + N), |
||||
chunkstart( |
||||
[](size_t N, const size_t* st, const size_t* cs) |
||||
{ |
||||
Vec out(N); |
||||
for(size_t i = 0; i < N; i++) out[i] = st[i] / cs[i]; |
||||
return out; |
||||
}(N, start, csize)), |
||||
mainind(Vec(count, count + N)), |
||||
chunkind( |
||||
[](size_t N, const size_t* st, const size_t* cn, const size_t* cs) |
||||
{ |
||||
Vec out(N); |
||||
for(size_t i = 0; i < N; i++) out[i] = (st[i] + cn[i]) / cs[i] - st[i] / cs[i] + 1; |
||||
return out; |
||||
}(N, start, count, csize)), |
||||
inchunkind(Vec(csize, csize + N)), |
||||
data(std::move(d)) |
||||
{ |
||||
} |
||||
|
||||
VType operator()(size_t lini) const |
||||
{ |
||||
Vec ind = mainind.Index(lini, mainind.Count()); |
||||
Vec cind(ind.size()), inind(ind.size()); |
||||
|
||||
for(size_t i = 0; i < ind.size(); i++) |
||||
{ |
||||
cind[i] = (ind[i] + start[i]) / inchunkind.Count(i) - chunkstart[i]; // indes of chunk
|
||||
inind[i] = (ind[i] + start[i]) % inchunkind.Count(i); // index inside chunk
|
||||
} |
||||
size_t chunk = chunkind.Index(cind); |
||||
size_t inside = inchunkind.Index(inind); |
||||
return data[chunk][inside]; |
||||
} |
||||
}; |
||||
|
||||
private: |
||||
// Create attribute from json value
|
||||
static AttVT CreateAtt(const Json::Value& val) |
||||
{ |
||||
if(val.type() == Json::intValue) return AttVT{std::in_place_type<int8>, val.asInt64()}; |
||||
if(val.type() == Json::uintValue) return AttVT(std::in_place_type<uint8>, val.asUInt64()); |
||||
if(val.type() == Json::realValue) return AttVT(std::in_place_type<double>, val.asDouble()); |
||||
if(val.type() == Json::stringValue) |
||||
{ |
||||
auto str = val.asString(); |
||||
return AttVT(std::in_place_type<MString>, MString(str.c_str(), str.size())); |
||||
} |
||||
if(val.type() == Json::booleanValue) return AttVT(std::in_place_type<bool>, val.asBool()); |
||||
return AttVT(); |
||||
} |
||||
|
||||
public: |
||||
// Read attributes from .zattrs
|
||||
static auto ReadAtts(const Json::Value& obj) |
||||
{ |
||||
std::vector<Attribute> out; |
||||
if(obj.type() != Json::objectValue) return out; |
||||
const auto keys = obj.getMemberNames(); |
||||
for(const auto& key: keys) |
||||
if(key != "_ARRAY_DIMENSIONS") out.emplace_back(key, CreateAtt(obj[key])); |
||||
return out; |
||||
} |
||||
}; |
||||
|
||||
class ZarrFunctions: public ZarrTypes |
||||
{ |
||||
std::unique_ptr<GenericCache> cache; |
||||
CURLRAII chandle; |
||||
MString url; |
||||
MString proxyurl; |
||||
|
||||
std::vector<std::vector<size_t>> chunks; |
||||
|
||||
// Find variable names in metadata
|
||||
static std::vector<MString> ReadVarNames(const Json::Value& meta); |
||||
|
||||
Error AddVar(const MString& name, const Json::Value& zattrs, const Json::Value& zarray); |
||||
|
||||
protected: |
||||
ZarrFunctions() |
||||
{ |
||||
auto oldprefix = michlib::GPL.UsePrefix("ZARR"); |
||||
cache.reset(CreateCache(michlib::GPL.ParameterSValue("Cache", ""))); |
||||
proxyurl = michlib::GPL.ParameterSValue("Proxy", ""); |
||||
if(proxyurl.Exist()) curl_easy_setopt(chandle, CURLOPT_PROXY, proxyurl.Buf()); |
||||
michlib::GPL.UsePrefix(oldprefix); |
||||
if(!cache) |
||||
{ |
||||
michlib::errmessage("Can't init data cache"); |
||||
cache.reset(new FakeCache); |
||||
} |
||||
} |
||||
|
||||
template<class VType> RetVal<ReadedData<VType>> Read(const MString& var, const size_t* start, const size_t* count) const |
||||
{ |
||||
using Vec = std::vector<size_t>; |
||||
|
||||
size_t ind = FindInd(var, vars); |
||||
const size_t N = vars[ind].NDim(); |
||||
const auto& csize = chunks[ind]; |
||||
|
||||
Vec chunkstart( |
||||
[](size_t N, const size_t* st, const size_t* cs) |
||||
{ |
||||
Vec out(N); |
||||
for(size_t i = 0; i < N; i++) out[i] = st[i] / cs[i]; |
||||
return out; |
||||
}(N, start, csize.data())); |
||||
ArrCounter chunkind( |
||||
[](size_t N, const size_t* st, const size_t* cn, const size_t* cs) |
||||
{ |
||||
Vec out(N); |
||||
for(size_t i = 0; i < N; i++) out[i] = (st[i] + cn[i] - 1) / cs[i] - st[i] / cs[i] + 1; |
||||
return out; |
||||
}(N, start, count, csize.data())); |
||||
|
||||
bool havefill = vars[ind].Fill().index() > 0; |
||||
VType fill = std::visit( |
||||
[](auto v) |
||||
{ |
||||
if constexpr(std::is_convertible_v<decltype(v), VType>) |
||||
return static_cast<VType>(v); |
||||
else |
||||
return std::numeric_limits<VType>::max(); |
||||
}, |
||||
vars[ind].Fill()); |
||||
|
||||
std::vector<std::unique_ptr<VType[]>> cdata; |
||||
|
||||
size_t chunksize = 1; |
||||
for(const auto c: csize) chunksize *= c; |
||||
|
||||
cdata.resize(chunkind.N()); |
||||
|
||||
for(; chunkind; ++chunkind) |
||||
{ |
||||
cdata[chunkind.Index()].reset(new VType[chunksize]); |
||||
auto res = GetChunk(var, chunkind.VIndex(chunkstart), chunksize, sizeof(VType), cdata[chunkind.Index()].get(), havefill ? &fill : nullptr); |
||||
if(!res) return res; |
||||
} |
||||
|
||||
return ReadedData<VType>(N, start, count, csize.data(), std::move(cdata)); |
||||
} |
||||
|
||||
Error GetChunk(const MString& var, const std::vector<size_t>& chunkind, size_t chunksize, size_t elsize, void* data, const void* fill) const; |
||||
|
||||
public: |
||||
Error Open(const MString& product, const MString& dataset, bool time = true); |
||||
}; |
||||
|
||||
using Zarr = NcZarrRead<ZarrFunctions>; |
||||
@ -1 +1 @@
|
||||
Subproject commit a2105413f8a10e8695c85706487d4db5d2f54a61 |
||||
Subproject commit f380988909fadd7f42c7ce09288c4ff3198ca318 |
||||
@ -0,0 +1,222 @@
|
||||
#define MICHLIB_NOSOURCE |
||||
#include "COPERNICUS.h" |
||||
#include "mirrorfuncs.h" |
||||
#include <libxml/parser.h> |
||||
#include <libxml/tree.h> |
||||
|
||||
using michlib::GPL; |
||||
|
||||
RetVal<std::vector<struct FileInfo>> COPERNICUSData::ReadRemoteFileList(const MString& url) const |
||||
{ |
||||
const static MString pref = "COPERNICUSData::ReadRemoteFileList"; |
||||
LIBXML_TEST_VERSION |
||||
|
||||
std::vector<struct FileInfo> out; |
||||
MString bucket, prefix; |
||||
|
||||
// Split url on prefix and bucket
|
||||
{ |
||||
size_t pos = url.Len(); |
||||
size_t count = 0; |
||||
for(size_t i = 0; i < url.Len(); i++) |
||||
{ |
||||
if(url[i] == '/') count++; |
||||
if(count == 4) |
||||
{ |
||||
pos = i; |
||||
break; |
||||
} |
||||
} |
||||
if(pos == url.Len()) return {pref, "Can't parse url: " + url}; |
||||
|
||||
bucket = url.SubStr(1, pos); |
||||
prefix = url.SubStr(pos + 2, url.Len() - pos - 1); |
||||
} |
||||
|
||||
MString cont; |
||||
bool next = true; |
||||
|
||||
CURLRAII chandle; |
||||
while(next) |
||||
{ |
||||
MString url = bucket + "?list-type=2&prefix=" + prefix; |
||||
if(cont.Exist()) url += "&continuation-token=" + cont; |
||||
cont = ""; |
||||
|
||||
auto [data, res] = GetUrl(chandle, url); |
||||
if(res != CURLE_OK) return {pref, MString("Can't download ") + url + ": " + chandle.Err()}; |
||||
|
||||
xmlDocPtr doc = xmlReadMemory(data.Buf(), data.Len(), "data.xml", nullptr, 0); |
||||
if(doc == nullptr) return {pref, MString("Can't download ") + url + ": XML parse error"}; |
||||
auto cur = xmlDocGetRootElement(doc); |
||||
if(cur == nullptr) |
||||
{ |
||||
xmlFreeDoc(doc); |
||||
return {pref, MString("Can't download ") + url + ": empty XML"}; |
||||
} |
||||
if(xmlStrEqual(cur->name, (const xmlChar*)"ListBucketResult") == 0) |
||||
{ |
||||
xmlFreeDoc(doc); |
||||
return {pref, MString("Can't download ") + url + ": unknown XML"}; |
||||
} |
||||
|
||||
for(const auto* n = cur->children; n; n = n->next) |
||||
{ |
||||
if(xmlStrEqual(n->name, (const xmlChar*)"NextContinuationToken") == 1) |
||||
{ |
||||
auto* content = xmlNodeGetContent(n); |
||||
cont = (char*)content; |
||||
xmlFree(content); |
||||
} |
||||
if(xmlStrEqual(n->name, (const xmlChar*)"Contents") == 1) |
||||
{ |
||||
MString fname; |
||||
MDateTime mtime; |
||||
size_t size = 0; |
||||
for(const auto* c = n->children; c; c = c->next) |
||||
{ |
||||
if(xmlStrEqual(c->name, (const xmlChar*)"Key") == 1) |
||||
{ |
||||
auto* content = xmlNodeGetContent(c); |
||||
fname = (char*)content; |
||||
xmlFree(content); |
||||
} |
||||
if(xmlStrEqual(c->name, (const xmlChar*)"LastModified") == 1) |
||||
{ |
||||
auto* content = xmlNodeGetContent(c); |
||||
mtime.FromString((char*)content); |
||||
xmlFree(content); |
||||
} |
||||
if(xmlStrEqual(c->name, (const xmlChar*)"Size") == 1) |
||||
{ |
||||
auto* content = xmlNodeGetContent(c); |
||||
size = MString((char*)content).ToInteger<size_t>(); |
||||
xmlFree(content); |
||||
} |
||||
} |
||||
out.emplace_back(bucket + "/" + fname, fname.SubStr(prefix.Len() + 2, fname.Len() - prefix.Len() - 1), mtime, size); |
||||
} |
||||
} |
||||
xmlFreeDoc(doc); |
||||
next = cont.Exist(); |
||||
} |
||||
|
||||
std::sort(out.begin(), out.end(), [](const struct FileInfo& a, const struct FileInfo& b) { return a.name < b.name; }); |
||||
return out; |
||||
} |
||||
|
||||
Error COPERNICUSData::Mirror(const CLArgs& args) const |
||||
{ |
||||
const static MString pref = "COPERNICUSData::Mirror"; |
||||
|
||||
GPL.UsePrefix("COPERNICUS"); |
||||
|
||||
// Local directory
|
||||
MString mirrorroot = GPL.ParameterSValue("MirrorTo", ""); |
||||
if(!mirrorroot.Exist()) return {pref, "Local mirror directory not specified"}; |
||||
|
||||
if(!args.contains("product")) return {pref, "Copernicus product not specified"}; |
||||
MString prod = args.at("product"); |
||||
CopernicusCatalog cat; |
||||
|
||||
std::vector<MString> dsets; |
||||
if(args.contains("dataset")) |
||||
dsets.push_back(args.at("dataset")); |
||||
else |
||||
{ |
||||
auto dlist = cat.DatasetList(prod); |
||||
if(!dlist) return dlist.Add(pref, "Can't get list of datasets"); |
||||
dsets = dlist.Value(); |
||||
} |
||||
|
||||
michlib::RegExpSimple filter((args.contains("filter") ? args.at("filter") : ".*").Buf()); |
||||
if(filter.Compile() != 0) return Error(pref, MString("Can't compile regular expression ") + filter.RegStr()); |
||||
|
||||
CURLRAII chandle; |
||||
for(const auto& dset: dsets) |
||||
{ |
||||
michlib::message("Mirroring " + dset); |
||||
auto url = cat.DatasetNativeURL(prod, dset); |
||||
if(!url) return {pref, "Can't find data for dataset " + dset + " from product " + prod}; |
||||
|
||||
MString locroot = mirrorroot + "/" + prod + "/" + dset; |
||||
|
||||
auto lfilesret = ReadLocalFileList(locroot); |
||||
if(!lfilesret) return lfilesret.Add(pref, "Can't get local file list"); |
||||
const auto& lfiles = lfilesret.Value(); |
||||
|
||||
auto rfilesret = ReadRemoteFileList(url.Value()); |
||||
if(!rfilesret) return rfilesret.Add(pref, "Can't get remote file list"); |
||||
const auto& rfiles = rfilesret.Value(); |
||||
|
||||
std::vector<size_t> down, rem; |
||||
std::vector<std::pair<size_t, size_t>> upd; |
||||
|
||||
{ |
||||
size_t rpos = 0, lpos = 0; |
||||
while(rpos != rfiles.size() || lpos != lfiles.size()) |
||||
{ |
||||
if(rpos == rfiles.size()) |
||||
while(lpos != lfiles.size()) |
||||
{ |
||||
if(filter.Match(lfiles[lpos].name.Buf())) rem.push_back(lpos); |
||||
lpos++; |
||||
} |
||||
if(lpos == lfiles.size()) |
||||
while(rpos != rfiles.size()) |
||||
{ |
||||
if(filter.Match(rfiles[rpos].name.Buf())) down.push_back(rpos); |
||||
rpos++; |
||||
} |
||||
if(rpos == rfiles.size() || lpos == lfiles.size()) continue; |
||||
|
||||
if(rfiles[rpos].name < lfiles[lpos].name) |
||||
{ |
||||
if(filter.Match(rfiles[rpos].name.Buf())) down.push_back(rpos); |
||||
rpos++; |
||||
} |
||||
else if(lfiles[lpos].name < rfiles[rpos].name) |
||||
{ |
||||
if(filter.Match(lfiles[lpos].name.Buf())) rem.push_back(lpos); |
||||
lpos++; |
||||
} |
||||
else |
||||
{ |
||||
auto delta = rfiles[rpos].mtime.Epoch() - lfiles[lpos].mtime.Epoch(); |
||||
if(delta < 0) delta = -delta; |
||||
if((delta > 0 || rfiles[rpos].size != lfiles[lpos].size) && filter.Match(lfiles[lpos].name.Buf())) upd.emplace_back(rpos, lpos); |
||||
lpos++; |
||||
rpos++; |
||||
} |
||||
} |
||||
} |
||||
|
||||
michlib::message(MString("New files: ") + down.size()); |
||||
michlib::message(MString("Obsolete files: ") + rem.size()); |
||||
michlib::message(MString("Modified files: ") + upd.size()); |
||||
|
||||
for(size_t i = 0; i < down.size(); i++) |
||||
{ |
||||
size_t ri = down[i]; |
||||
auto err = DownloadFile(chandle, rfiles[ri], locroot); |
||||
if(!err) return err.Add(pref, "Can't download file"); |
||||
} |
||||
|
||||
for(size_t i = 0; i < rem.size(); i++) |
||||
{ |
||||
size_t li = rem[i]; |
||||
auto err = RemoveFile(lfiles[li]); |
||||
if(!err) return err.Add(pref, "Can't remove file"); |
||||
} |
||||
|
||||
for(size_t i = 0; i < upd.size(); i++) |
||||
{ |
||||
size_t ri = upd[i].first; |
||||
size_t li = upd[i].second; |
||||
auto err = UpdateFile(chandle, rfiles[ri], lfiles[li], locroot); |
||||
if(!err) return err.Add(pref, "Can't update file"); |
||||
} |
||||
} |
||||
|
||||
return Error(); |
||||
} |
||||
@ -0,0 +1,21 @@
|
||||
#pragma once |
||||
#include "ParseArgs.h" |
||||
#include "copcat.h" |
||||
#include "mdatetime.h" |
||||
|
||||
using michlib::MDateTime; |
||||
using michlib::MString; |
||||
|
||||
class COPERNICUSData |
||||
{ |
||||
// Get remote file list from url
|
||||
RetVal<std::vector<struct FileInfo>> ReadRemoteFileList(const MString& url) const; |
||||
|
||||
public: |
||||
static constexpr const char* name = "COPERNICUS"; |
||||
|
||||
COPERNICUSData() = default; |
||||
|
||||
// Main mirror function
|
||||
Error Mirror(const CLArgs& args) const; |
||||
}; |
||||
@ -0,0 +1,86 @@
|
||||
#define MICHLIB_NOSOURCE |
||||
#include "GRIDVEL.h" |
||||
|
||||
MString GRIDVELData::Info() const |
||||
{ |
||||
if(!isOk()) return ""; |
||||
|
||||
// clang-format off
|
||||
return MString() + |
||||
" Region: (" + u.Xb() + " : " + u.Xe() + ") x (" + u.Yb() + " : " + u.Ye() + ")\n" + |
||||
" Grid: " + u.Nx() + "x" + u.Ny() + "\n" + |
||||
" Supported variables: u, v, U2"; |
||||
// clang-format on
|
||||
} |
||||
|
||||
MString GRIDVELData::Open(const CLArgs& args) |
||||
{ |
||||
MString uname = args.contains("u") ? args.at("u") : ""; |
||||
MString vname = args.contains("v") ? args.at("v") : ""; |
||||
if(!uname.Exist()) return "File with u component not specified"; |
||||
if(!vname.Exist()) return "File with v component not specified"; |
||||
|
||||
unit = args.contains("unit") ? args.at("unit") : "cm/s"; |
||||
fill = args.contains("fill") ? args.at("fill").ToReal() : 1e10; |
||||
|
||||
u.Open(uname); |
||||
if(!u.Opened()) return "Can't open file: " + uname; |
||||
v.Open(vname); |
||||
if(!v.Opened()) return "Can't open file: " + vname; |
||||
if(u != v) return "Files " + uname + " and " + vname + " have different grids"; |
||||
|
||||
return ""; |
||||
} |
||||
|
||||
bool GRIDVELData::Read(const MString& vname, std::map<MString, GRIDVELData::Data>& cache, size_t i) const |
||||
{ |
||||
if(cache.contains(vname)) return true; |
||||
if(!isOk()) return false; |
||||
|
||||
// Only rectangular grids are supported
|
||||
real xs = (u.Xe() - u.Xb()) / (u.Nx() - 1); |
||||
real ys = (u.Ye() - u.Yb()) / (u.Ny() - 1); |
||||
|
||||
Data out(u.Nx(), u.Ny(), u.Xb(), u.Yb(), xs, ys, MString(unit)); |
||||
|
||||
// U and U2 from u and v
|
||||
if(vname == "U" || vname == "U2") |
||||
{ |
||||
bool square = vname == "U2"; |
||||
if(!(Read("u", cache, i) && Read("v", cache, i))) return false; |
||||
cache[vname] = cache.at("u").CopyGrid(); |
||||
auto& U = cache.at(vname); |
||||
const auto& udata = cache.at("u"); |
||||
const auto& vdata = cache.at("v"); |
||||
if(udata.Unit().Exist()) U.SetUnit(square ? ("(" + udata.Unit() + ")2") : udata.Unit()); |
||||
for(size_t ind = 0; ind < U.N(); ind++) |
||||
{ |
||||
if(udata.IsFill(ind) || vdata.IsFill(ind)) |
||||
U.V(ind) = U.Fillval(); |
||||
else |
||||
U.V(ind) = square ? (udata(ind) * udata(ind) + vdata(ind) * vdata(ind)) : michlib::Hypot(udata(ind), vdata(ind)); |
||||
} |
||||
return true; |
||||
} |
||||
|
||||
if(vname == "u" || vname == "v") |
||||
{ |
||||
bool isu = vname == "u"; |
||||
for(IType ix = 0; ix < u.Nx(); ix++) |
||||
for(IType iy = 0; iy < u.Ny(); iy++) |
||||
{ |
||||
if(isu ? (u(ix, iy) == NAN || u(ix, iy) >= 1e10) : (v(ix, iy) == NAN || v(ix, iy) >= 1e10)) |
||||
out.V(ix, iy) = out.Fillval(); |
||||
else |
||||
out.V(ix, iy) = isu ? u(ix, iy) : v(ix, iy); |
||||
} |
||||
|
||||
if(out) |
||||
{ |
||||
cache[vname] = std::move(out); |
||||
return true; |
||||
} |
||||
} |
||||
|
||||
return false; |
||||
} |
||||
@ -0,0 +1,50 @@
|
||||
#pragma once |
||||
#include "DataAdapters/gridfile.h" |
||||
#include "mdatetime.h" |
||||
#include "simple2ddata.h" |
||||
|
||||
using michlib::MDateTime; |
||||
|
||||
class GRIDVELData |
||||
{ |
||||
michlib::GridFile u,v; |
||||
MString unit; |
||||
real fill; |
||||
|
||||
using IType=decltype(u.Nx()); |
||||
|
||||
public: |
||||
static constexpr const char* name = "GRIDVEL"; |
||||
|
||||
static constexpr const char* disabledactions = "genintfile"; |
||||
|
||||
using Data = Simple2DData; |
||||
|
||||
GRIDVELData() = default; |
||||
|
||||
MString Info() const; |
||||
// TODO: RetVal
|
||||
MString Open(const CLArgs& args); |
||||
|
||||
bool isOk() const { return u.Opened() &&v.Opened(); } |
||||
|
||||
size_t NTimes() const { return 1; } |
||||
|
||||
MDateTime Time(size_t i) const |
||||
{ |
||||
return MDateTime(); |
||||
} |
||||
|
||||
//time_t Timestep() const { return isOk() ? (times[1] - times[0]) : 0; }
|
||||
|
||||
explicit operator bool() const { return u.Opened() &&v.Opened(); } |
||||
|
||||
VarPresence CheckVar(const MString& vname) const |
||||
{ |
||||
if(vname == "u" || vname == "v") return VarPresence::INTERNAL; |
||||
if(vname == "U" || vname == "U2") return VarPresence::DERIVED; |
||||
return VarPresence::NONE; |
||||
} |
||||
|
||||
bool Read(const MString& vname, std::map<MString, Data>& cache, size_t i) const; |
||||
}; |
||||
@ -0,0 +1,129 @@
|
||||
#define MICHLIB_NOSOURCE |
||||
#include "TSCDATA.h" |
||||
|
||||
MString TSCDATAData::Open(const CLArgs& args) |
||||
{ |
||||
if(!args.contains("dataset")) return "path to data not specified"; |
||||
MString dataset = args.at("dataset"); |
||||
|
||||
michlib::NCFileA newnc; |
||||
std::vector<MString> newvnames, newlnames; |
||||
MString newhistory; |
||||
|
||||
newnc.Reset(dataset); |
||||
if(!newnc) return "Can't open file " + dataset; |
||||
|
||||
auto head = newnc.Header(); |
||||
if(head.Dimensions().size() != 2) return "Unsupported number of dimensions"; |
||||
if((head.Dimensions()[0].Name() != "longitude" || head.Dimensions()[1].Name() != "latitude") && |
||||
(head.Dimensions()[1].Name() != "longitude" || head.Dimensions()[0].Name() != "latitude")) |
||||
return "Unsupported dimensions names"; |
||||
|
||||
if(head.Dimensions()[0].Name() == "longitude") |
||||
{ |
||||
nx = head.Dimensions()[0].Len(); |
||||
ny = head.Dimensions()[1].Len(); |
||||
} |
||||
else |
||||
{ |
||||
ny = head.Dimensions()[0].Len(); |
||||
nx = head.Dimensions()[1].Len(); |
||||
} |
||||
|
||||
{ |
||||
bool lonfound = false, latfound = false; |
||||
for(const auto& v: head.Variables()) |
||||
if(v.Dimensions().size() == 1 && v.Type().Id() == NC_FLOAT) |
||||
{ |
||||
lonfound = lonfound || (v.Dimensions()[0].Name() == "longitude" && v.Name() == "longitude"); |
||||
latfound = latfound || (v.Dimensions()[0].Name() == "latitude" && v.Name() == "latitude"); |
||||
} |
||||
if(!lonfound) return "Longitude not found"; |
||||
if(!latfound) return "Latitude not found"; |
||||
} |
||||
|
||||
for(const auto& v: head.Variables()) |
||||
{ |
||||
if(v.Dimensions().size() != 2) continue; |
||||
if(v.Type().Id() != NC_FLOAT) continue; |
||||
if((v.Dimensions()[0].Name() != "longitude" || v.Dimensions()[1].Name() != "latitude") && (v.Dimensions()[1].Name() != "longitude" || v.Dimensions()[0].Name() != "latitude")) |
||||
continue; |
||||
newvnames.push_back(v.Name()); |
||||
auto lname = newnc.A<MString>(v.Name(), "long_name"); |
||||
newlnames.push_back(lname ? lname.Get() : ""); |
||||
} |
||||
if(newvnames.size() == 0) return "No variables found"; |
||||
|
||||
{ |
||||
auto his = newnc.A<MString>("history"); |
||||
if(his) newhistory = his; |
||||
} |
||||
|
||||
history = std::move(newhistory); |
||||
vnames = std::move(newvnames); |
||||
lnames = std::move(newlnames); |
||||
nc = std::move(newnc); |
||||
return ""; |
||||
} |
||||
|
||||
MString TSCDATAData::Info() const |
||||
{ |
||||
if(!nc) return ""; |
||||
|
||||
MString out; |
||||
|
||||
out += MString("Dimensions: ") + nx + " X " + ny + "\n"; |
||||
out += MString("Variables: "); |
||||
for(size_t i = 0; i < vnames.size(); i++) out += ((i == 0) ? "" : ", ") + vnames[i]; |
||||
out += "\n"; |
||||
auto his = nc.A<MString>("history"); |
||||
if(his) out += "Creator: " + his.Get() + "\n"; |
||||
return out; |
||||
} |
||||
|
||||
std::vector<TSCDATAData::DataType> TSCDATAData::ReadLons() const |
||||
{ |
||||
std::vector<DataType> out; |
||||
auto lons = nc.V<DataType>("longitude"); |
||||
if(lons) |
||||
{ |
||||
out.resize(lons.DimLen(0)); |
||||
for(size_t i = 0; i < out.size(); i++) out[i] = lons(i); |
||||
} |
||||
return out; |
||||
} |
||||
|
||||
std::vector<TSCDATAData::DataType> TSCDATAData::ReadLats() const |
||||
{ |
||||
std::vector<DataType> out; |
||||
auto lats = nc.V<DataType>("latitude"); |
||||
if(lats) |
||||
{ |
||||
out.resize(lats.DimLen(0)); |
||||
for(size_t i = 0; i < out.size(); i++) out[i] = lats(i); |
||||
} |
||||
return out; |
||||
} |
||||
|
||||
std::vector<TSCDATAData::DataType> TSCDATAData::ReadVar(const MString& name) const |
||||
{ |
||||
std::vector<DataType> out; |
||||
bool havevar = false; |
||||
for(size_t i = 0; i < vnames.size(); i++) havevar = havevar || (vnames[i] == name); |
||||
if(!havevar) return out; |
||||
|
||||
auto var = nc.V<DataType>(name, "latitude", "longitude"); |
||||
if(var) |
||||
{ |
||||
out.resize(nx * ny); |
||||
for(size_t iy = 0; iy < ny; iy++) |
||||
for(size_t ix = 0; ix < nx; ix++) out[iy * nx + ix] = var(iy, ix); |
||||
} |
||||
return out; |
||||
} |
||||
|
||||
TSCDATAData::DataType TSCDATAData::FillVal(const MString& name) const |
||||
{ |
||||
auto fill = nc.A<DataType>(name, "_FillValue"); |
||||
return fill ? fill.Get() : 0.0; |
||||
} |
||||
@ -0,0 +1,35 @@
|
||||
#pragma once |
||||
#include "DataAdapters/ncfilealt.h" |
||||
#include "ParseArgs.h" |
||||
#include "simple2ddata.h" |
||||
|
||||
class TSCDATAData |
||||
{ |
||||
michlib::NCFileA nc; |
||||
std::vector<MString> vnames, lnames; |
||||
size_t nx, ny; |
||||
MString history; |
||||
|
||||
using DataType = float; |
||||
|
||||
public: |
||||
static constexpr const char* name = "TSCDATA"; |
||||
|
||||
MString Info() const; |
||||
MString Open(const CLArgs& args); |
||||
|
||||
std::vector<DataType> ReadLons() const; |
||||
std::vector<DataType> ReadLats() const; |
||||
std::vector<DataType> ReadVar(const MString& name) const; |
||||
|
||||
const auto& VarNames() const { return vnames; } |
||||
const auto& LongNames() const { return lnames; } |
||||
|
||||
size_t Nx() const { return nx; } |
||||
size_t Ny() const { return ny; } |
||||
size_t NVar() const { return vnames.size(); } |
||||
|
||||
const auto& History() const { return history; } |
||||
|
||||
DataType FillVal(const MString& name) const; |
||||
}; |
||||
@ -0,0 +1,526 @@
|
||||
#define MICHLIB_NOSOURCE |
||||
#include "VYLET.h" |
||||
|
||||
MString VYLETData::Info() const |
||||
{ |
||||
if(!vylet) return ""; |
||||
MString out; |
||||
michlib::CompiledParser xy2lon, xy2lat; |
||||
real x, y; |
||||
michlib::ParserVars pv; |
||||
|
||||
pv["x"] = &x; |
||||
pv["y"] = &y; |
||||
|
||||
vylet->UsePrefix("Datafile_Info"); |
||||
MString xy2lonstr = vylet->ParameterSValue("xy2lon", "%y"); |
||||
MString xy2latstr = vylet->ParameterSValue("xy2lat", "%x"); |
||||
ArifmeticCompiler(xy2lonstr, xy2lon, &pv); |
||||
ArifmeticCompiler(xy2latstr, xy2lat, &pv); |
||||
|
||||
real lonb, latb, lone, late; |
||||
|
||||
vylet->UsePrefix(""); |
||||
x = vylet->ParameterRValue("x0", 0.0); |
||||
y = vylet->ParameterRValue("y0", 0.0); |
||||
xy2lon.Run(lonb); |
||||
xy2lat.Run(latb); |
||||
x = vylet->ParameterRValue("x1", 0.0); |
||||
y = vylet->ParameterRValue("y1", 0.0); |
||||
xy2lon.Run(lone); |
||||
xy2lat.Run(late); |
||||
|
||||
out += "Start time: " + start.ToString() + "\n"; |
||||
out += "End time: " + end.ToString() + " " + (invtime ? "(backward integration)" : "(forward integration)") + "\n"; |
||||
out += "Region: (" + MString(lonb) + " : " + lone + ") x (" + latb + " : " + late + ")\n"; |
||||
out += "Grid:" + MString(" ") + lons.size() + "x" + lats.size() + "\n"; |
||||
out += HasCoast() ? "Coast: checked\n" : "Coast: not checked\n"; |
||||
out += LengthFast() ? "Trajectory length: fast calculation\n" : "Trajectory length: exact calculation\n"; |
||||
|
||||
out += "Borders: "; |
||||
bool needcomma = false; |
||||
if(Left() >= 0.0) |
||||
{ |
||||
real lon; |
||||
if(needcomma) out += ", "; |
||||
out += "left="; |
||||
x = Left(); |
||||
xy2lon.Run(lon); |
||||
out += lon; |
||||
needcomma = true; |
||||
} |
||||
if(Right() <= maxx) |
||||
{ |
||||
real lon; |
||||
if(needcomma) out += ", "; |
||||
out += "right="; |
||||
x = Right(); |
||||
xy2lon.Run(lon); |
||||
out += lon; |
||||
needcomma = true; |
||||
} |
||||
if(Down() >= 0.0) |
||||
{ |
||||
real lat; |
||||
if(needcomma) out += ", "; |
||||
out += "bottom="; |
||||
y = Down(); |
||||
xy2lat.Run(lat); |
||||
out += lat; |
||||
needcomma = true; |
||||
} |
||||
if(Up() <= maxy) |
||||
{ |
||||
real lat; |
||||
if(needcomma) out += ", "; |
||||
out += "bottom="; |
||||
y = Up(); |
||||
xy2lat.Run(lat); |
||||
out += lat; |
||||
needcomma = true; |
||||
} |
||||
if(!needcomma) out += "no"; |
||||
out += "\n"; |
||||
|
||||
MString vars = "vylD"; |
||||
if(HasLyap()) vars += ", vylL"; |
||||
if(HasTime()) vars += ", vylT"; |
||||
vars += ", vylNx, vylNy, vylRx, vylRy, vylEx, vylEy, vylPhip, vylPhim, vylPhit"; |
||||
if(HasLyap()) vars += ", vyldS"; |
||||
vars += ", vylAngle, vylLen"; |
||||
if(HasTime()) vars += ", vylTmask"; |
||||
|
||||
out += "Supported variables: " + vars + "\n"; |
||||
|
||||
return out; |
||||
} |
||||
|
||||
MString VYLETData::Open(const CLArgs& args) |
||||
{ |
||||
if(!args.contains("dataset")) return "path to data not specified"; |
||||
MString dataset = args.at("dataset"); |
||||
decltype(vylet) nvylet; |
||||
|
||||
michlib::RegExpSimple havex("%x"), havey("%y"); |
||||
|
||||
nvylet.reset(new michlib::BFileR); |
||||
if(nvylet->Open(dataset) != ERR_NOERR) return "Can't open file " + dataset; |
||||
|
||||
nvylet->UsePrefix("ProgramInfo"); |
||||
if(nvylet->ParameterSValue("Task", "") != "AdvInt:Vylet") return "File " + dataset + " is not vylet file"; |
||||
|
||||
nvylet->UsePrefix(""); |
||||
if(nvylet->ParameterSValue("nettype", "") != "SQUARE") return "File " + dataset + " have unsupported net type"; |
||||
MString method = nvylet->ParameterSValue("Method", ""); |
||||
if(!(method == "Bicubic" || method == "BicubicL" || method == "BicubicI" || method == "BicubicIL")) return "File " + dataset + " have unsupported integration method"; |
||||
invtime = (method == "BicubicI" || method == "BicubicIL"); |
||||
|
||||
nvylet->UsePrefix("Datafile_Info"); |
||||
MString xy2lonstr = nvylet->ParameterSValue("xy2lon", "%y"); |
||||
MString xy2latstr = nvylet->ParameterSValue("xy2lat", "%x"); |
||||
if(havey.Match(xy2lonstr.Buf()) || havex.Match(xy2latstr.Buf())) return "File " + dataset + " have unsupported grid"; |
||||
auto res = ref.FromString(nvylet->ParameterSValue(invtime ? "EndDate" : "BeginDate", "")); |
||||
if(!res) return "Can't read reference time"; |
||||
|
||||
nvylet->UsePrefix(""); |
||||
start = R2Time(nvylet->ParameterRValue("tbeg", 0.0)); |
||||
end = R2Time(nvylet->ParameterRValue("tmax", 0.0)); |
||||
|
||||
nvylet->UsePrefix("Datafile"); |
||||
auto dx = nvylet->ParameterRValue("dx", 0.0); |
||||
auto dy = nvylet->ParameterRValue("dy", 0.0); |
||||
auto nx = nvylet->ParameterUValue("nx", 0); |
||||
auto ny = nvylet->ParameterUValue("ny", 0); |
||||
maxx = dx * (nx - 1); |
||||
maxy = dy * (ny - 1); |
||||
|
||||
michlib::CompiledParser xy2lon, xy2lat; |
||||
real x, y; |
||||
michlib::ParserVars pv; |
||||
|
||||
pv["x"] = &x; |
||||
pv["y"] = &y; |
||||
|
||||
ArifmeticCompiler(xy2lonstr, xy2lon, &pv); |
||||
ArifmeticCompiler(xy2latstr, xy2lat, &pv); |
||||
|
||||
nvylet->UsePrefix(""); |
||||
auto nlon = nvylet->ParameterUValue("Nx", 0) + 1; |
||||
auto nlat = nvylet->ParameterUValue("Ny", 0) + 1; |
||||
lons.resize(nlon); |
||||
lats.resize(nlat); |
||||
|
||||
elon.resize(nlon * nlat); |
||||
elat.resize(nlon * nlat); |
||||
|
||||
for(size_t iy = 0; iy < lats.size(); iy++) |
||||
for(size_t ix = 0; ix < lons.size(); ix++) |
||||
{ |
||||
x = (*nvylet)[2][iy * lons.size() + ix]; |
||||
y = (*nvylet)[3][iy * lons.size() + ix]; |
||||
xy2lon.Run(elon[iy * lons.size() + ix]); |
||||
xy2lat.Run(elat[iy * lons.size() + ix]); |
||||
} |
||||
|
||||
for(size_t ix = 0; ix < nlon; ix++) |
||||
{ |
||||
x = (*nvylet)[0][ix]; |
||||
y = (*nvylet)[1][ix]; |
||||
xy2lon.Run(lons[ix]); |
||||
} |
||||
|
||||
for(size_t iy = 0; iy < nlat; iy++) |
||||
{ |
||||
x = (*nvylet)[0][iy * nlon]; |
||||
y = (*nvylet)[1][iy * nlon]; |
||||
xy2lat.Run(lats[iy]); |
||||
} |
||||
|
||||
vylet = std::move(nvylet); |
||||
return ""; |
||||
} |
||||
|
||||
bool VYLETData::Read(const MString& vname, std::map<MString, VYLETData::Data>& cache, size_t tind) const |
||||
{ |
||||
if(tind != 0) return false; |
||||
if(cache.contains(vname)) return true; |
||||
Data out; |
||||
|
||||
if(vname == "vylD") out = ReadD(); |
||||
if(vname == "vylL") out = ReadL(); |
||||
if(vname == "vylT") out = ReadT(); |
||||
if(vname == "vylNx") out = ReadNx(); |
||||
if(vname == "vylNy") out = ReadNy(); |
||||
if(vname == "vylRx") out = ReadRx(); |
||||
if(vname == "vylRy") out = ReadRy(); |
||||
if(vname == "vylEx") out = ReadEx(); |
||||
if(vname == "vylEy") out = ReadEy(); |
||||
if(vname == "vylPhip") out = ReadPhip(); |
||||
if(vname == "vylPhim") out = ReadPhim(); |
||||
if(vname == "vylPhit") out = ReadPhit(); |
||||
if(vname == "vyldS") out = ReaddS(); |
||||
if(vname == "vylAngle") out = ReadAngle(); |
||||
if(vname == "vylLen") out = ReadLen(); |
||||
if(vname == "vylTmask") out = ReadTmask(); |
||||
|
||||
if(!out) return false; |
||||
cache[vname] = std::move(out); |
||||
return true; |
||||
} |
||||
|
||||
VYLETData::Data VYLETData::ReadD() const |
||||
{ |
||||
Data out(lons, lats); |
||||
if(!vylet) return Data(); |
||||
|
||||
const real xl = Left(); |
||||
const real xr = Right(); |
||||
const real yd = Down(); |
||||
const real yu = Up(); |
||||
|
||||
real x, y; |
||||
|
||||
for(size_t iy = 0; iy < lats.size(); iy++) |
||||
for(size_t ix = 0; ix < lons.size(); ix++) |
||||
{ |
||||
out(ix, iy) = 0.0; |
||||
x = (*vylet)[2][iy * lons.size() + ix]; |
||||
y = (*vylet)[3][iy * lons.size() + ix]; |
||||
if(yd > 0.0 && y < yd) out(ix, iy) = -1.0; |
||||
if(xl > 0.0 && x < xl) out(ix, iy) = -2.0; |
||||
if(yu < maxy && y > yu) out(ix, iy) = -3.0; |
||||
if(xr < maxx && x > xr) out(ix, iy) = -4.0; |
||||
if(out(ix, iy) >= 0.0) |
||||
out(ix, iy) = 6371.0 * michlib::GCD(M_PI * lons[ix] / 180.0, M_PI * lats[iy] / 180.0, M_PI * elon[iy * lons.size() + ix] / 180.0, M_PI * elat[iy * lons.size() + ix] / 180.0); |
||||
} |
||||
|
||||
out.SetUnit("km"); |
||||
out.SetLongName("Distance between start and end points"); |
||||
out.SetComment("Special values: -1.0 - trajectory leave region via south border, -2.0 - trajectory leave region via west border, -3.0 - trajectory leave region via north border, " |
||||
"-4.0 - trajectory leave region via east border"); |
||||
return out; |
||||
} |
||||
|
||||
VYLETData::Data VYLETData::ReadL() const |
||||
{ |
||||
Data out(lons, lats); |
||||
if(!vylet) return Data(); |
||||
|
||||
auto lcol = Name2ColNum("Log(G.sigma1)"); |
||||
auto tcol = Name2ColNum("time"); |
||||
if(lcol == 0 || tcol == 0) return Data(); |
||||
|
||||
MDateTime time; |
||||
real lambda, days; |
||||
|
||||
for(size_t iy = 0; iy < lats.size(); iy++) |
||||
for(size_t ix = 0; ix < lons.size(); ix++) |
||||
{ |
||||
time = R2Time((*vylet)[tcol - 1][iy * lons.size() + ix]); |
||||
lambda = (*vylet)[lcol - 1][iy * lons.size() + ix]; |
||||
days = (time - start).D(); |
||||
out(ix, iy) = lambda / days; |
||||
} |
||||
|
||||
out.SetUnit("days-1"); |
||||
out.SetLongName("Lyapunov exponent"); |
||||
return out; |
||||
} |
||||
|
||||
VYLETData::Data VYLETData::ReadT() const |
||||
{ |
||||
Data out(lons, lats); |
||||
if(!vylet) return Data(); |
||||
|
||||
auto tcol = Name2ColNum("time"); |
||||
if(tcol == 0) return Data(); |
||||
|
||||
const real xl = Left(); |
||||
const real xr = Right(); |
||||
const real yd = Down(); |
||||
const real yu = Up(); |
||||
|
||||
vylet->UsePrefix(""); |
||||
const real acc = vylet->ParameterRValue("accuracy", 1.0); |
||||
const real tstep = 2.0 * M_PI * 1000.0 / acc; |
||||
|
||||
const real maxdays = (end - start).D(); |
||||
|
||||
MDateTime time; |
||||
real days; |
||||
real x, y; |
||||
bool inside, maxtime; |
||||
|
||||
for(size_t iy = 0; iy < lats.size(); iy++) |
||||
for(size_t ix = 0; ix < lons.size(); ix++) |
||||
{ |
||||
x = (*vylet)[2][iy * lons.size() + ix]; |
||||
y = (*vylet)[3][iy * lons.size() + ix]; |
||||
time = R2Time((*vylet)[tcol - 1][iy * lons.size() + ix]); |
||||
days = (time - start).D(); |
||||
if(days <= tstep * 1.5) days = 0.0; |
||||
inside = x > xl && x < xr && y > yd && y < yu; |
||||
maxtime = days >= maxdays - tstep * 0.5; |
||||
if(maxtime) days = maxdays; |
||||
if(inside && !maxtime) days = -days; |
||||
out(ix, iy) = days; |
||||
} |
||||
|
||||
out.SetUnit("days"); |
||||
out.SetLongName("Time to reach border or coast"); |
||||
out.SetComment("Special values: 0.0 - initial point on the land, " + MString(maxdays) + " - trajectory don't reach border, negative values - trajectory beached on the coast"); |
||||
return out; |
||||
} |
||||
|
||||
VYLETData::Data VYLETData::ReadNx() const |
||||
{ |
||||
Data out(lons, lats); |
||||
if(!vylet) return Data(); |
||||
|
||||
auto ncol = Name2ColNum("nx=0"); |
||||
if(ncol == 0) return Data(); |
||||
|
||||
for(size_t iy = 0; iy < lats.size(); iy++) |
||||
for(size_t ix = 0; ix < lons.size(); ix++) out(ix, iy) = (*vylet)[ncol - 1][iy * lons.size() + ix]; |
||||
|
||||
out.SetLongName("Number of moments u=0"); |
||||
return out; |
||||
} |
||||
|
||||
VYLETData::Data VYLETData::ReadNy() const |
||||
{ |
||||
Data out(lons, lats); |
||||
if(!vylet) return Data(); |
||||
|
||||
auto ncol = Name2ColNum("ny=0"); |
||||
if(ncol == 0) return Data(); |
||||
|
||||
for(size_t iy = 0; iy < lats.size(); iy++) |
||||
for(size_t ix = 0; ix < lons.size(); ix++) out(ix, iy) = (*vylet)[ncol - 1][iy * lons.size() + ix]; |
||||
|
||||
out.SetLongName("Number of moments v=0"); |
||||
return out; |
||||
} |
||||
|
||||
VYLETData::Data VYLETData::ReadRx() const |
||||
{ |
||||
Data out(lons, lats); |
||||
if(!vylet) return Data(); |
||||
|
||||
for(size_t iy = 0; iy < lats.size(); iy++) |
||||
for(size_t ix = 0; ix < lons.size(); ix++) out(ix, iy) = elon[iy * lons.size() + ix] - lons[ix]; |
||||
|
||||
out.SetLongName("Displacement in zonal direction"); |
||||
|
||||
return out; |
||||
} |
||||
|
||||
VYLETData::Data VYLETData::ReadRy() const |
||||
{ |
||||
Data out(lons, lats); |
||||
if(!vylet) return Data(); |
||||
|
||||
for(size_t iy = 0; iy < lats.size(); iy++) |
||||
for(size_t ix = 0; ix < lons.size(); ix++) out(ix, iy) = elat[iy * lons.size() + ix] - lats[iy]; |
||||
|
||||
out.SetLongName("Displacement in meridional direction"); |
||||
|
||||
return out; |
||||
} |
||||
|
||||
VYLETData::Data VYLETData::ReadEx() const |
||||
{ |
||||
Data out(lons, lats); |
||||
if(!vylet) return Data(); |
||||
|
||||
for(size_t iy = 0; iy < lats.size(); iy++) |
||||
for(size_t ix = 0; ix < lons.size(); ix++) out(ix, iy) = elon[iy * lons.size() + ix]; |
||||
|
||||
out.SetLongName("Final longitude"); |
||||
|
||||
return out; |
||||
} |
||||
|
||||
VYLETData::Data VYLETData::ReadEy() const |
||||
{ |
||||
Data out(lons, lats); |
||||
if(!vylet) return Data(); |
||||
|
||||
for(size_t iy = 0; iy < lats.size(); iy++) |
||||
for(size_t ix = 0; ix < lons.size(); ix++) out(ix, iy) = elat[iy * lons.size() + ix]; |
||||
|
||||
out.SetLongName("Final latitude"); |
||||
|
||||
return out; |
||||
} |
||||
|
||||
VYLETData::Data VYLETData::ReadPhip() const |
||||
{ |
||||
Data out(lons, lats); |
||||
if(!vylet) return Data(); |
||||
|
||||
auto ncol = Name2ColNum("nphip"); |
||||
if(ncol == 0) return Data(); |
||||
|
||||
for(size_t iy = 0; iy < lats.size(); iy++) |
||||
for(size_t ix = 0; ix < lons.size(); ix++) out(ix, iy) = (*vylet)[ncol - 1][iy * lons.size() + ix]; |
||||
|
||||
out.SetLongName("Number of counterclockwise rotations"); |
||||
return out; |
||||
} |
||||
|
||||
VYLETData::Data VYLETData::ReadPhim() const |
||||
{ |
||||
Data out(lons, lats); |
||||
if(!vylet) return Data(); |
||||
|
||||
auto ncol = Name2ColNum("nphim"); |
||||
if(ncol == 0) return Data(); |
||||
|
||||
for(size_t iy = 0; iy < lats.size(); iy++) |
||||
for(size_t ix = 0; ix < lons.size(); ix++) out(ix, iy) = (*vylet)[ncol - 1][iy * lons.size() + ix]; |
||||
|
||||
out.SetLongName("Number of clockwise rotations"); |
||||
return out; |
||||
} |
||||
|
||||
VYLETData::Data VYLETData::ReadPhit() const |
||||
{ |
||||
Data out(lons, lats); |
||||
if(!vylet) return Data(); |
||||
|
||||
auto pcol = Name2ColNum("nphip"); |
||||
auto mcol = Name2ColNum("nphim"); |
||||
if(pcol == 0 || mcol == 0) return Data(); |
||||
|
||||
for(size_t iy = 0; iy < lats.size(); iy++) |
||||
for(size_t ix = 0; ix < lons.size(); ix++) out(ix, iy) = (*vylet)[pcol - 1][iy * lons.size() + ix] - (*vylet)[mcol - 1][iy * lons.size() + ix]; |
||||
|
||||
out.SetLongName("Difference between the number of rotations counterclockwise and clockwise"); |
||||
return out; |
||||
} |
||||
|
||||
VYLETData::Data VYLETData::ReaddS() const |
||||
{ |
||||
Data out(lons, lats); |
||||
if(!vylet) return Data(); |
||||
|
||||
auto l1col = Name2ColNum("Log(G.sigma1)"); |
||||
auto l2col = Name2ColNum("Log(G.sigma2)"); |
||||
if(l1col == 0 || l2col == 0) return Data(); |
||||
|
||||
for(size_t iy = 0; iy < lats.size(); iy++) |
||||
for(size_t ix = 0; ix < lons.size(); ix++) out(ix, iy) = (*vylet)[l1col - 1][iy * lons.size() + ix] + (*vylet)[l2col - 1][iy * lons.size() + ix]; |
||||
|
||||
out.SetLongName("Logarithm of area change multiplier"); |
||||
return out; |
||||
} |
||||
|
||||
VYLETData::Data VYLETData::ReadAngle() const |
||||
{ |
||||
Data out(lons, lats); |
||||
if(!vylet) return Data(); |
||||
|
||||
auto ncol = Name2ColNum("angle"); |
||||
if(ncol == 0) return Data(); |
||||
|
||||
for(size_t iy = 0; iy < lats.size(); iy++) |
||||
for(size_t ix = 0; ix < lons.size(); ix++) out(ix, iy) = (*vylet)[ncol - 1][iy * lons.size() + ix] / (2.0 * M_PI); |
||||
|
||||
out.SetLongName("Total rotation angle normalized to 2π"); |
||||
return out; |
||||
} |
||||
|
||||
VYLETData::Data VYLETData::ReadLen() const |
||||
{ |
||||
Data out(lons, lats); |
||||
if(!vylet) return Data(); |
||||
|
||||
auto ncol = Name2ColNum("length"); |
||||
if(ncol == 0) return Data(); |
||||
bool sisangle = !LengthFast(); |
||||
real mul = sisangle ? (6371.0 * M_PI / (60.0 * 180.0)) : 1.0; |
||||
|
||||
for(size_t iy = 0; iy < lats.size(); iy++) |
||||
for(size_t ix = 0; ix < lons.size(); ix++) out(ix, iy) = mul * (*vylet)[ncol - 1][iy * lons.size() + ix]; |
||||
|
||||
if(sisangle) |
||||
{ |
||||
out.SetUnit("km"); |
||||
out.SetLongName("Trajectory length"); |
||||
} |
||||
out.SetLongName("Trajectory length (in abstract units)"); |
||||
return out; |
||||
} |
||||
|
||||
VYLETData::Data VYLETData::ReadTmask() const |
||||
{ |
||||
Data out(lons, lats); |
||||
if(!vylet) return Data(); |
||||
|
||||
auto tcol = Name2ColNum("time"); |
||||
if(tcol == 0) return Data(); |
||||
|
||||
vylet->UsePrefix(""); |
||||
const real acc = vylet->ParameterRValue("accuracy", 1.0); |
||||
const real tstep = 2.0 * M_PI * 1000.0 / acc; |
||||
|
||||
const real maxdays = (end - start).D(); |
||||
|
||||
MDateTime time; |
||||
real days; |
||||
bool maxtime; |
||||
|
||||
for(size_t iy = 0; iy < lats.size(); iy++) |
||||
for(size_t ix = 0; ix < lons.size(); ix++) |
||||
{ |
||||
time = R2Time((*vylet)[tcol - 1][iy * lons.size() + ix]); |
||||
days = (time - start).D(); |
||||
maxtime = days >= maxdays - tstep * 0.5; |
||||
out(ix, iy) = maxtime ? 1.0 : NAN; |
||||
} |
||||
|
||||
out.SetLongName("Flag"); |
||||
out.SetComment("Values: 1.0 - the trajectory did not reach either the coast or the borders, NaN - overwise"); |
||||
return out; |
||||
} |
||||
@ -0,0 +1,145 @@
|
||||
#pragma once |
||||
#include "BFileR.h" |
||||
#include "ParseArgs.h" |
||||
#include "mdatetime.h" |
||||
#include "mregex.h" |
||||
#include "simple2ddata.h" |
||||
#include <memory> |
||||
#include <set> |
||||
|
||||
using michlib::M_PI; |
||||
using michlib::MDateTime; |
||||
using michlib::real; |
||||
using michlib::Round; |
||||
|
||||
class VYLETData |
||||
{ |
||||
std::unique_ptr<michlib::BFileR> vylet; |
||||
MDateTime start, end, ref; |
||||
bool invtime; |
||||
real maxx, maxy; |
||||
std::vector<real> lons, lats; |
||||
std::vector<real> elon, elat; |
||||
|
||||
MDateTime R2Time(real r) const |
||||
{ |
||||
auto sec = static_cast<time_t>(Round(r * MDateTime::secondsperday)); |
||||
return ref + (invtime ? -sec : sec); |
||||
}; |
||||
|
||||
public: |
||||
using Data = Rect2DData; |
||||
|
||||
private: |
||||
auto Name2ColNum(const MString& name) const |
||||
{ |
||||
auto nc = vylet->Columns(); |
||||
decltype(nc) i; |
||||
for(i = 1; i <= nc; i++) |
||||
if(vylet->ColumnName(i) == name) return i; |
||||
i = 0; |
||||
return i; |
||||
} |
||||
|
||||
Data ReadD() const; |
||||
Data ReadL() const; |
||||
Data ReadT() const; |
||||
Data ReadNx() const; |
||||
Data ReadNy() const; |
||||
Data ReadRx() const; |
||||
Data ReadRy() const; |
||||
Data ReadEx() const; |
||||
Data ReadEy() const; |
||||
Data ReadPhip() const; |
||||
Data ReadPhim() const; |
||||
Data ReadPhit() const; |
||||
Data ReaddS() const; |
||||
Data ReadAngle() const; |
||||
Data ReadLen() const; |
||||
Data ReadTmask() const; |
||||
|
||||
public: |
||||
static constexpr const char* name = "VYLET"; |
||||
|
||||
static constexpr const char* disabledactions = "genintfile uv"; |
||||
|
||||
MString DefaultVars() const |
||||
{ |
||||
MString vars = "vylD,vylRx,vylRy,vylEx,vylEy,vylAngle,vylLen"; |
||||
if(HasLyap()) vars += ",vylL"; |
||||
if(HasTime()) vars += ",vylT"; |
||||
return vars; |
||||
} |
||||
|
||||
VarPresence CheckVar(const MString& vname) const |
||||
{ |
||||
std::set<MString> vars{"vylD", "vylNx", "vylNy", "vylRx", "vylRy", "vylEx", "vylEy", "vylPhip", "vylPhim", "vylPhit", "vylAngle", "vylLen"}; |
||||
if(HasLyap()) vars.insert("vylL"); |
||||
if(HasLyap()) vars.insert("vyldS"); |
||||
if(HasTime()) vars.insert("vylT"); |
||||
if(HasTime()) vars.insert("vylTmask"); |
||||
return vars.contains(vname) ? VarPresence::DERIVED : VarPresence::NONE; |
||||
} |
||||
|
||||
bool Read(const MString& vname, std::map<MString, Data>& cache, size_t tind) const; |
||||
|
||||
size_t NTimes() const { return vylet ? 1 : 0; } |
||||
|
||||
MDateTime Time(size_t i) const |
||||
{ |
||||
if(i == 0 && vylet) return start; |
||||
return MDateTime(); |
||||
} |
||||
|
||||
bool HasLyap() const |
||||
{ |
||||
if(!vylet) return false; |
||||
vylet->UsePrefix(""); |
||||
MString method = vylet->ParameterSValue("Method", ""); |
||||
return method == "BicubicL" || method == "BicubicIL"; |
||||
} |
||||
|
||||
bool HasCoast() const |
||||
{ |
||||
if(!vylet) return false; |
||||
vylet->UsePrefix(""); |
||||
return vylet->ParameterBValue("checkcoast", true); |
||||
} |
||||
|
||||
real Left() const |
||||
{ |
||||
vylet->UsePrefix(""); |
||||
return vylet->ParameterRValue("xl", -1.0); |
||||
} |
||||
|
||||
real Right() const |
||||
{ |
||||
vylet->UsePrefix(""); |
||||
return vylet->ParameterRValue("xr", maxx + 1.0); |
||||
} |
||||
|
||||
real Down() const |
||||
{ |
||||
vylet->UsePrefix(""); |
||||
return vylet->ParameterRValue("yd", -1.0); |
||||
} |
||||
|
||||
real Up() const |
||||
{ |
||||
vylet->UsePrefix(""); |
||||
return vylet->ParameterRValue("yu", maxy + 1.0); |
||||
} |
||||
|
||||
bool HasBorders() const { return Left() >= 0.0 || Right() <= maxx || Down() >= 0.0 || Up() <= maxy; } |
||||
|
||||
bool HasTime() const { return HasCoast() || HasBorders(); } |
||||
|
||||
bool LengthFast() const |
||||
{ |
||||
vylet->UsePrefix(""); |
||||
return !vylet->ParameterBValue("sisangle", false); |
||||
} |
||||
|
||||
MString Info() const; |
||||
MString Open(const CLArgs& args); |
||||
}; |
||||
@ -0,0 +1,254 @@
|
||||
#define MICHLIB_NOSOURCE |
||||
#include "cache.h" |
||||
|
||||
SQLiteConnection::DBType SQLiteConnection::db = nullptr; |
||||
size_t SQLiteConnection::count = 0; |
||||
std::vector<SQLiteConnection::FuncType> SQLiteConnection::destructs = {}; |
||||
|
||||
PostgreSQLConnection::DBType PostgreSQLConnection::conn = nullptr; |
||||
size_t PostgreSQLConnection::count = 0; |
||||
std::vector<PostgreSQLConnection::FuncType> PostgreSQLConnection::destructs = {}; |
||||
|
||||
bool SQLiteCache::regdest = false; |
||||
|
||||
bool PostgreSQLCache::regdest = false; |
||||
|
||||
bool FileInfoCache::regdest = false; |
||||
|
||||
void FileInfoCache::GetDirId() |
||||
{ |
||||
if(dirid != 0) return; |
||||
const char* params[] = {dir.Buf()}; |
||||
int plens[] = {int_cast<int>(dir.Len())}; |
||||
int pfor[] = {0}; |
||||
|
||||
PGresultRAII res = PQexecParams(conn, "SELECT id FROM dirs WHERE name=$1::text;", 1, nullptr, params, plens, pfor, 1); |
||||
if(PQresultStatus(res) != PGRES_TUPLES_OK) |
||||
{ |
||||
michlib::errmessage(PQresStatus(PQresultStatus(res))); |
||||
michlib::errmessage(PQerrorMessage(conn)); |
||||
return; |
||||
} |
||||
else if(PQntuples(res) == 0) |
||||
{ |
||||
res = PQexecParams(conn, |
||||
"INSERT INTO dirs(name,id) VALUES ($1, (SELECT min(num.numid) FROM (SELECT generate_series(1, (SELECT COALESCE((SELECT max(id) FROM dirs), 1)) + " |
||||
"1, 1) AS numid) num LEFT JOIN dirs ON dirs.id=num.numid WHERE id IS NULL)) RETURNING id;", |
||||
1, nullptr, params, plens, pfor, 1); |
||||
|
||||
if(PQresultStatus(res) != PGRES_COMMAND_OK && PQresultStatus(res) != PGRES_TUPLES_OK) |
||||
{ |
||||
michlib::errmessage(PQresStatus(PQresultStatus(res))); |
||||
michlib::errmessage(PQerrorMessage(conn)); |
||||
} |
||||
if(PQntuples(res) == 0) return; |
||||
} |
||||
|
||||
if(PQgetlength(res, 0, 0) == sizeof(dirid)) dirid = *pointer_cast<const decltype(dirid)*>(PQgetvalue(res, 0, 0)); |
||||
michlib::message("Dirid: ", Invert(dirid)); |
||||
} |
||||
|
||||
FileInfoCache::FileInfoCache(FileInfoCache::CallbackType&& readfunc_, const MString& dir_): readfunc(std::move(readfunc_)), dir(dir_), dirid(0) |
||||
{ |
||||
if(!conn) return; |
||||
|
||||
if(!regdest) |
||||
{ |
||||
// Create table
|
||||
PGresultRAII res = PQexec(conn, "SET client_min_messages=WARNING;"); |
||||
|
||||
res = PQexec(conn, "BEGIN;" |
||||
"CREATE TABLE IF NOT EXISTS dirs(name TEXT PRIMARY KEY, id INTEGER UNIQUE NOT NULL CONSTRAINT id_is_positive CHECK(id>0));" |
||||
"CREATE TABLE IF NOT EXISTS files(name TEXT NOT NULL, size BIGINT NOT NULL CONSTRAINT size_is_positive CHECK(size>0), modtime TIMESTAMP(0) NOT NULL, " |
||||
"dirid INTEGER REFERENCES dirs(id) ON DELETE CASCADE, lastaccess TIMESTAMP(0) NOT NULL, data BYTEA NOT NULL, PRIMARY KEY(name,dirid));" |
||||
"COMMIT;"); |
||||
if(PQresultStatus(res) != PGRES_COMMAND_OK) |
||||
{ |
||||
michlib::errmessage(PQresStatus(PQresultStatus(res))); |
||||
michlib::errmessage(PQerrorMessage(conn)); |
||||
} |
||||
|
||||
res = PQexec(conn, "SET client_min_messages=NOTICE;"); |
||||
|
||||
conn.AddDestructor( |
||||
[](PostgreSQLConnection::DBType conn) |
||||
{ |
||||
PGresultRAII res = PQexec(conn, "BEGIN;" |
||||
"DELETE FROM files WHERE lastaccess+'100 days'::interval<localtimestamp;" |
||||
"DELETE FROM dirs WHERE id NOT IN (SELECT dirid FROM files);" |
||||
"COMMIT;"); |
||||
}); |
||||
|
||||
regdest = true; |
||||
} |
||||
|
||||
GetDirId(); |
||||
//UpdateCache();
|
||||
} |
||||
|
||||
Error FileInfoCache::UpdateCache(bool force) const |
||||
{ |
||||
const static MString pref = "FileInfoCache::UpdateCache"; |
||||
|
||||
DIRRAII dhandle; |
||||
|
||||
dhandle.reset(opendir(dir.Buf())); |
||||
|
||||
if(!dhandle) return {pref, "Can't open directory " + dir}; |
||||
|
||||
int dfd = dirfd(dhandle); |
||||
errno = 0; |
||||
struct dirent* dent = readdir(dhandle); |
||||
if(errno != 0) return {pref, "Can't read directory " + dir}; |
||||
struct stat st; |
||||
|
||||
do { |
||||
if(dent->d_name[0] != '.') |
||||
{ |
||||
int ret = fstatat(dfd, dent->d_name, &st, 0); |
||||
if(ret != 0) return {pref, "Can't stat " + dir + "/" + dent->d_name}; |
||||
if(S_ISREG(st.st_mode)) // Regular file
|
||||
{ |
||||
const char* params[] = {dent->d_name, pointer_cast<const char*>(&dirid)}; |
||||
int plens[] = {int_cast<int>(strlen(dent->d_name)), sizeof(dirid)}; |
||||
int pfor[] = {0, 1}; |
||||
bool querysucc = true; |
||||
time_t modtime; |
||||
size_t size; |
||||
|
||||
PGresultRAII res = PQexecParams(conn, "SELECT size,modtime FROM files WHERE name=$1::text AND dirid=$2::integer;", 2, nullptr, params, plens, pfor, 1); |
||||
|
||||
if(PQresultStatus(res) != PGRES_COMMAND_OK && PQresultStatus(res) != PGRES_TUPLES_OK) |
||||
{ |
||||
michlib::errmessage(PQresStatus(PQresultStatus(res))); |
||||
michlib::errmessage(PQerrorMessage(conn)); |
||||
querysucc = false; |
||||
} |
||||
else if(PQntuples(res) == 0 || PQntuples(res) > 1) |
||||
querysucc = false; |
||||
|
||||
if(querysucc) |
||||
{ |
||||
size = *pointer_cast<const decltype(size)*>(PQgetvalue(res, 0, 0)); |
||||
modtime = raw2epoch(*pointer_cast<const time_t*>(PQgetvalue(res, 0, 1))); |
||||
} |
||||
else |
||||
{ |
||||
size = int_cast<size_t>(st.st_size); |
||||
modtime = st.st_mtim.tv_sec; |
||||
} |
||||
|
||||
if(!querysucc || force || size != int_cast<size_t>(st.st_size) || modtime != st.st_mtim.tv_sec) |
||||
{ |
||||
auto ret = GetData(dent->d_name); |
||||
// Remove entry
|
||||
if(!ret && querysucc) |
||||
{ |
||||
PGresultRAII dres = PQexecParams(conn, "DELETE FROM files WHERE name=$1::text AND dirid=$2::integer;", 2, nullptr, params, plens, pfor, 1); |
||||
if(PQresultStatus(dres) != PGRES_COMMAND_OK) |
||||
{ |
||||
michlib::errmessage(PQresStatus(PQresultStatus(dres))); |
||||
michlib::errmessage(PQerrorMessage(conn)); |
||||
} |
||||
} |
||||
else // Update or insert
|
||||
{ |
||||
auto sizei = Invert(size); |
||||
auto modtimei = epoch2raw(modtime); |
||||
|
||||
const char* params[] = {dent->d_name, pointer_cast<const char*>(&sizei), pointer_cast<const char*>(&modtimei), pointer_cast<const char*>(&dirid), ret.value().Buf()}; |
||||
int plens[] = {int_cast<int>(strlen(dent->d_name)), sizeof(sizei), sizeof(modtimei), sizeof(dirid), int_cast<int>(ret.value().Len())}; |
||||
int pfor[] = {0, 1, 1, 1, 1}; |
||||
|
||||
PGresultRAII res = PQexecParams(conn, |
||||
"INSERT INTO files (name,size,modtime,dirid,lastaccess,data) VALUES($1::text, $2::bigint, $3::timestamp, $4::integer, localtimestamp, $5) " |
||||
"ON CONFLICT ON CONSTRAINT files_pkey DO UPDATE SET " |
||||
"size=EXCLUDED.size, modtime=EXCLUDED.modtime, lastaccess=EXCLUDED.lastaccess, data=EXCLUDED.data;", |
||||
5, nullptr, params, plens, pfor, 1); |
||||
if(PQresultStatus(res) != PGRES_COMMAND_OK) |
||||
{ |
||||
michlib::errmessage(PQresStatus(PQresultStatus(res))); |
||||
michlib::errmessage(PQerrorMessage(conn)); |
||||
} |
||||
} // Insert or update branch
|
||||
} // Need data update
|
||||
} // Regular file
|
||||
} // if(dent->d_name[0] != '.')
|
||||
dent = readdir(dhandle); |
||||
} while(dent != nullptr || errno != 0); |
||||
|
||||
return Error(); |
||||
} |
||||
|
||||
FileInfoCache::DataType FileInfoCache::GetInfo(const MString& name) const |
||||
{ |
||||
if(!*this) return GetData(name); |
||||
|
||||
bool querysucc = true; |
||||
MString data; |
||||
time_t modtime; |
||||
size_t size; |
||||
|
||||
{ |
||||
const char* params[] = {name.Buf(), pointer_cast<const char*>(&dirid)}; |
||||
int plens[] = {int_cast<int>(name.Len()), sizeof(dirid)}; |
||||
int pfor[] = {0, 1}; |
||||
|
||||
PGresultRAII res = |
||||
PQexecParams(conn, "UPDATE files SET lastaccess=localtimestamp WHERE name=$1::text AND dirid=$2::integer RETURNING data,size,modtime;", 2, nullptr, params, plens, pfor, 1); |
||||
if(PQresultStatus(res) != PGRES_COMMAND_OK && PQresultStatus(res) != PGRES_TUPLES_OK) |
||||
{ |
||||
michlib::errmessage(PQresStatus(PQresultStatus(res))); |
||||
michlib::errmessage(PQerrorMessage(conn)); |
||||
querysucc = false; |
||||
} |
||||
|
||||
if(PQntuples(res) == 0 || PQntuples(res) > 1) |
||||
{ |
||||
michlib::errmessage("Data for file ", dir + "/" + name, (PQntuples(res) == 0 ? " not found " : " duplicated "), "in cache"); |
||||
querysucc = false; |
||||
} |
||||
|
||||
if(querysucc) |
||||
{ |
||||
data = MString(PQgetvalue(res, 0, 0), PQgetlength(res, 0, 0)); |
||||
size = *pointer_cast<const decltype(size)*>(PQgetvalue(res, 0, 1)); |
||||
modtime = raw2epoch(*pointer_cast<const time_t*>(PQgetvalue(res, 0, 2))); |
||||
} |
||||
} |
||||
|
||||
{ |
||||
struct stat st; |
||||
|
||||
int ret = stat((dir + "/" + name).Buf(), &st); |
||||
if(ret != 0) return DataType(); |
||||
|
||||
if(querysucc && st.st_mtim.tv_sec == modtime && size == int_cast<size_t>(st.st_size)) return data; |
||||
modtime = st.st_mtim.tv_sec; |
||||
size = st.st_size; |
||||
} |
||||
|
||||
auto ret = GetData(name); |
||||
if(ret) |
||||
{ |
||||
auto sizei = Invert(size); |
||||
auto modtimei = epoch2raw(modtime); |
||||
|
||||
const char* params[] = {name.Buf(), pointer_cast<const char*>(&sizei), pointer_cast<const char*>(&modtimei), pointer_cast<const char*>(&dirid), ret.value().Buf()}; |
||||
int plens[] = {int_cast<int>(name.Len()), sizeof(sizei), sizeof(modtimei), sizeof(dirid), int_cast<int>(ret.value().Len())}; |
||||
int pfor[] = {0, 1, 1, 1, 1}; |
||||
|
||||
PGresultRAII res = PQexecParams(conn, |
||||
"INSERT INTO files (name,size,modtime,dirid,lastaccess,data) VALUES($1::text, $2::bigint, $3::timestamp, $4::integer, localtimestamp, $5) " |
||||
"ON CONFLICT ON CONSTRAINT files_pkey DO UPDATE SET " |
||||
"size=EXCLUDED.size, modtime=EXCLUDED.modtime, lastaccess=EXCLUDED.lastaccess, data=EXCLUDED.data;", |
||||
5, nullptr, params, plens, pfor, 1); |
||||
if(PQresultStatus(res) != PGRES_COMMAND_OK) |
||||
{ |
||||
michlib::errmessage(PQresStatus(PQresultStatus(res))); |
||||
michlib::errmessage(PQerrorMessage(conn)); |
||||
} |
||||
} |
||||
|
||||
return ret; |
||||
} |
||||
@ -0,0 +1,180 @@
|
||||
#define MICHLIB_NOSOURCE |
||||
#include "copcat.h" |
||||
#include "GPL.h" |
||||
#include "mirrorfuncs.h" |
||||
|
||||
const MString CopernicusCatalog::caturl = "https://stac.marine.copernicus.eu/metadata/catalog.stac.json"; |
||||
|
||||
CopernicusCatalog::CopernicusCatalog() |
||||
{ |
||||
auto oldprefix = michlib::GPL.UsePrefix("COPERNICUS"); |
||||
// Cache
|
||||
cache.reset(CreateCache(michlib::GPL.ParameterSValue("Cache", ""))); |
||||
if(!cache) |
||||
{ |
||||
michlib::errmessage("Can't init cache"); |
||||
cache.reset(new FakeCache); |
||||
} |
||||
|
||||
// Proxy
|
||||
auto proxyurl = michlib::GPL.ParameterSValue("Proxy", ""); |
||||
if(proxyurl.Exist()) curl_easy_setopt(chandle, CURLOPT_PROXY, proxyurl.Buf()); |
||||
|
||||
michlib::GPL.UsePrefix(oldprefix); |
||||
GetCatalog(); |
||||
} |
||||
|
||||
Error CopernicusCatalog::GetCatalog() |
||||
{ |
||||
if(Valid()) return Error(); |
||||
auto ret = GetJSON(caturl); |
||||
if(ret) |
||||
catalog = ret.Value(); |
||||
else |
||||
return ret.Add("CopernicusCatalog::GetCatalog", "can't download catalog"); |
||||
return Error(); |
||||
} |
||||
|
||||
RetVal<std::vector<MString>> CopernicusCatalog::ProductList() const |
||||
{ |
||||
static const MString pref = "CopernicusCatalog::ProductList"; |
||||
|
||||
if(!Valid()) return {pref, "no catalog"}; |
||||
|
||||
const auto& links = catalog["links"]; |
||||
if(links.type() != Json::arrayValue) return {pref, "no \"links\" section in the catalog"}; |
||||
|
||||
std::vector<MString> out; |
||||
for(Json::ArrayIndex i = 0; i < links.size(); i++) |
||||
{ |
||||
const auto& rel = links[i]["rel"]; |
||||
const auto& href = links[i]["href"]; |
||||
|
||||
if(rel.type() == Json::stringValue && href.type() == Json::stringValue && rel.asString() == "child") |
||||
{ |
||||
auto str = href.asString(); |
||||
str.erase(str.find('/')); |
||||
out.emplace_back(str.c_str()); |
||||
} |
||||
} |
||||
return out; |
||||
} |
||||
|
||||
RetVal<MString> CopernicusCatalog::ProductURL(const MString& prod) const |
||||
{ |
||||
static const MString pref = "CopernicusCatalog::ProductURL"; |
||||
|
||||
if(!Valid()) return {pref, "no catalog"}; |
||||
|
||||
const auto& links = catalog["links"]; |
||||
if(links.type() != Json::arrayValue) return {pref, "no \"links\" section in the catalog"}; |
||||
|
||||
for(Json::ArrayIndex i = 0; i < links.size(); i++) |
||||
{ |
||||
const auto& href = links[i]["href"]; |
||||
if(href.type() == Json::stringValue && href.asString() == (prod + "/product.stac.json").Buf()) return DirName(caturl) + "/" + MString(href.asString().c_str()); |
||||
} |
||||
return {pref, "unknown product: " + prod}; |
||||
} |
||||
|
||||
RetVal<std::vector<MString>> CopernicusCatalog::DatasetList(const MString& prod) const |
||||
{ |
||||
static const MString pref = "CopernicusCatalog::DatasetList"; |
||||
|
||||
MString url; |
||||
{ |
||||
auto ret = ProductURL(prod); |
||||
if(!ret) return ret.Add(pref, "Can't get url for the product " + prod); |
||||
url = ret.Value(); |
||||
} |
||||
|
||||
auto ret = GetJSON(url); |
||||
if(!ret) return ret.Add(pref, "Can't download product " + prod); |
||||
|
||||
const auto& links = ret.Value()["links"]; |
||||
if(links.type() != Json::arrayValue) return {pref, "no \"links\" section in the product " + prod + " description"}; |
||||
|
||||
std::vector<MString> out; |
||||
for(Json::ArrayIndex i = 0; i < links.size(); i++) |
||||
{ |
||||
const auto& rel = links[i]["rel"]; |
||||
const auto& href = links[i]["href"]; |
||||
|
||||
if(rel.type() == Json::stringValue && href.type() == Json::stringValue && rel.asString() == "item") |
||||
{ |
||||
auto str = href.asString(); |
||||
str.erase(str.find('/')); |
||||
out.emplace_back(str.c_str()); |
||||
} |
||||
} |
||||
return out; |
||||
} |
||||
|
||||
RetVal<MString> CopernicusCatalog::DatasetURL(const MString& prod, const MString& dataset) const |
||||
{ |
||||
static const MString pref = "CopernicusCatalog::DatasetURL"; |
||||
|
||||
MString url; |
||||
{ |
||||
auto ret = ProductURL(prod); |
||||
if(!ret) return ret.Add(pref, "Can't get url for the product " + prod); |
||||
url = ret.Value(); |
||||
} |
||||
|
||||
auto ret = GetJSON(url); |
||||
if(!ret) return ret.Add(pref, "Can't download product " + prod); |
||||
|
||||
const auto& links = ret.Value()["links"]; |
||||
if(links.type() != Json::arrayValue) return {pref, "no \"links\" section in the product " + prod + " description"}; |
||||
|
||||
for(Json::ArrayIndex i = 0; i < links.size(); i++) |
||||
{ |
||||
const auto& href = links[i]["href"]; |
||||
if(href.type() == Json::stringValue && href.asString() == (dataset + "/dataset.stac.json").Buf()) return DirName(url) + "/" + MString(href.asString().c_str()); |
||||
} |
||||
return {pref, "unknown dataset: " + dataset}; |
||||
} |
||||
|
||||
RetVal<MString> CopernicusCatalog::AssetURL(const MString& prod, const MString& dataset, const MString& asset) const |
||||
{ |
||||
static const MString pref = "CopernicusCatalog::AssetURL"; |
||||
|
||||
MString url; |
||||
{ |
||||
auto ret = DatasetURL(prod, dataset); |
||||
if(!ret) return ret.Add(pref, "Can't get url for the dataset " + dataset); |
||||
url = ret.Value(); |
||||
} |
||||
|
||||
auto ret = GetJSON(url); |
||||
if(!ret) return ret.Add(pref, "Can't download dataset " + dataset); |
||||
|
||||
const auto& href = ret.Value()["assets"][asset.Buf()]["href"]; |
||||
if(!href || href.type() != Json::stringValue) return {pref, "href for the asset " + asset + " not found"}; |
||||
return MString(href.asString().c_str()); |
||||
} |
||||
|
||||
RetVal<Json::Value> CopernicusCatalog::GetJSON(const MString& url) const |
||||
{ |
||||
const static MString pref = "CopernicusCatalog::GetJSON"; |
||||
|
||||
Json::Reader reader; |
||||
Json::Value obj; |
||||
MString content; |
||||
|
||||
auto [val, suc] = cache->Get(url); |
||||
if(suc) |
||||
content = std::move(val); |
||||
else |
||||
{ |
||||
michlib::message(url + " not found in cache, downloading"); |
||||
auto [out, res] = GetUrl(chandle, url); |
||||
if(res != CURLE_OK) return Error(pref, MString("can't download JSON: ") + chandle.Err()); |
||||
cache->Put(url, out, 3600); |
||||
content = std::move(out); |
||||
} |
||||
|
||||
reader.parse(content.Buf(), content.Buf() + content.Len(), obj, false); |
||||
|
||||
return obj; |
||||
} |
||||
@ -0,0 +1,41 @@
|
||||
#define MICHLIB_NOSOURCE |
||||
#include "curlfuncs.h" |
||||
#include <unistd.h> |
||||
|
||||
using michlib::pointer_cast; |
||||
using michlib::uint1; |
||||
|
||||
size_t Write2String(char* ptr, size_t size, size_t n, void* data) |
||||
{ |
||||
MString* out = pointer_cast<MString*>(data); |
||||
*out += MString(ptr, size * n); |
||||
return size * n; |
||||
} |
||||
|
||||
size_t Write2File(char* ptr, size_t size, size_t n, void* data) |
||||
{ |
||||
const int* fd = pointer_cast<const int*>(data); |
||||
size_t count = size * n; |
||||
const uint1* buf = pointer_cast<const uint1*>(ptr); |
||||
|
||||
while(count != 0) |
||||
{ |
||||
auto wr = write(*fd, buf, count); |
||||
if(wr == -1) return 0; |
||||
count -= wr; |
||||
buf += wr; |
||||
} |
||||
|
||||
return size * n; |
||||
} |
||||
|
||||
std::pair<MString, CURLcode> GetUrl(const CURLRAII& chandle, const MString& url) |
||||
{ |
||||
MString out; |
||||
|
||||
curl_easy_setopt(chandle, CURLOPT_URL, url.Buf()); |
||||
curl_easy_setopt(chandle, CURLOPT_WRITEFUNCTION, Write2String); |
||||
curl_easy_setopt(chandle, CURLOPT_WRITEDATA, &out); |
||||
auto res = curl_easy_perform(chandle); |
||||
return {out, res}; |
||||
} |
||||
@ -0,0 +1,295 @@
|
||||
#define MICHLIB_NOSOURCE |
||||
#include "layereddataz.h" |
||||
|
||||
MString LayeredDataZ::Info() const |
||||
{ |
||||
if(!isOk()) return ""; |
||||
MString d; |
||||
for(size_t i = 0; i < NDepths(); i++) d += MString(" ") + "(" + i + " " + Depth(i) + ")"; |
||||
|
||||
std::set<MString> vars; |
||||
for(const auto& f: nc) GetVars(f, vars); |
||||
|
||||
MString svars; |
||||
{ |
||||
bool first = true; |
||||
for(const auto& v: vars) |
||||
{ |
||||
svars += (first ? "" : ", ") + v; |
||||
first = false; |
||||
} |
||||
} |
||||
|
||||
// clang-format off
|
||||
return |
||||
"Dataset: " + Title() + "\n" + |
||||
" Begin date: " + Time(0).ToString() + "\n" + |
||||
" End date: " + Time(NTimes()-1).ToString() + "\n" + |
||||
" Time step: " + Timestep() + " seconds\n" + |
||||
" Time moments: " + NTimes() + "\n" + |
||||
" Region: (" + lonb + " : " + lone + ") x (" + latb + " : " + late + ")\n" + |
||||
" Grid: " + dname.nx + "x" + dname.ny + " (" + lonstep + " x " + latstep + ")\n" + |
||||
" Depths:" + d + "\n" + |
||||
" Supported variables: " + svars; |
||||
// clang-format on
|
||||
} |
||||
|
||||
MString LayeredDataZ::Open(const MString& dataset) |
||||
{ |
||||
nc.clear(); |
||||
MString proxyurl = GPL.ParameterSValue("USEPROXY", ""); |
||||
if(proxyurl.Exist()) proxy.Activate("all_proxy", proxyurl); |
||||
|
||||
nc.clear(); |
||||
size_t i = 1; |
||||
while(true) |
||||
{ |
||||
MString url = GPL.ParameterSValue(dataset + "_URL" + i, ""); |
||||
if(url.Exist()) |
||||
{ |
||||
// Split url on product and dataset
|
||||
auto words = url.Split(":"); |
||||
if(words.size() == 0 || words.size() > 2) |
||||
{ |
||||
nc.clear(); |
||||
return "Invalid url " + url; |
||||
} |
||||
MString product = words[0]; |
||||
MString dataset = words.size() == 2 ? words[1] : ""; |
||||
|
||||
nc.emplace_back(); |
||||
{ |
||||
auto ret = nc.back().OpenZarr(product, dataset); |
||||
if(!ret) |
||||
{ |
||||
nc.clear(); |
||||
return "Can't open " + dataset + " of " + product; |
||||
} |
||||
} |
||||
} |
||||
else |
||||
break; |
||||
i++; |
||||
} |
||||
if(nc.size() == 0) return "No urls for dataset " + dataset + " specified in config"; |
||||
|
||||
dname = GetDNames(nc[0]); |
||||
if(!(dname.lonname.Exist() && dname.latname.Exist())) |
||||
{ |
||||
nc.clear(); |
||||
return "Can't find longitude/latitude"; |
||||
} |
||||
if(!dname.timename.Exist()) |
||||
{ |
||||
nc.clear(); |
||||
return "Can't find time"; |
||||
} |
||||
|
||||
auto cn = GetCNames(nc[0]); |
||||
|
||||
// Read times
|
||||
for(auto& f: nc) |
||||
{ |
||||
auto ret = f.ReadTimes(cn.timename); |
||||
if(!ret) |
||||
{ |
||||
nc.clear(); |
||||
return "Can't read times"; |
||||
} |
||||
times.insert(times.end(), f.Times().begin(), f.Times().end()); |
||||
} |
||||
std::sort(times.begin(), times.end()); |
||||
auto last = std::unique(times.begin(), times.end()); |
||||
times.erase(last, times.end()); |
||||
|
||||
depthinv = false; |
||||
if(cn.depthname.Exist()) |
||||
{ |
||||
auto ret = nc[0].Read(cn.depthname, depths); |
||||
if(!ret) |
||||
{ |
||||
nc.clear(); |
||||
return "Can't read depths"; |
||||
} |
||||
if(depths.back() <= 0 && depths.front() <= 0) std::ranges::transform(depths, depths.begin(), std::negate{}); |
||||
if(depths.back() < depths.front() && depths.size() > 1) |
||||
{ |
||||
depthinv = true; |
||||
for(size_t i = 0; i < depths.size() - i - 1; i++) std::swap(depths[i], depths[depths.size() - i - 1]); |
||||
} |
||||
} |
||||
else // Surface only data
|
||||
{ |
||||
depths.resize(1); |
||||
depths[0] = 0; |
||||
} |
||||
|
||||
std::vector<double> lons, lats; |
||||
{ |
||||
auto ret = nc[0].Read(cn.lonname, lons); |
||||
if(!ret) |
||||
{ |
||||
nc.clear(); |
||||
return "Can't get longitudes"; |
||||
} |
||||
} |
||||
{ |
||||
auto ret = nc[0].Read(cn.latname, lats); |
||||
if(!ret) |
||||
{ |
||||
nc.clear(); |
||||
return "Can't get latitudes"; |
||||
} |
||||
} |
||||
|
||||
lonb = lons[0]; |
||||
latb = lats[0]; |
||||
lone = lons.back(); |
||||
late = lats.back(); |
||||
lonstep = (lone - lonb) / (dname.nx - 1); |
||||
latstep = (late - latb) / (dname.ny - 1); |
||||
|
||||
return ""; |
||||
} |
||||
|
||||
std::pair<const BaseParameters*, MString> LayeredDataZ::Parameters(michlib_internal::ParameterListEx& pars, const CLArgs& args, const struct Region& reg) const |
||||
{ |
||||
std::unique_ptr<struct Parameters> ppar{new struct Parameters}; |
||||
|
||||
ppar->layer = args.contains("layer") ? args.at("layer").ToInteger<size_t>() : 0; |
||||
if(!args.contains("depth") && ppar->layer >= NDepths()) return {nullptr, MString("Layer ") + ppar->layer + " is too deep!"}; |
||||
real depth = args.contains("depth") ? args.at("depth").ToReal() : Depth(ppar->layer); |
||||
|
||||
{ |
||||
auto dom = DetGeoDomain(lonb, lone); |
||||
real lon1 = ToGeoDomain(reg.lonb, dom); |
||||
real lon2 = ToGeoDomain(reg.lone, dom); |
||||
real lat1 = reg.latb; |
||||
real lat2 = reg.late; |
||||
bool global = lone - lonb + 1.5 * lonstep > 360.0; |
||||
|
||||
// Special case when the longitude lies in a small sector between the end and the start
|
||||
if(global) |
||||
{ |
||||
if(lon1 < lonb) lon1 = lone; |
||||
if(lon2 > lone) lon2 = lonb; |
||||
} |
||||
else |
||||
{ |
||||
if(lon1 < lonb) lon1 = lonb; |
||||
if(lon2 > lone) lon2 = lone; |
||||
} |
||||
|
||||
ppar->xb = static_cast<size_t>(Floor((lon1 - lonb) / lonstep)); |
||||
ppar->xe = static_cast<size_t>(Ceil((lon2 - lonb) / lonstep)); |
||||
|
||||
if(ppar->xb == ppar->xe) return {nullptr, "Lonb must be not equal late"}; |
||||
if(!global && ppar->xb > ppar->xe) return {nullptr, "Lonb must be lesser then lone"}; |
||||
|
||||
ppar->yb = static_cast<size_t>(Floor((lat1 - latb) / latstep)); |
||||
ppar->ye = static_cast<size_t>(Ceil((lat2 - latb) / latstep)); |
||||
if(ppar->ye > dname.ny - 1) ppar->ye = dname.ny - 1; |
||||
if(ppar->yb >= ppar->ye) return {nullptr, "Latb must be lesser then late"}; |
||||
|
||||
if(depth < 0.0 || depth > depths.back()) |
||||
ppar->layer = (depth < 0.0) ? 0 : (depths.size() - 1); |
||||
else |
||||
for(size_t i = 0; i < depths.size() - 1; i++) |
||||
{ |
||||
if(depth >= depths[i] && depth <= depths[i + 1]) |
||||
{ |
||||
ppar->layer = (depth - depths[i] <= depths[i + 1] - depth) ? i : (i + 1); |
||||
break; |
||||
} |
||||
} |
||||
} |
||||
|
||||
pars.SetParameter("depth", Depth(ppar->layer)); |
||||
pars.SetParameter("layer", ppar->layer); |
||||
pars.SetParameter("dataset", Title()); |
||||
pars.SetParameter("lonb", Lon(ppar->xb)); |
||||
pars.SetParameter("latb", Lat(ppar->yb)); |
||||
pars.SetParameter("lone", Lon(ppar->xe)); |
||||
pars.SetParameter("late", Lat(ppar->ye)); |
||||
|
||||
return {ppar.release(), ""}; |
||||
} |
||||
|
||||
bool LayeredDataZ::Read(const MString& vname, std::map<MString, LayeredDataZ::Data>& cache, const BaseParameters* ip, size_t i) const |
||||
{ |
||||
if(cache.contains(vname)) return true; |
||||
if(!isOk()) return false; |
||||
|
||||
auto p = dynamic_cast<const struct Parameters*>(ip); |
||||
auto [name, id, tid] = VarNameLoc(vname, times[i]); |
||||
if(!name.Exist()) // Conversion read
|
||||
return TransformationRead(this, vname, cache, ip, i); |
||||
|
||||
// Direct read
|
||||
bool nodepth = false; |
||||
Data data; |
||||
//auto head = nc[id]->Header();
|
||||
for(const auto& v: nc[id].VarNames()) |
||||
if(v == name) |
||||
{ |
||||
if(nc[id].NDim(v) == 3) nodepth = true; |
||||
data = ReadVarRaw(nc[id], name, tid, nodepth, p); |
||||
if(data) |
||||
{ |
||||
cache[vname] = std::move(data); |
||||
return true; |
||||
} |
||||
} |
||||
|
||||
return false; |
||||
} |
||||
|
||||
LayeredDataZ::Data LayeredDataZ::ReadVarRaw(const NC& f, const MString& name, size_t i, bool nodepth, const struct LayeredDataZ::Parameters* p) const |
||||
{ |
||||
real unitmul = 1.0; |
||||
real offset = 0.0, scale = 1.0; |
||||
|
||||
if(f.HasAtt(name, "add_offset")) offset = f.AttReal(name, "add_offset"); |
||||
if(f.HasAtt(name, "scale_factor")) scale = f.AttReal(name, "scale_factor"); |
||||
|
||||
MString unit; |
||||
if(f.HasAtt(name, "units")) unit = f.AttString(name, "units"); |
||||
if(unit == "m s-1" || unit == "m/s") |
||||
{ |
||||
unitmul = 100.0; |
||||
unit = "cm/s"; |
||||
} |
||||
|
||||
Data data((p->xb < p->xe) ? (p->xe - p->xb + 1) : (dname.nx + p->xe - p->xb + 1), p->ye - p->yb + 1, Lon(p->xb), Lat(p->yb), lonstep, latstep, std::move(unit)); |
||||
|
||||
auto trans = [scale, offset, unitmul](auto raw) -> real { return (raw * scale + offset) * unitmul; }; |
||||
|
||||
auto rlayer = depthinv ? depths.size() - p->layer - 1 : p->layer; |
||||
|
||||
if(p->xb < p->xe) |
||||
{ |
||||
auto ret = nodepth ? f.Read(name, data, trans, {dname.lonname, p->xb, p->xe - p->xb + 1}, {dname.latname, p->yb, p->ye - p->yb + 1}, {dname.timename, i, 1}) |
||||
: f.Read(name, data, trans, {dname.lonname, p->xb, p->xe - p->xb + 1}, {dname.latname, p->yb, p->ye - p->yb + 1}, {dname.timename, i, 1}, |
||||
{dname.depthname, rlayer, 1}); |
||||
if(!ret) return Data(); |
||||
} |
||||
else |
||||
{ |
||||
{ |
||||
auto ret = nodepth ? f.Read(name, data, trans, {dname.lonname, p->xb, dname.nx - p->xb}, {dname.latname, p->yb, p->ye - p->yb + 1}, {dname.timename, i, 1}) |
||||
: f.Read(name, data, trans, {dname.lonname, p->xb, dname.nx - p->xb}, {dname.latname, p->yb, p->ye - p->yb + 1}, {dname.timename, i, 1}, |
||||
{dname.depthname, rlayer, 1}); |
||||
if(!ret) return Data(); |
||||
} |
||||
{ |
||||
size_t shift = dname.nx - p->xb + 1; |
||||
auto shifteddata = [&data, shift](size_t ix, size_t iy) -> real& { return data(ix + shift, iy); }; |
||||
|
||||
auto ret = |
||||
nodepth ? f.Read(name, shifteddata, trans, {dname.lonname, 0, p->xe + 1}, {dname.latname, p->yb, p->ye - p->yb + 1}, {dname.timename, i, 1}) |
||||
: f.Read(name, shifteddata, trans, {dname.lonname, 0, p->xe + 1}, {dname.latname, p->yb, p->ye - p->yb + 1}, {dname.timename, i, 1}, {dname.depthname, rlayer, 1}); |
||||
if(!ret) return Data(); |
||||
} |
||||
} |
||||
return data; |
||||
} |
||||
@ -0,0 +1,134 @@
|
||||
#define MICHLIB_NOSOURCE |
||||
#include "mirrorfuncs.h" |
||||
#include "StringFunctions.h" |
||||
#include "filehelpers.h" |
||||
#include "merrors.h" |
||||
|
||||
using michlib::FD; |
||||
using michlib::message; |
||||
|
||||
bool MakePath(const MString& dname) |
||||
{ |
||||
struct stat st; |
||||
int ret = stat(dname.Buf(), &st); |
||||
if(ret == 0) return S_ISDIR(st.st_mode); |
||||
|
||||
auto dirs = michlib::Split_on_words(dname, "/", false); |
||||
MString cdir = ""; |
||||
for(const auto& dir: dirs) |
||||
{ |
||||
cdir += "/" + dir; |
||||
ret = stat(cdir.Buf(), &st); |
||||
if(ret == 0 && S_ISDIR(st.st_mode)) continue; |
||||
if(ret == 0 && !S_ISDIR(st.st_mode)) return false; |
||||
ret = mkdir(cdir.Buf(), 0755); |
||||
if(ret != 0) return false; |
||||
} |
||||
return true; |
||||
} |
||||
|
||||
RetVal<std::vector<struct FileInfo>> ReadLocalFileList(const MString& dir, const bool nofollow, const MString& path) |
||||
{ |
||||
const static MString pref = "ReadLocalFileList"; |
||||
|
||||
std::vector<struct FileInfo> out; |
||||
DIRRAII dhandle; |
||||
|
||||
MakePath(dir); |
||||
dhandle.reset(opendir(dir.Buf())); |
||||
|
||||
if(!dhandle) return {pref, "Can't open directory " + path + (path.Exist() ? "/" : "") + dir}; |
||||
|
||||
int dfd = dirfd(dhandle); |
||||
errno = 0; |
||||
struct dirent* dent = readdir(dhandle); |
||||
if(errno != 0) return {pref, "Can't read directory " + path + (path.Exist() ? "/" : "") + dir}; |
||||
struct stat st; |
||||
|
||||
do { |
||||
if(dent->d_name[0] != '.') |
||||
{ |
||||
int ret = fstatat(dfd, dent->d_name, &st, nofollow ? AT_SYMLINK_NOFOLLOW : 0); |
||||
if(ret != 0) return {pref, "Can't stat " + path + "/" + dir + "/" + dent->d_name}; |
||||
if(S_ISDIR(st.st_mode)) // Directory, recurse
|
||||
{ |
||||
auto list = ReadLocalFileList(dir + "/" + dent->d_name, nofollow, path + (path.Exist() ? "/" : "") + dent->d_name); |
||||
if(!list) return list; |
||||
out.insert(out.end(), list.Value().begin(), list.Value().end()); |
||||
} |
||||
if(S_ISREG(st.st_mode)) // Regular file
|
||||
{ |
||||
out.emplace_back(dir + "/" + dent->d_name, path + (path.Exist() ? "/" : "") + dent->d_name, MDateTime(st.st_mtim.tv_sec, st.st_mtim.tv_nsec), st.st_size); |
||||
} |
||||
// Ignore non-directories, non-files
|
||||
} |
||||
dent = readdir(dhandle); |
||||
} while(dent != nullptr || errno != 0); |
||||
|
||||
if(errno != 0) return {pref, "Can't read directory " + path + "/" + dir}; |
||||
std::sort(out.begin(), out.end(), [](const struct FileInfo& a, const struct FileInfo& b) { return a.name < b.name; }); |
||||
return out; |
||||
} |
||||
|
||||
Error DownloadFile(const CURLRAII& chandle, const struct FileInfo& rinfo, const MString& root) |
||||
{ |
||||
const static MString pref = "DownloadFile"; |
||||
|
||||
message("Downloading " + rinfo.url); |
||||
|
||||
MString dname = DirName(rinfo.name), fname = FileName(rinfo.name); |
||||
FD fd; |
||||
|
||||
if(!MakePath(root + "/" + dname)) return {pref, "Can't create directory " + root + "/" + dname}; |
||||
fd.Reset(creat((root + "/" + rinfo.name).Buf(), 0644)); |
||||
if(!fd) return {pref, "Can't create file " + root + "/" + rinfo.name}; |
||||
|
||||
int cfd = fd.Get(); |
||||
curl_easy_setopt(chandle, CURLOPT_WRITEFUNCTION, Write2File); |
||||
curl_easy_setopt(chandle, CURLOPT_WRITEDATA, &cfd); |
||||
curl_easy_setopt(chandle, CURLOPT_URL, rinfo.url.Buf()); |
||||
auto res = curl_easy_perform(chandle); |
||||
if(res != CURLE_OK) |
||||
{ |
||||
unlink((root + "/" + rinfo.name).Buf()); |
||||
return {pref, MString("Can't download file: ") + chandle.Err()}; |
||||
} |
||||
|
||||
{ |
||||
struct timespec times[2]; |
||||
times[0].tv_sec = times[1].tv_sec = rinfo.mtime.Epoch(); |
||||
times[0].tv_nsec = times[1].tv_nsec = 0; |
||||
|
||||
int ret = futimens(fd, times); |
||||
if(ret != 0) |
||||
{ |
||||
unlink((root + "/" + rinfo.name).Buf()); |
||||
return {pref, "Can't set mtime for file: " + root + "/" + rinfo.name}; |
||||
} |
||||
} |
||||
|
||||
return Error(); |
||||
} |
||||
|
||||
Error RemoveFile(const struct FileInfo& linfo) |
||||
{ |
||||
const static MString pref = "RemoveFile"; |
||||
|
||||
message("Remove " + linfo.url); |
||||
int ret = unlink(linfo.url.Buf()); |
||||
if(ret != 0) return {pref, "Can't remove file " + linfo.url}; |
||||
return Error(); |
||||
} |
||||
|
||||
Error UpdateFile(const CURLRAII& chandle, const struct FileInfo& rinfo, const struct FileInfo& linfo, const MString& root) |
||||
{ |
||||
const static MString pref = "UpdateFile"; |
||||
|
||||
message("Update " + linfo.url); |
||||
auto rm = RemoveFile(linfo); |
||||
if(!rm) return rm.Add(pref, "Can't remove file"); |
||||
auto df = DownloadFile(chandle, rinfo, root); |
||||
if(!df) return df.Add(pref, "Can't download file"); |
||||
|
||||
return Error(); |
||||
} |
||||
@ -1,105 +1,143 @@
|
||||
#define MICHLIB_NOSOURCE |
||||
#include "ncfilew.h" |
||||
|
||||
MString NCFileW::CreateFile(NCFileW::Type stype, const MString& name, const MString& history, int compression, size_t nx, size_t ny) |
||||
MString NCFileW::CreateFile(NCFileW::Type stype, const MString& name, int compression, size_t nx, size_t ny) |
||||
{ |
||||
if(stype == UNKNOWN) return "Can't determine file type"; |
||||
|
||||
compress = compression; |
||||
|
||||
int ret; |
||||
MString text; |
||||
const float node_offset = 0.0; |
||||
|
||||
ret = nc_create(name.Buf(), NC_CLOBBER | NC_NETCDF4, &ncid); |
||||
if(ret != NC_NOERR) return "Can't create netcdf file: " + name; |
||||
Open(name); |
||||
if(!*this) return "Can't create netcdf file " + name + ": " + ErrMessage(); |
||||
|
||||
ret = nc_put_att_text(ncid, NC_GLOBAL, "history", history.Len() + 1, history.Buf()); |
||||
if(ret != NC_NOERR) return "Can't write history attribute in the netcdf file"; |
||||
AddAtt("node_offset", node_offset); |
||||
|
||||
switch(stype) |
||||
{ |
||||
case(G1V1): |
||||
case(PSET): |
||||
{ |
||||
ret = nc_def_dim(ncid, "i", nx, &xdimid); |
||||
if(ret != NC_NOERR) return "Can't create dimension in the netcdf file"; |
||||
ret = nc_def_var(ncid, "longitude", NC_FLOAT, 1, &xdimid, &xid); |
||||
if(ret != NC_NOERR) return "Can't create longitude variable in the netcdf file"; |
||||
ret = nc_def_var(ncid, "latitude", NC_FLOAT, 1, &xdimid, &yid); |
||||
if(ret != NC_NOERR) return "Can't create latitude variable in the netcdf file"; |
||||
AddDim("i", nx); |
||||
AddVar("x", NC_FLOAT, "i"); |
||||
AddVar("y", NC_FLOAT, "i"); |
||||
break; |
||||
} |
||||
case(G1V2): |
||||
case(GPSET): |
||||
{ |
||||
ret = nc_def_dim(ncid, "longitude", nx, &xdimid); |
||||
if(ret != NC_NOERR) return "Can't create x-dimension in the netcdf file"; |
||||
ret = nc_def_dim(ncid, "latitude", ny, &ydimid); |
||||
if(ret != NC_NOERR) return "Can't create y-dimension in the netcdf file"; |
||||
ret = nc_def_var(ncid, "longitude", NC_FLOAT, 1, &xdimid, &xid); |
||||
if(ret != NC_NOERR) return "Can't create longitude variable in the netcdf file"; |
||||
ret = nc_def_var(ncid, "latitude", NC_FLOAT, 1, &ydimid, &yid); |
||||
if(ret != NC_NOERR) return "Can't create latitude variable in the netcdf file"; |
||||
AddDim("i", nx); |
||||
AddVar("longitude", NC_FLOAT, "i"); |
||||
AddVar("latitude", NC_FLOAT, "i"); |
||||
break; |
||||
} |
||||
case(G2V2): |
||||
case(RGRID): |
||||
{ |
||||
ret = nc_def_dim(ncid, "longitude", nx, &xdimid); |
||||
if(ret != NC_NOERR) return "Can't create x-dimension in the netcdf file"; |
||||
ret = nc_def_dim(ncid, "latitude", ny, &ydimid); |
||||
if(ret != NC_NOERR) return "Can't create y-dimension in the netcdf file"; |
||||
ret = nc_def_var(ncid, "longitude", NC_FLOAT, 2, dimid, &xid); |
||||
if(ret != NC_NOERR) return "Can't create longitude variable in the netcdf file"; |
||||
ret = nc_def_var(ncid, "latitude", NC_FLOAT, 2, dimid, &yid); |
||||
if(ret != NC_NOERR) return "Can't create latitude variable in the netcdf file"; |
||||
AddDim("x", nx); |
||||
AddDim("y", ny); |
||||
AddVar("x", NC_FLOAT, "x"); |
||||
AddVar("y", NC_FLOAT, "y"); |
||||
break; |
||||
} |
||||
case(GRGRID): |
||||
{ |
||||
AddDim("longitude", nx); |
||||
AddDim("latitude", ny); |
||||
AddVar("longitude", NC_FLOAT, "longitude"); |
||||
AddVar("latitude", NC_FLOAT, "latitude"); |
||||
break; |
||||
} |
||||
case(GRID): |
||||
{ |
||||
AddDim("x", nx); |
||||
AddDim("y", ny); |
||||
AddVar("x", NC_FLOAT, "y", "x"); |
||||
AddVar("y", NC_FLOAT, "y", "x"); |
||||
break; |
||||
} |
||||
case(GGRID): |
||||
{ |
||||
AddDim("longitude", nx); |
||||
AddDim("latitude", ny); |
||||
AddVar("longitude", NC_FLOAT, "latitude", "longitude"); |
||||
AddVar("latitude", NC_FLOAT, "latitude", "longitude"); |
||||
break; |
||||
} |
||||
case(UNKNOWN): return "Can't determine file type"; |
||||
} |
||||
|
||||
ret = nc_def_var_deflate(ncid, xid, 1, 1, compress); |
||||
if(ret != NC_NOERR) return "Can't set deflate parameters for longitude variable in the netcdf file"; |
||||
ret = nc_def_var_deflate(ncid, yid, 1, 1, compress); |
||||
if(ret != NC_NOERR) return "Can't set deflate parameters for latitude variable in the netcdf file"; |
||||
|
||||
text = "longitude"; |
||||
ret = nc_put_att_text(ncid, xid, "standard_name", text.Len() + 1, text.Buf()); |
||||
if(ret != NC_NOERR) return "Can't write standard_name attribute of longitude variable in the netcdf file"; |
||||
text = "latitude"; |
||||
ret = nc_put_att_text(ncid, yid, "standard_name", text.Len() + 1, text.Buf()); |
||||
if(ret != NC_NOERR) return "Can't write standard_name attribute of latitude variable in the netcdf file"; |
||||
if(IsGeoType(stype)) |
||||
{ |
||||
SetComp("longitude", compress); |
||||
SetComp("latitude", compress); |
||||
AddAtt("longitude", "standard_name", "longitude"); |
||||
AddAtt("longitude", "long_name", "Longitude"); |
||||
AddAtt("latitude", "standard_name", "latitude"); |
||||
AddAtt("latitude", "long_name", "Latitude"); |
||||
} |
||||
else |
||||
{ |
||||
SetComp("x", compress); |
||||
SetComp("y", compress); |
||||
AddAtt("x", "long_name", "x-coordinate"); |
||||
AddAtt("y", "long_name", "y-coordinate"); |
||||
} |
||||
if(!*this) return "Can't set grid in the netcdf file " + name + ": " + ErrMessage(); |
||||
|
||||
type = stype; |
||||
return ""; |
||||
} |
||||
|
||||
MString NCFileW::AddTimeData(const TimeData& tdata, bool tisindex) |
||||
{ |
||||
tdep = tisindex; |
||||
|
||||
AddDim("time", tdata.steps.size()); |
||||
AddVar("time", Type2NCType<decltype(tdata.steps)::value_type>, "time"); |
||||
SetComp("time", compress); |
||||
AddAtt("time", "standard_name", "time"); |
||||
AddAtt("time", "long_name", "time"); |
||||
AddAtt("time", "units", tdata.UnitName()); |
||||
WriteVar("time", tdata.steps.data()); |
||||
|
||||
if(!*this) return MString("Can't add time data to the netcdf file: ") + ErrMessage(); |
||||
|
||||
return ""; |
||||
} |
||||
|
||||
MString NCFileW::AddVariable(const MString& name, const MString& stname, const MString& lname, const MString& units, const MString& comment) |
||||
{ |
||||
if(type == UNKNOWN) return "File not opened"; |
||||
if(HaveVar(name)) return "Variable " + name + " already defined"; |
||||
|
||||
struct Var v(name, 0); |
||||
|
||||
int ret; |
||||
if(type == G1V1) |
||||
ret = nc_def_var(ncid, v.name.Buf(), NC_FLOAT, 1, &xdimid, &v.id); |
||||
if(Is1DType(type)) |
||||
{ |
||||
if(tdep) |
||||
AddVar(name, NC_FLOAT, "time", "i"); |
||||
else |
||||
ret = nc_def_var(ncid, v.name.Buf(), NC_FLOAT, 2, dimid, &v.id); |
||||
if(ret != NC_NOERR) return "Can't create " + v.name + " variable in the netcdf file"; |
||||
|
||||
ret = nc_def_var_deflate(ncid, v.id, 1, 1, compress); |
||||
if(ret != NC_NOERR) return "Can't set deflate parameters for " + v.name + " variable in the netcdf file"; |
||||
AddVar(name, NC_FLOAT, "i"); |
||||
} |
||||
else if(IsGeoType(type)) |
||||
{ |
||||
if(tdep) |
||||
AddVar(name, NC_FLOAT, "time", "latitude", "longitude"); |
||||
else |
||||
AddVar(name, NC_FLOAT, "latitude", "longitude"); |
||||
} |
||||
else |
||||
{ |
||||
if(tdep) |
||||
AddVar(name, NC_FLOAT, "time", "y", "x"); |
||||
else |
||||
AddVar(name, NC_FLOAT, "y", "x"); |
||||
} |
||||
|
||||
if(stname.Exist()) ret = nc_put_att_text(ncid, v.id, "standard_name", stname.Len() + 1, stname.Buf()); |
||||
if(ret != NC_NOERR) return "Can't write standard_name attribute of " + v.name + " variable in the netcdf file"; |
||||
if(lname.Exist()) ret = nc_put_att_text(ncid, v.id, "long_name", lname.Len() + 1, lname.Buf()); |
||||
if(ret != NC_NOERR) return "Can't write long_name attribute of " + v.name + " variable in the netcdf file"; |
||||
if(units.Exist()) ret = nc_put_att_text(ncid, v.id, "units", units.Len() + 1, units.Buf()); |
||||
if(ret != NC_NOERR) return "Can't write units attribute of " + v.name + " variable in the netcdf file"; |
||||
if(comment.Exist()) ret = nc_put_att_text(ncid, v.id, "comment", comment.Len() + 1, comment.Buf()); |
||||
if(ret != NC_NOERR) return "Can't write comment attribute of " + v.name + " variable in the netcdf file"; |
||||
SetComp(name, compress); |
||||
if(stname.Exist()) AddAtt(name, "standard_name", stname); |
||||
if(lname.Exist()) AddAtt(name, "long_name", lname); |
||||
if(units.Exist()) AddAtt(name, "units", units); |
||||
if(comment.Exist()) AddAtt(name, "comment", comment); |
||||
AddAtt(name, "_FillValue", fill); |
||||
|
||||
ret = nc_put_att_float(ncid, v.id, "_FillValue", NC_FLOAT, 1, &fill); |
||||
if(ret != NC_NOERR) return "Can't write _FillValue attribute of " + v.name + " variable in the netcdf file"; |
||||
if(!*this) return "Can't add variable " + name + ": " + ErrMessage(); |
||||
|
||||
vars.push_back(v); |
||||
return ""; |
||||
} |
||||
|
||||
@ -0,0 +1,231 @@
|
||||
#define MICHLIB_NOSOURCE |
||||
#include "ncsimple.h" |
||||
|
||||
RetVal<std::vector<NCSimpleFunctions::Attribute>> NCSimpleFunctions::ReadAtts(int vid) const |
||||
{ |
||||
static const MString pref = "NCSimple::ReadAtts"; |
||||
|
||||
int natt; |
||||
int ret; |
||||
if(vid == NC_GLOBAL) |
||||
ret = nc_inq_natts(ncid, &natt); |
||||
else |
||||
ret = nc_inq_var(ncid, vid, nullptr, nullptr, nullptr, nullptr, &natt); |
||||
if(ret != NC_NOERR) return Error(pref, MString("Can't inquire number of attributes: ") + nc_strerror(ret)); |
||||
|
||||
std::vector<Attribute> out; |
||||
char name[NC_MAX_NAME + 1]; |
||||
|
||||
for(int aid = 0; aid < natt; aid++) |
||||
{ |
||||
nc_type type; |
||||
size_t len; |
||||
|
||||
ret = nc_inq_attname(ncid, vid, aid, name); |
||||
if(ret != NC_NOERR) return Error(pref, MString("Can't inquire attribute name: ") + nc_strerror(ret)); |
||||
ret = nc_inq_atttype(ncid, vid, name, &type); |
||||
if(ret != NC_NOERR) return Error(pref, MString("Can't inquire attribute type: ") + nc_strerror(ret)); |
||||
ret = nc_inq_attlen(ncid, vid, name, &len); |
||||
if(ret != NC_NOERR) return Error(pref, MString("Can't inquire attribute length: ") + nc_strerror(ret)); |
||||
|
||||
if(type == NC_DOUBLE || type == NC_FLOAT) |
||||
{ |
||||
if(len == 1) |
||||
{ |
||||
double d; |
||||
ret = nc_get_att_double(ncid, vid, name, &d); |
||||
if(ret != NC_NOERR) return Error(pref, MString("Can't read attribute ") + name + ": " + nc_strerror(ret)); |
||||
out.emplace_back(MString(name), d); |
||||
} |
||||
else |
||||
{ |
||||
std::vector<double> dd(len); |
||||
ret = nc_get_att_double(ncid, vid, name, dd.data()); |
||||
if(ret != NC_NOERR) return Error(pref, MString("Can't read attribute ") + name + ": " + nc_strerror(ret)); |
||||
for(size_t i = 0; i < dd.size(); i++) out.emplace_back(MString(name) + "[" + i + "]", dd[i]); |
||||
} |
||||
} |
||||
else if(type == NC_BYTE || type == NC_SHORT || type == NC_INT || type == NC_INT64) |
||||
{ |
||||
if(len == 1) |
||||
{ |
||||
long long i; |
||||
ret = nc_get_att_longlong(ncid, vid, name, &i); |
||||
if(ret != NC_NOERR) return Error(pref, MString("Can't read attribute ") + name + ": " + nc_strerror(ret)); |
||||
out.emplace_back(MString(name), int_cast<int8>(i)); |
||||
} |
||||
else |
||||
{ |
||||
std::vector<long long> ii(len); |
||||
ret = nc_get_att_longlong(ncid, vid, name, ii.data()); |
||||
if(ret != NC_NOERR) return Error(pref, MString("Can't read attribute ") + name + ": " + nc_strerror(ret)); |
||||
for(size_t i = 0; i < ii.size(); i++) out.emplace_back(MString(name) + "[" + i + "]", int_cast<int8>(ii[i])); |
||||
} |
||||
} |
||||
else if(type == NC_UBYTE || type == NC_USHORT || type == NC_UINT || type == NC_UINT64) |
||||
{ |
||||
if(len == 1) |
||||
{ |
||||
unsigned long long u; |
||||
ret = nc_get_att_ulonglong(ncid, vid, name, &u); |
||||
if(ret != NC_NOERR) return Error(pref, MString("Can't read attribute ") + name + ": " + nc_strerror(ret)); |
||||
out.emplace_back(MString(name), int_cast<uint8>(u)); |
||||
} |
||||
else |
||||
{ |
||||
std::vector<unsigned long long> uu(len); |
||||
ret = nc_get_att_ulonglong(ncid, vid, name, uu.data()); |
||||
if(ret != NC_NOERR) return Error(pref, MString("Can't read attribute ") + name + ": " + nc_strerror(ret)); |
||||
for(size_t i = 0; i < uu.size(); i++) out.emplace_back(MString(name) + "[" + i + "]", int_cast<uint8>(uu[i])); |
||||
} |
||||
} |
||||
else if(type == NC_CHAR) |
||||
{ |
||||
std::vector<char> ss(len + 1, 0); |
||||
ret = nc_get_att_text(ncid, vid, name, ss.data()); |
||||
if(ret != NC_NOERR) return Error(pref, MString("Can't read attribute ") + name + ": " + nc_strerror(ret)); |
||||
out.emplace_back(MString(name), MString(ss.data())); |
||||
} |
||||
else |
||||
return Error(pref, MString("Unsupported type of attribute ") + name); |
||||
|
||||
// Ignore all other types
|
||||
} |
||||
|
||||
return out; |
||||
} |
||||
|
||||
Error NCSimpleFunctions::Open(const MString& filename) |
||||
{ |
||||
static const MString pref = "NCSimple::Open"; |
||||
|
||||
// Cleanup
|
||||
gats.clear(); |
||||
dims.clear(); |
||||
vars.clear(); |
||||
|
||||
nc_close(ncid); |
||||
|
||||
std::vector<Attribute> newgats; |
||||
std::vector<Dimension> newdims; |
||||
std::vector<Variable> newvars; |
||||
|
||||
char name[NC_MAX_NAME + 1]; |
||||
|
||||
// Open
|
||||
{ |
||||
auto ret = nc_open(filename.Buf(), 0, &ncid); |
||||
if(ret != NC_NOERR) return Error(pref, "Can't open file " + filename + ": " + nc_strerror(ret)); |
||||
} |
||||
|
||||
// Dimensions
|
||||
{ |
||||
int ndim; |
||||
auto ret = nc_inq_dimids(ncid, &ndim, nullptr, 1); |
||||
if(ret != NC_NOERR) return Error(pref, "Can't inquire number of dimensions in file " + filename + ": " + nc_strerror(ret)); |
||||
std::vector<int> dimids(ndim); |
||||
ret = nc_inq_dimids(ncid, nullptr, dimids.data(), 1); |
||||
if(ret != NC_NOERR) return Error(pref, "Can't inquire dimension ids in file " + filename + ": " + nc_strerror(ret)); |
||||
|
||||
size_t len; |
||||
for(const auto id: dimids) |
||||
{ |
||||
ret = nc_inq_dim(ncid, id, name, &len); |
||||
if(ret != NC_NOERR) return Error(pref, "Can't inquire dimension name and size in file " + filename + ": " + nc_strerror(ret)); |
||||
newdims.emplace_back(name, len); |
||||
} |
||||
} |
||||
|
||||
// Global attributes
|
||||
{ |
||||
auto ret = ReadAtts(NC_GLOBAL); |
||||
if(!ret) return ret.Add(pref, "Can't read global attributes in file " + filename); |
||||
newgats = std::move(ret.Value()); |
||||
} |
||||
|
||||
// Variables
|
||||
{ |
||||
int nvar; |
||||
auto ret = nc_inq_varids(ncid, &nvar, nullptr); |
||||
if(ret != NC_NOERR) return Error(pref, "Can't inquire number of variables in file " + filename + ": " + nc_strerror(ret)); |
||||
std::vector<int> varids(nvar); |
||||
ret = nc_inq_varids(ncid, nullptr, varids.data()); |
||||
if(ret != NC_NOERR) return Error(pref, "Can't inquire variables ids in file " + filename + ": " + nc_strerror(ret)); |
||||
|
||||
for(const auto vid: varids) |
||||
{ |
||||
nc_type nctype; |
||||
int ndim; |
||||
|
||||
auto ret = nc_inq_var(ncid, vid, name, &nctype, &ndim, nullptr, nullptr); |
||||
if(ret != NC_NOERR) return Error(pref, "Can't inquire variable info in file " + filename + ": " + nc_strerror(ret)); |
||||
|
||||
VarType vt = VarType::UNDEF; |
||||
|
||||
if(nctype == NC_FLOAT) vt = VarType::FLOAT; |
||||
if(nctype == NC_DOUBLE) vt = VarType::DOUBLE; |
||||
if(nctype == NC_BYTE) vt = VarType::INT1; |
||||
if(nctype == NC_SHORT) vt = VarType::INT2; |
||||
if(nctype == NC_INT) vt = VarType::INT4; |
||||
if(nctype == NC_INT64) vt = VarType::INT8; |
||||
if(nctype == NC_UBYTE) vt = VarType::UINT1; |
||||
|
||||
if(vt == VarType::UNDEF) return Error(pref, "Unsupported type of variable " + MString(name) + " in file " + filename); |
||||
|
||||
std::vector<int> dimids(ndim); |
||||
ret = nc_inq_vardimid(ncid, vid, dimids.data()); |
||||
if(ret != NC_NOERR) return Error(pref, "Can't inquire variable dimensions in file " + filename + ": " + nc_strerror(ret)); |
||||
|
||||
std::vector<size_t> dims; |
||||
char dname[NC_MAX_NAME + 1]; |
||||
for(const auto did: dimids) |
||||
{ |
||||
auto ret = nc_inq_dimname(ncid, did, dname); |
||||
if(ret != NC_NOERR) return Error(pref, "Can't inquire dimension name in file " + filename + ": " + nc_strerror(ret)); |
||||
size_t ind = newdims.size(); |
||||
for(size_t i = 0; i < newdims.size() && ind == newdims.size(); i++) |
||||
if(dname == newdims[i].Name()) ind = i; |
||||
if(ind == newdims.size()) return Error(pref, "Can't find dimension " + MString(dname) + " of variable " + name + " in file " + filename); |
||||
dims.push_back(ind); |
||||
} |
||||
|
||||
auto atts = ReadAtts(vid); |
||||
const char fname[] = "_FillValue"; |
||||
if(!atts) return Error(pref, "Can't get attributes of variable " + MString(name) + " in file " + filename); |
||||
std::vector<Attribute> vatts = std::move(atts.Value()); |
||||
|
||||
Variable::FillType fill; |
||||
if(FindInd(fname, vatts) != vatts.size()) |
||||
{ |
||||
if(nctype == NC_FLOAT || nctype == NC_DOUBLE) |
||||
{ |
||||
double d; |
||||
auto ret = nc_get_att_double(ncid, vid, fname, &d); |
||||
if(ret != NC_NOERR) return Error(pref, "Can't get fill value for the variable " + MString(name) + " in file " + filename + ": " + nc_strerror(ret)); |
||||
fill = d; |
||||
} |
||||
if(nctype == NC_BYTE || nctype == NC_SHORT || nctype == NC_INT || nctype == NC_INT64) |
||||
{ |
||||
long long l; |
||||
auto ret = nc_get_att_longlong(ncid, vid, fname, &l); |
||||
if(ret != NC_NOERR) return Error(pref, "Can't get fill value for the variable " + MString(name) + " in file " + filename + ": " + nc_strerror(ret)); |
||||
fill = int_cast<int8>(l); |
||||
} |
||||
if(nctype == NC_UBYTE || nctype == NC_USHORT || nctype == NC_UINT || nctype == NC_UINT64) |
||||
{ |
||||
unsigned long long u; |
||||
auto ret = nc_get_att_ulonglong(ncid, vid, fname, &u); |
||||
if(ret != NC_NOERR) return Error(pref, "Can't get fill value for the variable " + MString(name) + " in file " + filename + ": " + nc_strerror(ret)); |
||||
fill = int_cast<uint8>(u); |
||||
} |
||||
} |
||||
|
||||
newvars.emplace_back(name, vt, std::move(dims), std::move(vatts), fill); |
||||
} |
||||
} |
||||
|
||||
gats = std::move(newgats); |
||||
dims = std::move(newdims); |
||||
vars = std::move(newvars); |
||||
return Error(); |
||||
} |
||||
@ -0,0 +1,250 @@
|
||||
#define MICHLIB_NOSOURCE |
||||
#include "zarr.h" |
||||
#include "copcat.h" |
||||
#include <blosc.h> |
||||
|
||||
std::vector<MString> ZarrFunctions::ReadVarNames(const Json::Value& meta) |
||||
{ |
||||
std::vector<MString> out; |
||||
if(meta.type() != Json::objectValue) return out; |
||||
const auto keys = meta.getMemberNames(); |
||||
for(const auto& key: keys) |
||||
{ |
||||
if(!key.ends_with("/.zarray")) continue; |
||||
const auto vname = key.substr(0, key.size() - 8); |
||||
const auto& zattr = meta[vname + "/.zattrs"]; |
||||
if(!(zattr && zattr.type() == Json::objectValue)) continue; |
||||
|
||||
MString name(vname.c_str(), vname.size()); |
||||
bool found = false; |
||||
for(size_t id = 0; id < out.size(); id++) |
||||
if(out[id] == name) |
||||
{ |
||||
found = true; |
||||
break; |
||||
} |
||||
if(!found) out.emplace_back(std::move(name)); |
||||
} |
||||
return out; |
||||
} |
||||
|
||||
Error ZarrFunctions::AddVar(const MString& name, const Json::Value& zattrs, const Json::Value& zarray) |
||||
{ |
||||
static const MString pref = "Zarr::AddVar"; |
||||
|
||||
VarType newtype; |
||||
|
||||
Variable::FillType fill; |
||||
|
||||
// Checks for parameters in zarray
|
||||
{ |
||||
const auto& cid = zarray["compressor"]["id"]; |
||||
if(!cid || cid.type() != Json::stringValue || cid.asString() != "blosc") return {pref, "Unsupported compressor: " + MString(cid.asString().c_str())}; |
||||
} |
||||
{ |
||||
const auto& zf = zarray["zarr_format"]; |
||||
if(!zf || (zf.type() != Json::uintValue && zf.type() != Json::intValue) || zf.asUInt() != 2) return {pref, "Unsupported format version: " + MString(zf.asUInt())}; |
||||
} |
||||
{ |
||||
const auto& ord = zarray["order"]; |
||||
if(!ord || ord.type() != Json::stringValue || ord.asString() != "C") return {pref, "Order in not C"}; |
||||
} |
||||
{ |
||||
const auto& f = zarray["filters"]; |
||||
if(f.type() != Json::nullValue) return {pref, "Filters is not null"}; |
||||
} |
||||
|
||||
// Read dtype
|
||||
{ |
||||
const auto& dtype = zarray["dtype"]; |
||||
if(!dtype || dtype.type() != Json::stringValue) return {pref, "No datatype"}; |
||||
const auto str = dtype.asString(); |
||||
if(str == "<f4") |
||||
newtype = VarType::FLOAT; |
||||
else if(str == "<f8") |
||||
newtype = VarType::DOUBLE; |
||||
else if(str == "|i1") |
||||
newtype = VarType::INT1; |
||||
else if(str == "|u1") |
||||
newtype = VarType::UINT1; |
||||
else if(str == "<i2") |
||||
newtype = VarType::INT2; |
||||
else if(str == "<i4") |
||||
newtype = VarType::INT4; |
||||
else if(str == "<i8") |
||||
newtype = VarType::INT8; |
||||
else |
||||
return {pref, "Unsupported datatype: " + MString(str.c_str())}; |
||||
} |
||||
|
||||
// Read fill_value
|
||||
{ |
||||
const auto& fillval = zarray["fill_value"]; |
||||
if(!fillval) // return {pref, "No fillval"};
|
||||
fill = 0; |
||||
else if(fillval.type() == Json::uintValue) |
||||
fill = fillval.asUInt64(); |
||||
else if(fillval.type() == Json::intValue) |
||||
fill = fillval.asInt64(); |
||||
else if(fillval.type() == Json::realValue) |
||||
fill = fillval.asDouble(); |
||||
else if(fillval.type() == Json::stringValue && fillval.asString() == "NaN") |
||||
fill = NAN; |
||||
} |
||||
|
||||
// Read attributes
|
||||
auto atts = ReadAtts(zattrs); |
||||
|
||||
std::vector<MString> dnames; |
||||
std::vector<size_t> dsizes; |
||||
std::vector<size_t> csizes; |
||||
std::vector<size_t> dids; |
||||
|
||||
// Read dimensions names
|
||||
{ |
||||
const auto& arrdim = zattrs["_ARRAY_DIMENSIONS"]; |
||||
if(!(arrdim && arrdim.type() == Json::arrayValue)) return {pref, "_ARRAY_DIMENSIONS not found"}; |
||||
for(Json::ArrayIndex i = 0; i < arrdim.size(); i++) |
||||
if(const auto& dim = arrdim[i]; dim.type() == Json::stringValue) |
||||
{ |
||||
const auto val = dim.asString(); |
||||
dnames.emplace_back(val.c_str(), val.size()); |
||||
} |
||||
} |
||||
|
||||
// Read dimensions sizes
|
||||
{ |
||||
const auto& shape = zarray["shape"]; |
||||
if(!(shape && shape.type() == Json::arrayValue)) return {pref, "shape not found"}; |
||||
for(Json::ArrayIndex i = 0; i < shape.size(); i++) |
||||
if(const auto& s = shape[i]; s.type() == Json::uintValue || s.type() == Json::intValue) dsizes.push_back(s.asUInt()); |
||||
} |
||||
|
||||
// Read chunk sizes
|
||||
{ |
||||
const auto& chunk = zarray["chunks"]; |
||||
if(!(chunk && chunk.type() == Json::arrayValue)) return {pref, "chunks not found"}; |
||||
for(Json::ArrayIndex i = 0; i < chunk.size(); i++) |
||||
if(const auto& c = chunk[i]; c.type() == Json::uintValue || c.type() == Json::intValue) csizes.push_back(c.asUInt()); |
||||
} |
||||
|
||||
if(dnames.size() != dsizes.size() || dnames.size() != csizes.size()) return {pref, "shape and chunks are in contradiction"}; |
||||
|
||||
dids.resize(dnames.size()); |
||||
|
||||
// Check dimensions names and sizes
|
||||
for(size_t i = 0; i < dnames.size(); i++) |
||||
{ |
||||
bool found = false; |
||||
for(size_t id = 0; id < dims.size(); id++) |
||||
if(dims[id].Name() == dnames[i]) |
||||
{ |
||||
found = true; |
||||
if(dims[id].Size() != dsizes[i]) |
||||
return {pref, "According to previous data, the dimension " + dnames[i] + " has a size of " + dims[id].Size() + ", but here it is defined as " + dsizes[i]}; |
||||
dids[i] = id; |
||||
break; |
||||
} |
||||
if(!found) |
||||
{ |
||||
dids[i] = dims.size(); |
||||
dims.emplace_back(dnames[i], dsizes[i]); |
||||
} |
||||
} |
||||
|
||||
vars.emplace_back(name, newtype, std::move(dids), std::move(atts), fill); |
||||
chunks.push_back(std::move(csizes)); |
||||
|
||||
return Error(); |
||||
} |
||||
|
||||
Error ZarrFunctions::GetChunk(const MString& var, const std::vector<size_t>& chunkind, size_t chunksize, size_t elsize, void* data, const void* fill) const |
||||
{ |
||||
static const MString pref = "Zarr::GetChunk"; |
||||
|
||||
MString str = url + "/" + var + "/"; |
||||
for(size_t i = 0; i < chunkind.size(); i++) str += (i == 0 ? "" : ".") + MString(chunkind[i]); |
||||
|
||||
auto [content, suc] = cache->Get(str); |
||||
|
||||
if(!suc) |
||||
{ |
||||
michlib::message(str + " not found in cache, downloading"); |
||||
CURLRAII myhandle; // TODO: remove this workaround of unknown bug
|
||||
if(proxyurl.Exist()) curl_easy_setopt(myhandle, CURLOPT_PROXY, proxyurl.Buf()); |
||||
//auto [out, res] = GetUrl(chandle, str);
|
||||
auto [out, res] = GetUrl(myhandle, str); |
||||
if(res != CURLE_OK) return Error(pref, MString("can't download chunk: ") + chandle.Err()); |
||||
long respcode; |
||||
//curl_easy_getinfo(chandle, CURLINFO_RESPONSE_CODE, &respcode);
|
||||
curl_easy_getinfo(myhandle, CURLINFO_RESPONSE_CODE, &respcode); |
||||
if(respcode == 403) out = ""; // Failed chunk download mean that this chunk contains only fill
|
||||
cache->Put(str, out, 3600); |
||||
content = std::move(out); |
||||
} |
||||
|
||||
if(content.Exist()) |
||||
{ |
||||
size_t nb, cb, bs; |
||||
blosc_cbuffer_sizes(content.Buf(), &nb, &cb, &bs); |
||||
if(cb != content.Len()) return Error(pref, MString("bytes download: ") + content.Len() + ", but compressed bytes " + cb); |
||||
if(nb != chunksize * elsize) return Error(pref, MString("decompressed bytes: ") + nb + ", but buffer size " + chunksize * elsize); |
||||
auto res = blosc_decompress_ctx(content.Buf(), data, chunksize * elsize, 1); |
||||
if(int_cast<size_t>(res) != chunksize * elsize) return Error(pref, MString("decompress only ") + res + " bytes of " + chunksize * elsize); |
||||
} |
||||
else |
||||
{ |
||||
if(fill == nullptr) return Error(pref, MString("can't download chunk: ") + chandle.Err()); |
||||
for(size_t i = 0; i < chunksize; i++) memcpy(michlib::P1(data) + i * elsize, fill, elsize); |
||||
} |
||||
|
||||
return Error(); |
||||
} |
||||
|
||||
Error ZarrFunctions::Open(const MString& product, const MString& dataset, bool time) |
||||
{ |
||||
static const MString pref = "Zarr::Open"; |
||||
|
||||
gats.clear(); |
||||
dims.clear(); |
||||
vars.clear(); |
||||
|
||||
CopernicusCatalog cat; |
||||
Json::Value json; |
||||
|
||||
MString realdataset; |
||||
if(!dataset.Exist()) |
||||
{ |
||||
auto dsets = cat.DatasetList(product); |
||||
if(!dsets) return dsets.Add(pref, "Can't get default dataset of product " + product); |
||||
realdataset = dsets.Value()[0]; |
||||
} |
||||
else |
||||
realdataset = dataset; |
||||
|
||||
{ |
||||
auto urlret = time ? cat.DatasetTimeURL(product, realdataset) : cat.DatasetGeoURL(product, realdataset); |
||||
if(!urlret) return urlret.Add(pref, "Can't get url for the dataset " + realdataset + " of product " + product); |
||||
url = urlret.Value(); |
||||
|
||||
auto ret = cat.GetJSON(url + "/.zmetadata"); |
||||
if(ret) |
||||
json = ret.Value(); |
||||
else |
||||
return ret.Add(pref, "can't download .zmetadata"); |
||||
} |
||||
|
||||
const auto& meta = json["metadata"]; |
||||
if(!meta) return {pref, "No \"metadata\" key in JSON data"}; |
||||
|
||||
if(meta[".zattrs"]) gats = ReadAtts(meta[".zattrs"]); |
||||
auto vnames = ReadVarNames(meta); |
||||
|
||||
for(size_t i = 0; i < vnames.size(); i++) |
||||
{ |
||||
auto err = AddVar(vnames[i], meta[(vnames[i] + "/.zattrs").Buf()], meta[(vnames[i] + "/.zarray").Buf()]); |
||||
if(!err) return err.Add(pref, "Can't init variable " + vnames[i]); |
||||
} |
||||
|
||||
return Error(); |
||||
} |
||||
Loading…
Reference in new issue