Browse Source

Making averaging more generic

interpolate
Michael Uleysky 2 years ago
parent
commit
14118a1d53
  1. 10
      include/NEMO.h
  2. 43
      include/basedata.h
  3. 48
      include/simple2ddata.h

10
include/NEMO.h

@ -92,7 +92,7 @@ class NEMOData
auto unit = f.A<MString>(name, "units");
if(unit && unit.Get() == "m s-1") unitmul = 100.0;
Data data((xb < xe) ? (xe - xb + 1) : (nx + xe - xb + 1), ye - yb + 1, Lonb(), Latb());
Data data((xb < xe) ? (xe - xb + 1) : (nx + xe - xb + 1), ye - yb + 1, Lonb(), Latb(), 1.0 / 12.0, 1.0 / 12.0);
if(xb < xe)
{
@ -154,17 +154,16 @@ class NEMOData
Data ReadVar(const MString& name, const std::vector<size_t>& tindex) const
{
Data out;
Averager<Data> out;
if(tindex.size() == 0 || !isOk()) return out;
std::vector<size_t> count;
for(size_t i = 0; i < tindex.size(); i++)
{
Data dat = ReadVar(name, tindex[i]);
if(!dat) return Data();
out.Add(dat, count);
out.Add(dat);
}
out.Div(count);
out.Div();
return out;
}
@ -391,4 +390,3 @@ template<> inline NEMOData::Data NEMOData::Read<vartype::Vartype::SAL>(const std
template<> inline NEMOData::Data NEMOData::Read<vartype::Vartype::MLD>(const std::vector<size_t>& tindex) const { return ReadVar("mlotst", tindex); }
template<> inline NEMOData::Data NEMOData::Read<vartype::Vartype::SSH>(const std::vector<size_t>& tindex) const { return ReadVar("zos", tindex); }
template<> inline NEMOData::Data NEMOData::Read<vartype::Vartype::W>(const std::vector<size_t>& tindex) const { return ReadVar("wo", tindex); }

43
include/basedata.h

@ -24,6 +24,47 @@ class BaseData
static real Fillval() { return fillval; }
explicit operator bool() const { return data.size() != 0; }
explicit operator bool() const { return N() != 0; }
};
template<class Data> class Averager: public Data
{
std::vector<size_t> count;
public:
void Add(const Data& d)
{
if(!d) return;
if(!*this)
{
*static_cast<Data*>(this) = d;
count.resize(Data::N());
for(size_t i = 0; i < Data::N(); i++) count[i] = (Data::V(i) == Data::Fillval()) ? 0 : 1;
return;
}
if(Data::N() != d.N()) return;
for(size_t i = 0; i < Data::N(); i++)
if(d(i) != Data::Fillval())
{
if(Data::V(i) == Data::Fillval())
{
Data::V(i) = d(i);
count[i] = 1;
}
else
{
Data::V(i) += d(i);
count[i]++;
}
}
}
Averager& Div()
{
if(!*this) return *this;
for(size_t i = 0; i < Data::N(); i++)
if(count[i] != 0) Data::V(i) /= count[i];
return *this;
}
};

48
include/simple2ddata.h

@ -5,11 +5,12 @@ class Simple2DData: public BaseData
{
real x0 = 0.0, y0 = 0.0;
size_t nx = 0, ny = 0;
real xstep = 0.0, ystep = 0.0;
public:
Simple2DData() = default;
Simple2DData(size_t nx_, size_t ny_, real x0_, real y0_): BaseData(nx_ * ny_), x0(x0_), y0(y0_), nx(nx_), ny(ny_) {}
Simple2DData(size_t nx_, size_t ny_, real x0_, real y0_, real xs_, real ys_): BaseData(nx_ * ny_), x0(x0_), y0(y0_), nx(nx_), ny(ny_), xstep(xs_), ystep(ys_) {}
const real& V(size_t i) const { return BaseData::V(i); }
real& V(size_t i) { return BaseData::V(i); }
@ -26,50 +27,9 @@ class Simple2DData: public BaseData
size_t Nx() const { return nx; }
size_t Ny() const { return ny; }
real Lon(size_t ix, size_t iy) const { return x0 + ix / 12.0; }
real Lat(size_t ix, size_t iy) const { return y0 + iy / 12.0; }
real Lon(size_t ix, [[maybe_unused]] size_t iy) const { return x0 + ix * xstep; }
real Lat([[maybe_unused]] size_t ix, size_t iy) const { return y0 + iy * ystep; }
real Lon(size_t i) const { return Lon(i % nx, i / nx); }
real Lat(size_t i) const { return Lat(i % nx, i / nx); }
void InitCount(std::vector<size_t>& count) const
{
count.resize(data.size());
for(size_t i = 0; i < N(); i++) count[i] = (V(i) == fillval) ? 0 : 1;
}
Simple2DData& Add(const Simple2DData& d, std::vector<size_t>& count)
{
if(!d) return *this;
if(!*this)
{
*this = d;
InitCount(count);
}
if(N() != d.N() || N() != count.size()) return *this;
for(size_t i = 0; i < N(); i++)
if(d.V(i) != fillval)
{
if(V(i) == fillval)
{
V(i) = d(i);
count[i] = 1;
}
else
{
V(i) += d(i);
count[i]++;
}
}
return *this;
}
Simple2DData& Div(const std::vector<size_t>& count)
{
if(N() != count.size()) return *this;
for(size_t i = 0; i < N(); i++)
if(count[i] != 0) V(i) /= count[i];
return *this;
}
};

Loading…
Cancel
Save