Line data    Source code

       1             : /******************************************************************************
       2             :  *
       3             :  * Project:  netCDF read/write Driver
       4             :  * Purpose:  GDAL bindings over netCDF library.
       5             :  * Author:   Winor Chen <wchen329 at wisc.edu>
       6             :  *
       7             :  ******************************************************************************
       8             :  * Copyright (c) 2019, Winor Chen <wchen329 at wisc.edu>
       9             :  *
      10             :  * SPDX-License-Identifier: MIT
      11             :  ****************************************************************************/
      12             : #include "netcdfsg.h"
      13             : #include "netcdfdataset.h"
      14             : #include "ogr_core.h"
      15             : 
      16             : namespace nccfdriver
      17             : {
      18          67 : OGRwkbGeometryType RawToOGR(geom_t type, int axis_count)
      19             : {
      20          67 :     OGRwkbGeometryType ret = wkbNone;
      21             : 
      22          67 :     switch (type)
      23             :     {
      24           0 :         case NONE:
      25           0 :             break;
      26           4 :         case LINE:
      27           6 :             ret = axis_count == 2   ? wkbLineString
      28           2 :                   : axis_count == 3 ? wkbLineString25D
      29             :                                     : wkbNone;
      30           4 :             break;
      31           6 :         case MULTILINE:
      32           8 :             ret = axis_count == 2   ? wkbMultiLineString
      33           2 :                   : axis_count == 3 ? wkbMultiLineString25D
      34             :                                     : wkbNone;
      35           6 :             break;
      36           8 :         case POLYGON:
      37          11 :             ret = axis_count == 2   ? wkbPolygon
      38           3 :                   : axis_count == 3 ? wkbPolygon25D
      39             :                                     : wkbNone;
      40           8 :             break;
      41          30 :         case MULTIPOLYGON:
      42          36 :             ret = axis_count == 2   ? wkbMultiPolygon
      43           6 :                   : axis_count == 3 ? wkbMultiPolygon25D
      44             :                                     : wkbNone;
      45          30 :             break;
      46          12 :         case POINT:
      47          21 :             ret = axis_count == 2   ? wkbPoint
      48           9 :                   : axis_count == 3 ? wkbPoint25D
      49             :                                     : wkbNone;
      50          12 :             break;
      51           6 :         case MULTIPOINT:
      52           8 :             ret = axis_count == 2   ? wkbMultiPoint
      53           2 :                   : axis_count == 3 ? wkbMultiPoint25D
      54             :                                     : wkbNone;
      55           6 :             break;
      56           1 :         case UNSUPPORTED:
      57           1 :             break;
      58             :     }
      59             : 
      60          67 :     return ret;
      61             : }
      62             : 
      63         620 : geom_t OGRtoRaw(OGRwkbGeometryType type)
      64             : {
      65         620 :     geom_t ret = NONE;
      66         620 :     auto eFlatType = wkbFlatten(type);
      67             : 
      68         620 :     if (eFlatType == wkbPoint)
      69             :     {
      70         102 :         ret = POINT;
      71             :     }
      72             : 
      73         518 :     else if (eFlatType == wkbLineString)
      74             :     {
      75          12 :         ret = LINE;
      76             :     }
      77             : 
      78         506 :     else if (eFlatType == wkbPolygon)
      79             :     {
      80          42 :         ret = POLYGON;
      81             :     }
      82             : 
      83         464 :     else if (eFlatType == wkbMultiPoint)
      84             :     {
      85          19 :         ret = MULTIPOINT;
      86             :     }
      87             : 
      88         445 :     else if (eFlatType == wkbMultiLineString)
      89             :     {
      90          16 :         ret = MULTILINE;
      91             :     }
      92             : 
      93         429 :     else if (eFlatType == wkbMultiPolygon)
      94             :     {
      95         412 :         ret = MULTIPOLYGON;
      96             :     }
      97             : 
      98             :     // if the feature type isn't NONE potentially give a warning about measures
      99         620 :     if (ret != NONE && wkbHasM(type))
     100             :     {
     101           0 :         CPLError(CE_Warning, CPLE_NotSupported,
     102             :                  "A partially supported measured feature type was detected. X, "
     103             :                  "Y, Z Geometry will be preserved but the measure axis and "
     104             :                  "related information will be removed.");
     105             :     }
     106             : 
     107         620 :     return ret;
     108             : }
     109             : 
     110          42 : bool OGRHasZandSupported(OGRwkbGeometryType type)
     111             : {
     112          34 :     return type == wkbPoint25D || type == wkbLineString25D ||
     113          31 :            type == wkbPolygon25D || type == wkbMultiPoint25D ||
     114          76 :            type == wkbMultiLineString25D || type == wkbMultiPolygon25D;
     115             : }
     116             : 
     117             : }  // namespace nccfdriver
     118             : 
     119          75 : bool netCDFDataset::DetectAndFillSGLayers(int ncid)
     120             : {
     121             :     // Discover simple geometry variables
     122             :     int var_count;
     123          75 :     nc_inq_nvars(ncid, &var_count);
     124          75 :     std::set<int> vidList;
     125             : 
     126          75 :     nccfdriver::scanForGeometryContainers(ncid, vidList);
     127             : 
     128          75 :     if (!vidList.empty())
     129             :     {
     130         140 :         for (auto vid : vidList)
     131             :         {
     132             :             try
     133             :             {
     134          73 :                 LoadSGVarIntoLayer(ncid, vid);
     135             :             }
     136             : 
     137          14 :             catch (nccfdriver::SG_Exception &e)
     138             :             {
     139           7 :                 CPLError(CE_Warning, CPLE_AppDefined,
     140             :                          "Translation of a simple geometry layer has been "
     141             :                          "terminated prematurely due to an error.\n%s",
     142           7 :                          e.get_err_msg());
     143             :             }
     144             :         }
     145             :     }
     146             : 
     147         150 :     return !vidList.empty();
     148             : }
     149             : 
     150          73 : CPLErr netCDFDataset::LoadSGVarIntoLayer(int ncid, int nc_basevarId)
     151             : {
     152             :     std::shared_ptr<nccfdriver::SGeometry_Reader> sg(
     153         140 :         new nccfdriver::SGeometry_Reader(ncid, nc_basevarId));
     154          67 :     int cont_id = sg->getContainerId();
     155         134 :     nccfdriver::SGeometry_PropertyScanner pr(ncid, cont_id);
     156             :     OGRwkbGeometryType owgt =
     157          67 :         nccfdriver::RawToOGR(sg->getGeometryType(), sg->get_axisCount());
     158             : 
     159         134 :     std::string return_gm = "";
     160             : 
     161          67 :     if (sg->getGridMappingVarID() != nccfdriver::INVALID_VAR_ID)
     162          86 :         SetProjectionFromVar(ncid, nc_basevarId, true,
     163          43 :                              sg->getGridMappingName().c_str(), &return_gm,
     164             :                              sg.get(), /*paosRemovedMDItems=*/nullptr);
     165             : 
     166             :     // Geometry Type invalid, avoid further processing
     167          67 :     if (owgt == wkbNone)
     168             :     {
     169           1 :         throw nccfdriver::SG_Exception_BadFeature();
     170             :     }
     171             : 
     172             :     char baseName[NC_MAX_CHAR + 1];
     173          66 :     memset(baseName, 0, NC_MAX_CHAR + 1);
     174          66 :     nc_inq_varname(ncid, nc_basevarId, baseName);
     175             : 
     176          66 :     OGRSpatialReference *poSRS = nullptr;
     177          66 :     if (return_gm != "")
     178             :     {
     179          43 :         poSRS = new OGRSpatialReference();
     180          43 :         if (poSRS->importFromWkt(return_gm.c_str()) != OGRERR_NONE)
     181             :         {
     182           0 :             delete poSRS;
     183           0 :             throw nccfdriver::SG_Exception_General_Malformed("SRS settings");
     184             :         }
     185             : 
     186          43 :         poSRS->SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
     187             :     }
     188             : 
     189             :     std::shared_ptr<netCDFLayer> poL(
     190         132 :         new netCDFLayer(this, ncid, baseName, owgt, poSRS));
     191             : 
     192          66 :     if (poSRS != nullptr)
     193             :     {
     194          43 :         poSRS->Release();
     195             :     }
     196             : 
     197          66 :     poL->EnableSGBypass();
     198          66 :     OGRFeatureDefn *defn = poL->GetLayerDefn();
     199          66 :     defn->SetGeomType(owgt);
     200             : 
     201             :     // Add properties
     202          66 :     std::vector<int> props = pr.ids();
     203         174 :     for (size_t itr = 0; itr < props.size(); itr++)
     204             :     {
     205         108 :         poL->AddField(props[itr]);
     206             :     }
     207             : 
     208             :     // Set simple geometry object
     209          66 :     poL->SetSGeometryRepresentation(sg);
     210             : 
     211             :     // Create layer
     212          66 :     papoLayers.push_back(poL);
     213             : 
     214         132 :     return CE_None;
     215             : }
     216             : 
     217             : /* Creates and fills any needed variables that haven't already been created
     218             :  */
     219         314 : bool netCDFDataset::SGCommitPendingTransaction()
     220             : {
     221             :     try
     222             :     {
     223         314 :         if (bSGSupport)
     224             :         {
     225             :             // Go through all the layers and resize dimensions accordingly
     226          85 :             for (size_t layerInd = 0; layerInd < papoLayers.size(); layerInd++)
     227             :             {
     228             :                 auto poLayer =
     229          45 :                     dynamic_cast<netCDFLayer *>(papoLayers[layerInd].get());
     230          45 :                 if (!poLayer)
     231           0 :                     continue;
     232             :                 nccfdriver::ncLayer_SG_Metadata &layerMD =
     233          45 :                     poLayer->getLayerSGMetadata();
     234          45 :                 nccfdriver::geom_t wType = layerMD.getWritableType();
     235             : 
     236             :                 // Resize node coordinates
     237          45 :                 int ncoord_did = layerMD.get_node_coord_dimID();
     238          45 :                 if (ncoord_did != nccfdriver::INVALID_DIM_ID)
     239             :                 {
     240          42 :                     vcdf.nc_resize_vdim(
     241             :                         ncoord_did, layerMD.get_next_write_pos_node_coord());
     242             :                 }
     243             : 
     244             :                 // Resize node count (for all except POINT)
     245          45 :                 if (wType != nccfdriver::POINT)
     246             :                 {
     247          35 :                     int ncount_did = layerMD.get_node_count_dimID();
     248          35 :                     if (ncount_did != nccfdriver::INVALID_DIM_ID)
     249             :                     {
     250          32 :                         vcdf.nc_resize_vdim(
     251             :                             ncount_did,
     252             :                             layerMD.get_next_write_pos_node_count());
     253             :                     }
     254             :                 }
     255             : 
     256             :                 // Resize part node count (for MULTILINE, POLYGON, MULTIPOLYGON)
     257          45 :                 if (wType == nccfdriver::MULTILINE ||
     258          36 :                     wType == nccfdriver::POLYGON ||
     259             :                     wType == nccfdriver::MULTIPOLYGON)
     260             :                 {
     261          30 :                     int pnc_did = layerMD.get_pnc_dimID();
     262          30 :                     if (pnc_did != nccfdriver::INVALID_DIM_ID)
     263             :                     {
     264          27 :                         vcdf.nc_resize_vdim(pnc_did,
     265             :                                             layerMD.get_next_write_pos_pnc());
     266             :                     }
     267             :                 }
     268             : 
     269          45 :                 nccfdriver::geom_t geometry_type = layerMD.getWritableType();
     270             : 
     271             :                 /* Delete interior ring stuff if not detected
     272             :                  */
     273             : 
     274          83 :                 if (!layerMD.getInteriorRingDetected() &&
     275          22 :                     (geometry_type == nccfdriver::MULTIPOLYGON ||
     276          83 :                      geometry_type == nccfdriver::POLYGON) &&
     277          20 :                     layerMD.get_containerRealID() != nccfdriver::INVALID_VAR_ID)
     278             :                 {
     279          17 :                     SetDefineMode(true);
     280             : 
     281             :                     int err_code =
     282          17 :                         nc_del_att(cdfid, layerMD.get_containerRealID(),
     283             :                                    CF_SG_INTERIOR_RING);
     284          17 :                     NCDF_ERR(err_code);
     285          17 :                     if (err_code != NC_NOERR)
     286             :                     {
     287           0 :                         std::string frmt = std::string("attribute: ") +
     288           0 :                                            std::string(CF_SG_INTERIOR_RING);
     289             :                         throw nccfdriver::SGWriter_Exception_NCDelFailure(
     290           0 :                             layerMD.get_containerName().c_str(), frmt.c_str());
     291             :                     }
     292             : 
     293             :                     // Invalidate variable writes as well - Interior Ring
     294          17 :                     vcdf.nc_del_vvar(layerMD.get_intring_varID());
     295             : 
     296          17 :                     if (geometry_type == nccfdriver::POLYGON)
     297             :                     {
     298             :                         err_code =
     299           4 :                             nc_del_att(cdfid, layerMD.get_containerRealID(),
     300             :                                        CF_SG_PART_NODE_COUNT);
     301           4 :                         NCDF_ERR(err_code);
     302           4 :                         if (err_code != NC_NOERR)
     303             :                         {
     304             :                             std::string frmt =
     305           0 :                                 std::string("attribute: ") +
     306           0 :                                 std::string(CF_SG_PART_NODE_COUNT);
     307             :                             throw nccfdriver::SGWriter_Exception_NCDelFailure(
     308           0 :                                 layerMD.get_containerName().c_str(),
     309           0 :                                 frmt.c_str());
     310             :                         }
     311             : 
     312             :                         // Invalidate variable writes as well - Part Node Count
     313           4 :                         vcdf.nc_del_vvar(layerMD.get_pnc_varID());
     314             : 
     315             :                         // Invalidate dimension as well - Part Node Count
     316           4 :                         vcdf.nc_del_vdim(layerMD.get_pnc_dimID());
     317             :                     }
     318             : 
     319          17 :                     SetDefineMode(false);
     320             :                 }
     321             :             }
     322             : 
     323          40 :             vcdf.nc_vmap();
     324          40 :             this->FieldScribe.commit_transaction();
     325          40 :             this->GeometryScribe.commit_transaction();
     326             :         }
     327             :     }
     328             : 
     329           0 :     catch (nccfdriver::SG_Exception &sge)
     330             :     {
     331           0 :         CPLError(CE_Fatal, CPLE_FileIO,
     332             :                  "An error occurred while writing the target netCDF File. "
     333             :                  "Translation will be terminated.\n%s",
     334           0 :                  sge.get_err_msg());
     335           0 :         return false;
     336             :     }
     337         314 :     return true;
     338             : }
     339             : 
     340         180 : void netCDFDataset::SGLogPendingTransaction()
     341             : {
     342         180 :     GeometryScribe.log_transaction();
     343         180 :     FieldScribe.log_transaction();
     344         180 : }
     345             : 
     346             : /* Takes an index and using the layer geometry builds the equivalent
     347             :  * OGRFeature.
     348             :  */
     349         492 : OGRFeature *netCDFLayer::buildSGeometryFeature(size_t featureInd)
     350             : {
     351             :     OGRGeometry *geometry;
     352             : 
     353         492 :     switch (m_simpleGeometryReader->getGeometryType())
     354             :     {
     355          39 :         case nccfdriver::POINT:
     356          39 :             geometry = new OGRPoint;
     357          39 :             break;
     358          16 :         case nccfdriver::LINE:
     359          16 :             geometry = new OGRLineString;
     360          16 :             break;
     361          13 :         case nccfdriver::POLYGON:
     362          13 :             geometry = new OGRPolygon;
     363          13 :             break;
     364          25 :         case nccfdriver::MULTIPOINT:
     365          25 :             geometry = new OGRMultiPoint;
     366          25 :             break;
     367          18 :         case nccfdriver::MULTILINE:
     368          18 :             geometry = new OGRMultiLineString;
     369          18 :             break;
     370         381 :         case nccfdriver::MULTIPOLYGON:
     371         381 :             geometry = new OGRMultiPolygon;
     372         381 :             break;
     373           0 :         default:
     374           0 :             throw nccfdriver::SG_Exception_BadFeature();
     375             :             break;
     376             :     }
     377             : 
     378         492 :     const auto wkb = m_simpleGeometryReader->serializeToWKB(featureInd);
     379         492 :     geometry->importFromWkb(wkb.data(), wkb.size(), wkbVariantIso);
     380         492 :     geometry->assignSpatialReference(this->GetSpatialRef());
     381             : 
     382         492 :     OGRFeatureDefn *defn = this->GetLayerDefn();
     383         492 :     OGRFeature *feat = new OGRFeature(defn);
     384         492 :     feat->SetGeometryDirectly(geometry);
     385             : 
     386         492 :     int dimId = m_simpleGeometryReader->getInstDim();
     387             : 
     388         492 :     this->FillFeatureFromVar(feat, dimId, featureInd);
     389             : 
     390         492 :     feat->SetFID(featureInd);
     391         984 :     return feat;
     392             : }
     393             : 
     394        1932 : std::string netCDFDataset::generateLogName()
     395             : {
     396        1932 :     return std::string(CPLGenerateTempFilename(nullptr));
     397             : }