5 changed files with 334 additions and 2 deletions
@ -0,0 +1,60 @@
|
||||
#pragma once |
||||
#include "ParseArgs.h" |
||||
#include "basedata.h" |
||||
#include "mdatetime.h" |
||||
#include "mregex.h" |
||||
#include "DataAdapters/ncfilealt.h" |
||||
#include <dirent.h> |
||||
#include <set> |
||||
#include <utility> |
||||
|
||||
using michlib::GPL; |
||||
using michlib::MDateTime; |
||||
using michlib::RegExp; |
||||
using michlib::uint2; |
||||
using michlib::uint4; |
||||
using michlib::Floor; |
||||
using michlib::Ceil; |
||||
using michlib::pointer_cast; |
||||
|
||||
class MODISBINLOCALData |
||||
{ |
||||
MString datapath; |
||||
std::vector<MString> dataset; |
||||
std::vector<MDateTime> times; |
||||
|
||||
bool SufFromDataset(const MString& datasetstr); |
||||
|
||||
struct Parameters: public BaseParameters |
||||
{ |
||||
real lonb, latb, lone, late; |
||||
virtual ~Parameters() override = default; |
||||
}; |
||||
|
||||
public: |
||||
using Data = UngriddedData; |
||||
|
||||
MString Info() const; |
||||
MString Open(const CLArgs& args); |
||||
Data Read(const MString& var, const BaseParameters* ip, size_t tind) const; |
||||
Data ReadFile(const MString& suf, const struct Parameters* p, size_t tind) const; |
||||
|
||||
std::pair<const BaseParameters*, MString> Parameters(michlib_internal::ParameterListEx& pars, const CLArgs& args) const; |
||||
|
||||
bool CheckVar(const MString& vname) const |
||||
{ |
||||
if(dataset.empty()) return false; |
||||
bool ischl = dataset.front() == "A_CHL" || dataset.front() == "T_CHL"; |
||||
return (vname == "chl" && ischl) || (vname == "temp" && !ischl); |
||||
} |
||||
|
||||
bool isOk() const { return times.size() > 0; } |
||||
|
||||
size_t NTimes() const { return times.size(); } |
||||
|
||||
MDateTime Time(size_t i) const |
||||
{ |
||||
if(!isOk() || i >= times.size()) return MDateTime(); |
||||
return times[i]; |
||||
} |
||||
}; |
||||
@ -0,0 +1,262 @@
|
||||
#define MICHLIB_NOSOURCE |
||||
#include "MODISBINLOCAL.h" |
||||
|
||||
MString MODISBINLOCALData::Info() const |
||||
{ |
||||
MString ds = ""; |
||||
for(size_t i = 0; i < dataset.size(); i++) |
||||
{ |
||||
if(i != 0) ds += " "; |
||||
ds += dataset[i]; |
||||
} |
||||
|
||||
// clang-format off
|
||||
return |
||||
"Dataset: MODIS Level-3 Binned Data 4.6 km - " + ds + "\n" + |
||||
" Begin date: " + Time(0).ToString() + "\n" + |
||||
" End date: " + Time(NTimes()-1).ToString() + "\n" + |
||||
" Time moments: " + NTimes() + "\n" + |
||||
" Supported variables: temp, chl"; |
||||
// clang-format on
|
||||
} |
||||
|
||||
MString MODISBINLOCALData::Open(const CLArgs& args) |
||||
{ |
||||
GPL.UsePrefix("MODISBINLOCAL"); |
||||
datapath = GPL.ParameterSValue("Datapath", ""); |
||||
MString var, datasetstr; |
||||
if(!args.contains("var") && !args.contains("dataset")) var = "temp"; |
||||
if(args.contains("var")) |
||||
{ |
||||
var = args.at("var"); |
||||
if(var != "temp" && var != "chl") return "Unknown variable: " + var; |
||||
} |
||||
if(var.Exist() && !args.contains("dataset")) datasetstr = (var == "temp") ? "ALL" : "CHL"; |
||||
if(args.contains("dataset")) datasetstr = args.at("dataset"); |
||||
|
||||
bool ischl = SufFromDataset(datasetstr); |
||||
if(dataset.size() == 0) return "Incorrect dataset: " + datasetstr; |
||||
if(!var.Exist()) var = ischl ? "chl" : "temp"; |
||||
if((var == "temp" && ischl) || (var == "chl" && !ischl)) return "Arguments dataset and var are in contradiction"; |
||||
|
||||
MString regstring = "([0-9]{4}-[0-9]{2}-[0-9]{2})_("; |
||||
for(size_t i = 0; i < dataset.size(); i++) |
||||
{ |
||||
if(i != 0) regstring += "|"; |
||||
regstring += dataset[i]; |
||||
} |
||||
regstring += ").nc"; |
||||
|
||||
RegExp regex(regstring.Buf()); |
||||
regex.Compile(); |
||||
|
||||
DIR* dir = opendir(datapath.Buf()); |
||||
struct dirent* de; |
||||
if(nullptr == dir) return "Can't open directory " + datapath; |
||||
while((de = readdir(dir))) |
||||
{ |
||||
if(!regex.Match(de->d_name)) continue; |
||||
times.emplace_back(MString(de->d_name + regex.Off(1), regex.Len(1))); |
||||
} |
||||
closedir(dir); |
||||
std::sort(times.begin(), times.end()); |
||||
auto last = std::unique(times.begin(), times.end()); |
||||
times.erase(last, times.end()); |
||||
times.shrink_to_fit(); |
||||
return ""; |
||||
} |
||||
|
||||
std::pair<const BaseParameters*, MString> MODISBINLOCALData::Parameters(michlib_internal::ParameterListEx& pars, const CLArgs& args) const |
||||
{ |
||||
std::unique_ptr<struct Parameters> ppar{new struct Parameters}; |
||||
if(!(args.contains("lonb") && args.contains("lone") && args.contains("latb") && args.contains("late"))) return {nullptr, "Region not specified (lonb, lone, latb, late)"}; |
||||
ppar->lonb = args.at("lonb").ToReal(); |
||||
ppar->lone = args.at("lone").ToReal(); |
||||
ppar->latb = args.at("latb").ToReal(); |
||||
ppar->late = args.at("late").ToReal(); |
||||
|
||||
pars.SetParameter("lonb", ppar->lonb); |
||||
pars.SetParameter("latb", ppar->latb); |
||||
pars.SetParameter("lone", ppar->lone); |
||||
pars.SetParameter("late", ppar->late); |
||||
|
||||
return {ppar.release(), ""}; |
||||
} |
||||
|
||||
MODISBINLOCALData::Data MODISBINLOCALData::Read(const MString& var, const BaseParameters* ip, size_t tind) const |
||||
{ |
||||
auto p = pointer_cast<const struct Parameters*>(ip); |
||||
|
||||
Averager<Data> out; |
||||
|
||||
for(const auto& suf: dataset) |
||||
{ |
||||
Data dat = ReadFile(suf, p, tind); |
||||
if(dat) out.Add(dat); |
||||
} |
||||
return out.Div(); |
||||
} |
||||
|
||||
MODISBINLOCALData::Data MODISBINLOCALData::ReadFile(const MString& suf, const struct MODISBINLOCALData::Parameters* p, size_t tind) const |
||||
{ |
||||
michlib::NCFileA nc; |
||||
michlib::message("Reading " + datapath + "/" + times[tind].ToString() + "_" + suf + ".nc"); |
||||
nc.Reset(datapath + "/" + times[tind].ToString() + "_" + suf + ".nc"); |
||||
if(!nc) return Data(); |
||||
MString dname; |
||||
if(suf == "A_CHL" || suf == "T_CHL") dname = "chlor_a"; |
||||
if(suf == "A_SST" || suf == "T_SST") dname = "sst"; |
||||
if(suf == "A_NSST" || suf == "T_NSST") dname = "sst"; |
||||
if(suf == "A_SST4" || suf == "T_SST4") dname = "sst4"; |
||||
|
||||
struct BinIndexType |
||||
{ |
||||
uint4 start_num; |
||||
uint4 begin; |
||||
uint4 extent; |
||||
uint4 max; |
||||
}; |
||||
struct BinListType |
||||
{ |
||||
uint4 bin_num; |
||||
uint2 nobs; |
||||
uint2 nscenes; |
||||
float weights; |
||||
float time_rec; |
||||
}; |
||||
struct BinDataType |
||||
{ |
||||
float sum; |
||||
float sum_squared; |
||||
}; |
||||
|
||||
std::vector<uint4> beg; |
||||
std::vector<uint4> loc; |
||||
std::vector<float> val; |
||||
|
||||
std::vector<uint4> bins; |
||||
|
||||
// Read index
|
||||
{ |
||||
auto ind = nc.V<struct BinIndexType>("/level-3_binned_data/BinIndex"); |
||||
if(!ind) return Data(); |
||||
beg.resize(ind.DimLen(0) + 1); |
||||
beg[0] = 1; |
||||
for(size_t i = 1; i <= ind.DimLen(0); i++) beg[i] = beg[i - 1] + ind(i - 1).max; |
||||
} |
||||
|
||||
size_t ilatb, ilate; |
||||
{ |
||||
size_t nr = beg.size() - 2; |
||||
ilatb = static_cast<size_t>(Ceil((p->latb + 90.0) * nr / 180.0 - 0.5)); |
||||
ilate = static_cast<size_t>(Floor((p->late + 90.0) * nr / 180.0 - 0.5)); |
||||
} |
||||
|
||||
for(size_t i = ilatb; i <= ilate; i++) |
||||
{ |
||||
real lb = michlib::ToGeoDomainNeg(p->lonb); |
||||
real le = michlib::ToGeoDomainNeg(p->lone); |
||||
size_t nl = beg[i + 1] - beg[i]; |
||||
size_t ilonb, ilone, il; |
||||
ilonb = static_cast<size_t>(Ceil((lb + 180.0) * nl / 360.0 - 0.5)); |
||||
ilone = static_cast<size_t>(Floor((le + 180.0) * nl / 360.0 - 0.5)); |
||||
il = ilonb; |
||||
do bins.push_back(beg[i] + il); |
||||
while(++il % nl != ilone + 1); |
||||
} |
||||
std::sort(bins.begin(), bins.end()); |
||||
|
||||
auto InReg = [&bins = std::as_const(bins)](uint4 n) -> size_t |
||||
{ |
||||
size_t a = 0, b = bins.size() - 1; |
||||
if(bins[a] == n) return a; |
||||
if(bins[b] == n) return b; |
||||
while(true) |
||||
{ |
||||
size_t c = (a + b) / 2; |
||||
if(bins[c] == n) return c; |
||||
if(n < bins[c]) |
||||
b = c; |
||||
else |
||||
a = c; |
||||
if(b - a == 1) return bins.size(); |
||||
} |
||||
}; |
||||
|
||||
auto Row = [&beg = std::as_const(beg)](uint4 n) -> size_t |
||||
{ |
||||
if(n >= beg.back() || n == 0) return beg.size() - 1; |
||||
size_t a = 0, b = beg.size() - 2; |
||||
if(n >= beg[a] && n < beg[a + 1]) return a; |
||||
if(n >= beg[b] && n < beg[b + 1]) return b; |
||||
while(true) |
||||
{ |
||||
size_t c = (a + b) / 2; |
||||
if(n >= beg[c] && n < beg[c + 1]) return c; |
||||
if(n < beg[c]) |
||||
b = c; |
||||
else |
||||
a = c; |
||||
} |
||||
}; |
||||
auto Lat = [&Row = Row](uint4 n) |
||||
{ |
||||
auto r = Row(n); |
||||
auto nr = Row(0) - 1; // Number of rows
|
||||
return -90.0 + (r + 0.5) * (180.0 / nr); |
||||
}; |
||||
auto Lon = [&Row = Row, &beg = std::as_const(beg)](uint4 n) -> real |
||||
{ |
||||
auto r = Row(n); |
||||
if(r == beg.size() - 1) return NAN; |
||||
return -180.0 + (n - beg[r] + 0.5) * (360.0 / (beg[r + 1] - beg[r])); |
||||
}; |
||||
|
||||
{ |
||||
auto lst = nc.V<struct BinListType>("/level-3_binned_data/BinList"); |
||||
if(!lst) return Data(); |
||||
auto dat = nc.V<struct BinDataType>("/level-3_binned_data/" + dname); |
||||
if(!dat) return Data(); |
||||
if(dat.DimLen(0) != lst.DimLen(0)) return Data(); |
||||
val.resize(dat.DimLen(0)); |
||||
loc.resize(lst.DimLen(0)); |
||||
for(size_t i = 0; i < dat.DimLen(0); i++) |
||||
{ |
||||
loc[i] = lst(i).bin_num; |
||||
val[i] = dat(i).sum / lst(i).weights; |
||||
} |
||||
} |
||||
if(loc.size() != val.size()) return Data(); |
||||
|
||||
Data out( |
||||
bins.size(), [&bins = std::as_const(bins), &Lon = Lon](size_t i) { return Lon(bins[i]); }, [&bins = std::as_const(bins), &Lat = Lat](size_t i) { return Lat(bins[i]); }); |
||||
for(size_t i = 0; i < out.N(); i++) out(i) = out.Fillval(); |
||||
|
||||
for(size_t i = 0; i < loc.size(); i++) |
||||
{ |
||||
auto ind = InReg(loc[i]); |
||||
if(ind < bins.size()) out.V(ind) = val[i]; |
||||
} |
||||
return out; |
||||
} |
||||
|
||||
bool MODISBINLOCALData::SufFromDataset(const MString& datasetstr) |
||||
{ |
||||
std::set<MString> sufs; |
||||
auto words = michlib::Split_on_words(datasetstr, ", \t\x00", false); |
||||
bool ischl = words.front() == "A_CHL" || words.front() == "T_CHL" || words.front() == "CHL"; |
||||
|
||||
for(auto ci = words.cbegin(); ci != words.cend(); ++ci) |
||||
{ |
||||
if(ischl && (*ci == "A_CHL" || *ci == "CHL")) sufs.insert("A_CHL"); |
||||
if(ischl && (*ci == "T_CHL" || *ci == "CHL")) sufs.insert("T_CHL"); |
||||
if(!ischl && (*ci == "A_SST" || *ci == "SST" || *ci == "DSST" || *ci == "A_DSST" || *ci == "A_ALL" || *ci == "ALL")) sufs.insert("A_SST"); |
||||
if(!ischl && (*ci == "T_SST" || *ci == "SST" || *ci == "DSST" || *ci == "T_DSST" || *ci == "T_ALL" || *ci == "ALL")) sufs.insert("T_SST"); |
||||
if(!ischl && (*ci == "A_SST4" || *ci == "SST4" || *ci == "DSST" || *ci == "A_DSST" || *ci == "A_ALL" || *ci == "ALL")) sufs.insert("A_SST4"); |
||||
if(!ischl && (*ci == "T_SST4" || *ci == "SST4" || *ci == "DSST" || *ci == "T_DSST" || *ci == "T_ALL" || *ci == "ALL")) sufs.insert("T_SST4"); |
||||
if(!ischl && (*ci == "A_NSST" || *ci == "NSST" || *ci == "A_ALL" || *ci == "ALL")) sufs.insert("A_NSST"); |
||||
if(!ischl && (*ci == "T_NSST" || *ci == "NSST" || *ci == "T_ALL" || *ci == "ALL")) sufs.insert("T_NSST"); |
||||
} |
||||
for(const auto& w: sufs) dataset.push_back(w); |
||||
return ischl; |
||||
} |
||||
Loading…
Reference in new issue