michael_ul
/
odm
1
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity
Browse Source

Added new action GENINTFILE

interpolate
Michael Uleysky 2 years ago
parent
d844b0d014
commit
8ca3fa7bb2
6 changed files with 167 additions and 40 deletions
Whitespace
Show all changes Ignore whitespace when comparing lines Ignore changes in amount of whitespace Ignore changes in whitespace at EOL
Unified View
Diff Options
Show Stats Download Patch File Download Diff File
  1. 101
      include/actiongenintfile.h
  2. 50
      include/actions.h
  3. 1
      include/odm.h
  4. 30
      include/traits.h
  5. 17
      include/uvdata.h
  6. 8
      src/odm.cpp

101
include/actiongenintfile.h
Unescape View File

@ -0,0 +1,101 @@
#pragma once
#include "actions.h"
using michlib::Cos;
template<class D> struct DoAction_<D, ActionID::GENINTFILE>
{
static MString DoAction(const CLArgs& args, D& data);
};
template<class D> MString DoAction_<D, ActionID::GENINTFILE>::DoAction(const CLArgs& args, D& ds)
{
michlib_internal::ParameterListEx pars;
pars.UsePrefix("");
pars.SetParameter("source", args.at("source"));
if(!args.contains("out")) return "Output file name not specified";
auto [tindexes, err] = GetTIndexes(ds, args, pars);
if(err.Exist()) return err;
if(tindexes.size() < 10) return "Too few time moments";
// Check isotropy
for(size_t it = 0; it < tindexes.size() - 1; it++)
if(ds.Time(tindexes[it + 1]) - ds.Time(tindexes[it]) != ds.Time(tindexes[1]) - ds.Time(tindexes[0])) return "The time step is not isotropic";
std::unique_ptr<const BaseParameters> sourcepars;
if constexpr(ParametersSupported<D>)
{
auto [p, err] = ds.Parameters(pars, args);
if(err.Exist()) return err;
sourcepars.reset(p);
}
auto p = sourcepars.get();
MString name = args.at("out");
BFileW fw;
fw.Create(name, 2);
fw.SetColumnName(1, "u");
fw.SetColumnName(2, "v");
for(size_t it = 0; it < tindexes.size(); it++)
{
auto data = ReadUV(ds, p, tindexes[it]);
if(!data) return "Can't read data";
if(it == 0)
{
real lonb = data.Lon(0, 0), latb = data.Lat(0, 0), lone = data.Lon(data.Nx() - 1, data.Ny() - 1), late = data.Lat(data.Nx() - 1, data.Ny() - 1);
fw.UsePrefix("");
fw.SetParameter("nx", data.Nx());
fw.SetParameter("ny", data.Ny());
fw.SetParameter("dx", data.XStep() * 60.0);
fw.SetParameter("dy", data.YStep() * 60.0);
fw.SetParameter("nt", tindexes.size());
fw.SetParameter("dt", (ds.Time(tindexes[1]) - ds.Time(tindexes[0])) / 86400.0);
fw.SetParameter("FillValue", data.Fillval());
fw.UsePrefix("Info");
fw.SetParameter("lonb", lonb);
fw.SetParameter("latb", latb);
fw.SetParameter("lone", lone);
fw.SetParameter("late", late);
fw.SetParameter("DataSource", ds.DataTitle());
fw.SetParameter("xy2lon", MString(lonb) + "+%x/60");
fw.SetParameter("xy2lon_rp", MString("POP 60 DIV ") + lonb + " ADD");
fw.SetParameter("lonlat2x", MString("(%lon-") + lonb + ")*60");
fw.SetParameter("xy2lat", MString(latb) + "+%y/60");
fw.SetParameter("xy2lat_rp", MString("EXCH POP 60 DIV ") + latb + " ADD");
fw.SetParameter("lonlat2y", MString("(%lat-") + latb + ")*60");
fw.SetParameter("lonlatuv2dotx", "%u*(0.864/1.852)/Cos(%lat*M_PI/180)");
fw.SetParameter("lonlatuv2doty", "%v*(0.864/1.852)");
fw.SetParameter("dotxdotylonlat2u_rp", "EXCH POP EXCH POP COSD 1.852 0.864 DIV MUL MUL");
fw.SetParameter("dotxdotylonlat2v_rp", "POP POP 1.852 0.864 DIV MUL EXCH POP");
fw.SetParameter("dotxdotylonlat2u", "(1.852/0.864)*Cos(%lat*M_PI/180)*%dotx");
fw.SetParameter("dotxdotylonlat2v", "(1.852/0.864)*%doty");
fw.SetParameter("BeginDate", ds.Time(tindexes.front()).ToTString());
fw.SetParameter("EndDate", ds.Time(tindexes.back()).ToTString());
fw.SetParameter("Timestep", ds.Time(tindexes[1]) - ds.Time(tindexes[0]));
fw.SetParameter("DataID", michlib::UniqueId());
}
for(size_t ilat = 0; ilat < data.Ny(); ilat++)
{
for(size_t ilon = 0; ilon < data.Nx(); ilon++)
{
if(data.IsCoast(ilon, ilat))
{
fw.Write(data.Fillval());
fw.Write(data.Fillval());
}
else
{
fw.Write(data.U(ilon, ilat) * (0.864 / 1.852) / Cos(data.Lat(ilon, ilat) * M_PI / 180.0));
fw.Write(data.V(ilon, ilat) * (0.864 / 1.852));
}
}
}
}
fw.Finalize();
fw.Close();
return "";
};

50
include/actions.h
Unescape View File

@ -14,7 +14,8 @@ enum class ActionID
{ {
INFO, INFO,
TSC, TSC,
UV UV,
GENINTFILE
}; };
template<ActionID id> struct ActionCarrier template<ActionID id> struct ActionCarrier
@ -39,10 +40,14 @@ template<class T> struct IsActionSupportedTest<T, ActionID::UV>
{ {
static constexpr bool ans = ReadPSupported<T> || ReadSupported<T>; static constexpr bool ans = ReadPSupported<T> || ReadSupported<T>;
}; };
template<class T> struct IsActionSupportedTest<T, ActionID::GENINTFILE>
{
static constexpr bool ans = ReadIsGrid<T>;
};
template<class T, ActionID id> constexpr bool IsActionSupported = IsActionSupportedTest<T, id>::ans; template<class T, ActionID id> constexpr bool IsActionSupported = IsActionSupportedTest<T, id>::ans;
using ActionVariants = std::variant<ActionCarrier<ActionID::INFO>, ActionCarrier<ActionID::TSC>, ActionCarrier<ActionID::UV>>; using ActionVariants = std::variant<ActionCarrier<ActionID::INFO>, ActionCarrier<ActionID::TSC>, ActionCarrier<ActionID::UV>, ActionCarrier<ActionID::GENINTFILE>>;
class Action: public ActionVariants class Action: public ActionVariants
{ {
@ -58,6 +63,8 @@ class Action: public ActionVariants
*this = ActionVariants(ActionCarrier<ActionID::TSC>()); *this = ActionVariants(ActionCarrier<ActionID::TSC>());
else if(act == "uv") else if(act == "uv")
*this = ActionVariants(ActionCarrier<ActionID::UV>()); *this = ActionVariants(ActionCarrier<ActionID::UV>());
else if(act == "genintfile")
*this = ActionVariants(ActionCarrier<ActionID::GENINTFILE>());
else else
return "Unknown action: " + act; return "Unknown action: " + act;
return ""; return "";
@ -192,32 +199,37 @@ template<class D> BaseData Read(const D& data, const MString& vname, const BaseP
return BaseData(); return BaseData();
} }
template<class D> UVData<ReadType<D>> ReadUV(const D& data, const BaseParameters* p, const TIndex& tindex) template<class D> UVData<ReadType<D>> ReadUV(const D& data, const BaseParameters* p, size_t ind)
{ {
using RT = ReadType<D>; using RT = ReadType<D>;
using UV = UVData<RT>; using UV = UVData<RT>;
size_t ind; michlib::message("Time: " + data.Time(ind).ToTString());
if(tindex.size() == 1) RT u, v;
if constexpr(ReadPSupported<D>)
{ {
ind = tindex[0]; u = data.Read("u", p, ind);
michlib::message("Time: " + data.Time(ind).ToTString()); v = data.Read("v", p, ind);
RT u, v; }
if constexpr(ReadPSupported<D>) else if constexpr(ReadSupported<D>)
{ {
u = data.Read("u", p, ind); u = data.Read("u", ind);
v = data.Read("v", p, ind); v = data.Read("v", ind);
}
else if constexpr(ReadSupported<D>)
{
u = data.Read("u", ind);
v = data.Read("v", ind);
}
return UV(u, v);
} }
return UV(u, v);
}
template<class D> UVData<ReadType<D>> ReadUV(const D& data, const BaseParameters* p, const TIndex& tindex)
{
using RT = ReadType<D>;
using UV = UVData<RT>;
if(tindex.size() == 1)
return ReadUV(data, p, tindex[0]);
else else
{ {
Averager<UV> out; Averager<UV> out;
bool ok = true; bool ok = true;
size_t ind;
for(size_t i = 0; i < tindex.size(); i++) for(size_t i = 0; i < tindex.size(); i++)
{ {
if(!ok) break; if(!ok) break;

1
include/odm.h
Unescape View File

@ -1,4 +1,5 @@
#pragma once #pragma once
#include "actiongenintfile.h"
#include "actioninfo.h" #include "actioninfo.h"
#include "actiontsc.h" #include "actiontsc.h"
#include "actionuv.h" #include "actionuv.h"

30
include/traits.h
Unescape View File

@ -32,13 +32,6 @@ concept ReadSupported = requires(T t, const MString& vname, size_t i) {
} -> std::convertible_to<real>; } -> std::convertible_to<real>;
}; };
template<class T>
concept ReadUVPSupported = requires(T t, const BaseParameters* ip, size_t i) {
{
t.ReadUV(ip, i).U(0, 0) + t.template ReadUV(ip, i).V(0, 0)
} -> std::convertible_to<real>;
};
template<class T> template<class T>
concept IsUVData = requires(T t) { concept IsUVData = requires(T t) {
{ {
@ -52,6 +45,19 @@ concept IsUVData = requires(T t) {
} -> std::convertible_to<real>; } -> std::convertible_to<real>;
}; };
namespace internal
{
template<class T>
concept HaveXYStep = requires(T t) {
{
t.XStep()
} -> std::convertible_to<real>;
{
t.YStep()
} -> std::convertible_to<real>;
};
}
namespace internal namespace internal
{ {
template<typename...> using void_ = void; template<typename...> using void_ = void;
@ -72,3 +78,13 @@ template<class D> struct GetReadType_<D, false, true>
} // namespace internal } // namespace internal
template<class D> using ReadType = internal::GetReadType_<D, ReadPSupported<D>, ReadSupported<D>>::type; template<class D> using ReadType = internal::GetReadType_<D, ReadPSupported<D>, ReadSupported<D>>::type;
template<class T>
concept ReadIsGrid = requires {
{
ReadType<T>().XStep()
} -> std::convertible_to<real>;
{
ReadType<T>().YStep()
} -> std::convertible_to<real>;
};

17
include/uvdata.h
Unescape View File

@ -1,24 +1,12 @@
#pragma once #pragma once
#include "geohelpers.h" #include "geohelpers.h"
#include "traits.h"
#include <vector> #include <vector>
using michlib::M_PI; using michlib::M_PI;
using michlib::real; using michlib::real;
using michlib::Sqrt; using michlib::Sqrt;
namespace internal
{
template<class T>
concept HaveXYStep = requires(T t) {
{
t.XStep()
} -> std::convertible_to<real>;
{
t.YStep()
} -> std::convertible_to<real>;
};
}
enum class Metric enum class Metric
{ {
EUCLID, EUCLID,
@ -147,6 +135,9 @@ template<class OneVarData> class UVDataDims<OneVarData, true>: public UVDataStor
real Lat([[maybe_unused]] size_t ix, size_t iy) const { return y0 + iy * ystep; } 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 Lon(size_t i) const { return Lon(i % nx, i / nx); }
real Lat(size_t i) const { return Lat(i % nx, i / nx); } real Lat(size_t i) const { return Lat(i % nx, i / nx); }
real XStep() const { return xstep; }
real YStep() const { return ystep; }
}; };
template<class OneVarData> class UVDataDims<OneVarData, false>: public UVDataStorage template<class OneVarData> class UVDataDims<OneVarData, false>: public UVDataStorage

8
src/odm.cpp
Unescape View File

@ -4,7 +4,7 @@ inline void Usage(const MString& arg0)
{ {
message(arg0 + " (<key>=<value>...)"); message(arg0 + " (<key>=<value>...)");
message("Keys are:"); message("Keys are:");
message(" action. What the program should do. May be: info, tsc, uv. Default: info."); message(" action. What the program should do. May be: info, tsc, uv, genintfile. Default: info.");
message(" Keys for action=info. Print some information about dataset."); message(" Keys for action=info. Print some information about dataset.");
message(" source. Required. May be: NEMO, HYCOM, AVISO, AVISOLOCAL, BINFILE"); message(" source. Required. May be: NEMO, HYCOM, AVISO, AVISOLOCAL, BINFILE");
message(" Keys for source=NEMO"); message(" Keys for source=NEMO");
@ -41,6 +41,12 @@ inline void Usage(const MString& arg0)
message(" skipx, skipy. How many poinst skipped in the sparsed grid. Used only if velout parameter is present. Default: 1, output each point"); message(" skipx, skipy. How many poinst skipped in the sparsed grid. Used only if velout parameter is present. Default: 1, output each point");
message(" stpout. Output file for stationary points. If absent, this data not calculated."); message(" stpout. Output file for stationary points. If absent, this data not calculated.");
message(" Keys for specific sources see in the section action=tsc."); message(" Keys for specific sources see in the section action=tsc.");
message(" Keys for action=genintfile. Create file with velocity field for interpolation.");
message(" source. Required. May be: NEMO, HYCOM, AVISO, AVISOLOCAL");
message(" time. Time moment or regular expression. If present, timeb and timee must be absent");
message(" timeb, timee. Time interval. If present, time must be absent");
message(" out. Output file, required.");
message(" Keys for specific sources see in the section action=tsc.");
} }
int main(int argc, char** argv) int main(int argc, char** argv)

Write Preview
Loading…
Cancel
Save