You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.

118 lines
3.8 KiB

#pragma once
#include "BFileW.h"
#include "actiondep.h"
#include <memory>
using michlib::BFileW;
using michlib::Cos;
using michlib::M_PI;
ADD_ACTION(GenIntFile, genintfile, ReadIsGrid<Source>);
template<class D> MString ActionGenIntFile::DoAction(const CLArgs& args, D& ds)
{
auto [reg, regerr] = GetRegion<D>(args);
if(regerr.Exist()) return regerr;
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"));
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>)
{
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.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 "";
};