odm/actions/actiontsc.h

#pragma once
#include "BFileW.h"
#include "actiondep.h"
#include "ncfuncs.h"
#include <memory>

using michlib::BFileW;

ADD_ACTION(TSC, tsc, ReadPSupported<Source> || ReadSupported<Source>);

template<class D> MString ActionTSC::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"));

 auto [tindexes, err] = GetTIndexes(ds, args, pars);
 if(err.Exist()) return err;

 if(!args.contains("var")) return "Variable not specified";
 MString vname = args.at("var");
 if(!ds.CheckVar(vname)) return "Variable " + vname + " not exists in this dataset";
 pars.SetParameter("variable", vname);

 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();

 auto data = Read(ds, vname, p, tindexes);
 if(!data) return "Can't read data";
 if(!data.Unit().Exist()) michlib::errmessage("Unknown measurement unit!");

 MString name   = args.contains("out") ? args.at("out") : "out.bin";
 MString outfmt = args.contains("format") ? args.at("format") : (GetExt(name) == "nc" ? "nc" : "bin");

 if(outfmt == "bin")
 {
  BFileW fw;

  fw.Create(name, 3);
  fw.SetColumnName(1, "Longitude");
  fw.SetColumnName(2, "Latitude");
  fw.SetColumnName(3, vname + ", " + (data.Unit().Exist() ? data.Unit() : "unknown"));
  fw.SetParameters(pars);
  for(size_t i = 0; i < data.N(); i++)
  {
   fw.Write(data.Lon(i));
   fw.Write(data.Lat(i));
   fw.Write(data.IsFill(i) ? NAN : data(i));
  }
  fw.Finalize();
  fw.Close();
  return "";
 }

 if(outfmt == "nc" || outfmt == "netcdf")
 {
  int     ret;
  int     ncid, dimid, dimidxy[2];
  int     lonid, latid, vid;
  MString text;
  auto    fill = std::numeric_limits<float>::max();
  float   fval;
  int     compress = 3;

  if(args.contains("compress")) compress = args.at("compress").ToInt();

  ret = nc_create(name.Buf(), NC_CLOBBER | NC_NETCDF4, &ncid);
  if(ret != NC_NOERR) return "Can't create netcdf file: " + name;

  ret = nc_put_att_text(ncid, NC_GLOBAL, "history", args.at("_cmdline").Len() + 1, args.at("_cmdline").Buf());
  if(ret != NC_NOERR) return "Can't write history attribute in the netcdf file";

  if constexpr(ReadIs2DGeoRectArray<D>)
  {
   ret = nc_def_dim(ncid, "longitude", data.Nx(), &dimidxy[0]);
   if(ret != NC_NOERR) return "Can't create x-dimension in the netcdf file";
   ret = nc_def_dim(ncid, "latitude", data.Ny(), &dimidxy[1]);
   if(ret != NC_NOERR) return "Can't create y-dimension in the netcdf file";
   ret = nc_def_var(ncid, "longitude", NC_FLOAT, 1, &dimidxy[0], &lonid);
   if(ret != NC_NOERR) return "Can't create longitude variable in the netcdf file";
   ret = nc_def_var(ncid, "latitude", NC_FLOAT, 1, &dimidxy[1], &latid);
   if(ret != NC_NOERR) return "Can't create latitude variable in the netcdf file";
   ret = nc_def_var(ncid, vname.Buf(), NC_FLOAT, 2, dimidxy, &vid);
   if(ret != NC_NOERR) return "Can't create " + vname + " variable in the netcdf file";
  }
  else
  {
   ret = nc_def_dim(ncid, "i", data.N(), &dimid);
   if(ret != NC_NOERR) return "Can't create dimension in the netcdf file";
   ret = nc_def_var(ncid, "longitude", NC_FLOAT, 1, &dimid, &lonid);
   if(ret != NC_NOERR) return "Can't create longitude variable in the netcdf file";
   ret = nc_def_var(ncid, "latitude", NC_FLOAT, 1, &dimid, &latid);
   if(ret != NC_NOERR) return "Can't create latitude variable in the netcdf file";
   ret = nc_def_var(ncid, vname.Buf(), NC_FLOAT, 1, &dimid, &vid);
   if(ret != NC_NOERR) return "Can't create " + vname + " variable in the netcdf file";
  }

  ret = nc_def_var_deflate(ncid, lonid, 1, 1, compress);
  if(ret != NC_NOERR) return "Can't set deflate parameters for longitude variable in the netcdf file";
  ret = nc_def_var_deflate(ncid, latid, 1, 1, compress);
  if(ret != NC_NOERR) return "Can't set deflate parameters for latitude variable in the netcdf file";
  ret = nc_def_var_deflate(ncid, vid, 1, 1, compress);
  if(ret != NC_NOERR) return "Can't set deflate parameters for " + vname + " variable in the netcdf file";

  text = "longitude";
  ret  = nc_put_att_text(ncid, lonid, "standard_name", text.Len() + 1, text.Buf());
  if(ret != NC_NOERR) return "Can't write standard_name attribute of longitude variable in the netcdf file";
  text = "latitude";
  ret  = nc_put_att_text(ncid, latid, "standard_name", text.Len() + 1, text.Buf());
  if(ret != NC_NOERR) return "Can't write standard_name attribute of latitude variable in the netcdf file";

  text = NCFuncs::Name2StName(vname);
  if(text.Exist())
  {
   ret = nc_put_att_text(ncid, vid, "standard_name", text.Len() + 1, text.Buf());
   if(ret != NC_NOERR) return "Can't write standard_name attribute of " + vname + " variable in the netcdf file";
  }
  else
  {
   text = NCFuncs::Name2LongName(vname);
   ret  = nc_put_att_text(ncid, vid, "long_name", text.Len() + 1, text.Buf());
   if(ret != NC_NOERR) return "Can't write long_name attribute of " + vname + " variable in the netcdf file";
  }

  text = data.Unit();
  if(text.Exist())
  {
   ret = nc_put_att_text(ncid, vid, "units", text.Len() + 1, text.Buf());
   if(ret != NC_NOERR) return "Can't write units attribute of " + vname + " variable in the netcdf file";
  }

  ret = nc_put_att_float(ncid, vid, "_FillValue", NC_FLOAT, 1, &fill);
  if(ret != NC_NOERR) return "Can't write _FillValue attribute of " + vname + " variable in the netcdf file";

  ret = nc_enddef(ncid);
  if(ret != NC_NOERR) return "Can't finish definition of the netcdf file";

  if constexpr(ReadIs2DGeoRectArray<D>)
  {
   size_t i[2];
   for(i[0] = 0; i[0] < data.Nx(); i[0]++)
   {
    fval = data.Ix2Lon(i[0]);
    ret  = nc_put_var1_float(ncid, lonid, &i[0], &fval);
    if(ret != NC_NOERR) return "Can't write longitude variable in the netcdf file: " + MString(i[0]);
   }
   for(i[1] = 0; i[1] < data.Ny(); i[1]++)
   {
    fval = data.Iy2Lat(i[1]);
    ret  = nc_put_var1_float(ncid, latid, &i[1], &fval);
    if(ret != NC_NOERR) return "Can't write latitude variable in the netcdf file: " + MString(i[1]);
   }
   for(i[0] = 0; i[0] < data.Nx(); i[0]++)
    for(i[1] = 0; i[1] < data.Ny(); i[1]++)
    {
     fval = data.IsFill(i[0], i[1]) ? fill : data(i[0], i[1]);
     ret  = nc_put_var1_float(ncid, vid, i, &fval);
     if(ret != NC_NOERR) return "Can't write " + vname + " variable in the netcdf file: " + MString(i[0]) + ", " + i[1];
    }
  }
  else
  {
   for(size_t i = 0; i < data.N(); i++)
   {
    fval = data.Lon(i);
    ret  = nc_put_var1_float(ncid, lonid, &i, &fval);
    if(ret != NC_NOERR) return "Can't write longitude variable in the netcdf file: " + MString(i);
    fval = data.Lat(i);
    ret  = nc_put_var1_float(ncid, latid, &i, &fval);
    if(ret != NC_NOERR) return "Can't write latitude variable in the netcdf file: " + MString(i);
    fval = data.IsFill(i) ? fill : data(i);
    ret  = nc_put_var1_float(ncid, vid, &i, &fval);
    if(ret != NC_NOERR) return "Can't write " + vname + " variable in the netcdf file: " + MString(i);
   }
  }

  ret = nc_close(ncid);
  if(ret != NC_NOERR) return "Can't close netcdf file";

  return "";
 }

 return "Unknown format: " + outfmt;
};