Adapter support in the NCFileW class

include/ncfilew.h

@@ -1,4 +1,5 @@
 #pragma once
+#include "Adapter.h"
 #include "MString.h"
 #include "actiondep.h"
 #include "traits.h"
@@ -386,6 +387,14 @@ class NCFileW: public NCFileWBase
  public:
  static constexpr auto Fill() { return fill; }
 
+ michlib::Error AddTimeData(const std::vector<MDateTime>& times);
+
+ michlib::Error Create(const MString& name, const Adapter& ad, const michlib_internal::ParameterListEx& pars, int compression);
+
+ template<class D>
+ requires (std::is_same_v<D, Data2D> || std::is_same_v<D, Data3D>)
+ michlib::Error WriteVariable(const MString& name, const D& data, size_t it);
+
  template<class D> MString Create(const D& data, const MString& name, int compression)
  {
   if(type != UNKNOWN) return "File already created";

src/ncfilew.cpp

@@ -141,3 +141,227 @@ MString NCFileW::AddVariable(const MString& name, const MString& stname, const M
 
  return "";
 }
+
+michlib::Error NCFileW::AddTimeData(const std::vector<MDateTime>& times)
+{
+ static const MString pref = "NCFileW::AddTimeData";
+
+ static constexpr auto                                            stepunits = std::to_array<uint>({3600 * 24, 3600, 60, 1});
+ static constexpr std::array<const std::string, stepunits.size()> stepnames{"days", "hours", "minutes", "seconds"};
+
+ MDateTime                  refdate;
+ size_t                     unitind(stepunits.size());
+ std::vector<michlib::int4> steps(times.size());
+
+ if(steps.size() == 0) return {pref, "No time data available"};
+ refdate = times[0];
+ unitind = 0;
+
+ for(size_t i = 0; i < steps.size(); i++)
+ {
+  auto delta = times[i] - refdate;
+  while(delta.Seconds() % stepunits[unitind] != 0) unitind++;
+ }
+
+ for(size_t i = 0; i < steps.size(); i++)
+ {
+  auto delta = times[i] - refdate;
+  steps[i]   = michlib::int_cast<decltype(steps)::value_type>(delta.Seconds() / stepunits[unitind]);
+ }
+
+ MString unit = stepnames[unitind].c_str() + (" since " + refdate.ToString());
+
+ AddDim("time", steps.size());
+ AddVar("time", Type2NCType<decltype(steps)::value_type>, "time");
+ SetComp("time", compress);
+ AddAtt("time", "standard_name", "time");
+ AddAtt("time", "long_name", "time");
+ AddAtt("time", "units", unit);
+ WriteVar("time", steps.data());
+
+ if(!*this) return {pref, MString("Can't add time data to the netcdf file: ") + ErrMessage()};
+
+ return {};
+}
+
+michlib::Error NCFileW::Create(const MString& name, const Adapter& ad, const michlib_internal::ParameterListEx& pars, int compression)
+{
+ static const MString pref = "NCFileW::Create";
+ if(Ok()) return {pref, "File already opened"};
+
+ const float node_offset = 0.0;
+
+ Open(name);
+ if(!*this) return {pref, "Can't create netcdf file " + name + ": " + ErrMessage()};
+
+ MString rets;
+
+ rets = AddAtt("node_offset", node_offset);
+ if(rets.Exist()) return {pref, rets};
+
+ rets = AddAtts(pars);
+ if(rets.Exist()) return {pref, rets};
+
+ compress = compression;
+
+ {
+  auto ret = AddTimeData(ad.IndexedTimes());
+  if(!ret) return ret.Add(pref, "Can't write time data");
+ }
+
+ for(const auto& v: ad.Vars())
+ {
+  bool is2D = v.second.Read2D ? true : false;
+  bool is3D = v.second.Read3D ? true : false;
+
+  const auto& vname   = v.first;                     // Variable name
+  const auto& proj    = *v.second.proj;              // Horizontal projection
+  const auto& vert    = *v.second.vert;              // Vertical projection
+  const auto& units   = v.second.info.targetunit;    // Unit
+  const auto& stname  = v.second.info.stname;        // Standard name
+  const auto& lname   = v.second.info.lname;         // Long name
+  const auto& comment = v.second.info.comment;       // Comment
+  const auto& xdname  = v.second.info.rinfo->xdname; // X-dimension name
+  const auto& ydname  = v.second.info.rinfo->ydname; // Y-dimension name
+  const auto& zdname  = v.second.info.rinfo->zdname; // Z-dimension name
+  const auto& fillval = v.second.info.fillval;       // _FillValue
+
+  if((is2D && is3D) || (!is2D && !is3D)) return {pref, "Can't determine type of variable " + vname};
+
+  // Z-dimension
+  if(is3D && vert.IsValid())
+  {
+   // TODO Поддержка неизотропных сеток
+   if(!vert.IsIsotropic()) return {pref, "Non-isotropic vertical grids still unsupported"};
+   // Check if z-dimension already exist
+   int  did;
+   auto st = nc_inq_dimid(*this, zdname.Buf(), &did);
+   if(st != NC_NOERR) // Add zdata
+   {
+    AddDim(zdname, vert.Nz());
+    AddVar(zdname, Type2NCType<float>, zdname);
+    SetComp(zdname, compress);
+    AddAtt(zdname, "standard_name", "depth");
+    AddAtt(zdname, "long_name", "Depth");
+    AddAtt(zdname, "units", "m");
+    std::vector<float> data(vert.Nz());
+    for(size_t i = 0; i < vert.Nz(); i++) data[i] = static_cast<float>(vert.Depth(i));
+    WriteVar(zdname, data.data());
+    if(!*this) return {pref, MString("Can't add depth data to the netcdf file: ") + ErrMessage()};
+   }
+   else // Check compatibility with zdata in file
+   {
+    size_t len;
+    st = nc_inq_dimlen(*this, did, &len);
+    if(st != NC_NOERR) return {pref, "Can't get size of depth dimension"};
+    if(len != vert.Nz()) return {pref, "Length of depth dimension does'nt corrrespond to number of vertical layers of variable " + vname};
+   }
+  }
+
+  // X- and Y-dimensions
+  // TODO Поддержка непрямоугольных сеток
+  if(!proj.IsRectangular()) return {pref, "Non-rectangular horizontal grids still unsupported"};
+  {
+   int  lonid, latid;
+   auto stlon = nc_inq_dimid(*this, xdname.Buf(), &lonid);
+   auto stlat = nc_inq_dimid(*this, ydname.Buf(), &latid);
+
+   if(stlon != NC_NOERR) // Add longitudes
+   {
+    AddDim(xdname, proj.Nx());
+    AddVar(xdname, Type2NCType<float>, xdname);
+    SetComp(xdname, compress);
+    AddAtt(xdname, "standard_name", "longitude");
+    AddAtt(xdname, "long_name", "Longitude");
+    AddAtt(xdname, "units", "degrees_east");
+    std::vector<float> data(proj.Nx());
+    for(size_t i = 0; i < proj.Nx(); i++) data[i] = static_cast<float>(proj.Lon(i));
+    WriteVar(xdname, data.data());
+    if(!*this) return {pref, MString("Can't add longitude data to the netcdf file: ") + ErrMessage()};
+   }
+   else // Check compatibility with longitudes in file
+   {
+    size_t len;
+    auto   st = nc_inq_dimlen(*this, lonid, &len);
+    if(st != NC_NOERR) return {pref, "Can't get size of longitude dimension"};
+    if(len != proj.Nx()) return {pref, "Length of longitude dimension does'nt corrrespond to number of longitudes for variable " + vname};
+   }
+
+   if(stlat != NC_NOERR) // Add latitudes
+   {
+    AddDim(ydname, proj.Ny());
+    AddVar(ydname, Type2NCType<float>, ydname);
+    SetComp(ydname, compress);
+    AddAtt(ydname, "standard_name", "latitude");
+    AddAtt(ydname, "long_name", "Latitude");
+    AddAtt(ydname, "units", "degrees_north");
+    std::vector<float> data(proj.Ny());
+    for(size_t i = 0; i < proj.Ny(); i++) data[i] = static_cast<float>(proj.Lat(i));
+    WriteVar(ydname, data.data());
+    if(!*this) return {pref, MString("Can't add latitude data to the netcdf file: ") + ErrMessage()};
+   }
+   else // Check compatibility with longitudes in file
+   {
+    size_t len;
+    auto   st = nc_inq_dimlen(*this, latid, &len);
+    if(st != NC_NOERR) return {pref, "Can't get size of latitude dimension"};
+    if(len != proj.Ny()) return {pref, "Length of latitude dimension does'nt corrrespond to number of latitudes for variable " + vname};
+   }
+  }
+
+  // Variable
+  if(is2D) AddVar(vname, Type2NCType<float>, "time", ydname, xdname);
+  if(is3D) AddVar(vname, Type2NCType<float>, "time", zdname, ydname, xdname);
+  SetComp(vname, compress);
+  if(stname.Exist()) AddAtt(vname, "standard_name", stname);
+  if(lname.Exist()) AddAtt(vname, "long_name", lname);
+  if(units.Exist()) AddAtt(vname, "units", units);
+  if(comment.Exist()) AddAtt(vname, "comment", comment);
+  AddAtt(vname, "_FillValue", static_cast<float>(fillval));
+  if(!*this) return {pref, MString("Can't add variable " + vname + " to the netcdf file: ") + ErrMessage()};
+ } // end vars for
+
+ return {};
+}
+
+template<class D>
+requires(std::is_same_v<D, Data2D> || std::is_same_v<D, Data3D>)
+michlib::Error NCFileW::WriteVariable(const MString& name, const D& data, size_t it)
+{
+ static constexpr const bool   is2D = std::is_same_v<D, Data2D>;
+ static const MString          pref = "NCFileW::WriteVariable";
+ static constexpr const size_t ND   = is2D ? 3 : 4;
+
+ // Check variable presence
+ int varid;
+ {
+  auto ret = nc_inq_varid(*this, name.Buf(), &varid);
+  if(ret != NC_NOERR) return {pref, "Can't find variable " + name};
+ }
+
+ // Check is variable 2D
+ {
+  int  ndim;
+  auto ret = nc_inq_var(*this, varid, nullptr, nullptr, &ndim, nullptr, nullptr);
+  if(ret != NC_NOERR) return {pref, "Can't get number of dimensions for the variable " + name};
+  if(ndim != ND) return {pref, "Number of dimensions for the variable " + name + " is " + ndim + ", but must be " + ND};
+ }
+
+ // Check dimensions?
+
+ // Write
+ const double* pdata;
+ if constexpr(is2D)
+  pdata = &data(0, 0);
+ else
+  pdata = &data(0, 0, 0);
+
+ //WriteVar(name, it, michlib::pointer_cast<const float*>(pdata));
+ WriteVar(name, it, pdata);
+ if(!*this) return {pref, MString("Failed to write variable " + name + ": ") + ErrMessage()};
+
+ return {};
+}
+
+template michlib::Error NCFileW::WriteVariable<Data2D>(const MString& name, const Data2D& data, size_t it);
+template michlib::Error NCFileW::WriteVariable<Data3D>(const MString& name, const Data3D& data, size_t it);