|
|
|
|
@ -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; |
|
|
|
|
} |
Michael Uleysky