1 #include <cpl_string.h>
2 #include "ogrsf_frmts.h"
3 
4 // R headers moved outside extern "C" 070808 RSB re. note from BDR
5 // #ifdef __cplusplus
6 // extern "C" {
7 // #endif
8 
9 #include <Rdefines.h>
10 #include <R.h>
11 #include "rgdal.h"
12 
13 #ifdef __cplusplus
14 extern "C" {
15 #endif
16 
17 // RSB 081009
18 void wrtDF(int, int, SEXP, SEXP, SEXP, OGRFeature*, int);
19 
OGR_write(SEXP inp)20 SEXP OGR_write(SEXP inp)
21 {
22 
23 //  SEXP inp is an input list built in ogr_write() and documented
24 //  in code there
25 
26 // poFeature->SetFID((long) INTEGER_POINTER(VECTOR_ELT(inp, 12))[i]) 130502
27 
28 #ifdef GDALV2
29     GDALDriver *poDriver;
30     GDALDataset *poDS;
31 #else
32     OGRSFDriver *poDriver;
33     OGRDataSource *poDS;
34 #endif
35     OGRLayer *poLayer;
36     OGRFeatureDefn *poDefn;
37     char **papszCreateOptions = NULL;
38     char **papszCreateOptionsLayer = NULL;
39     SEXP ans, wkbtype_attr, comms, ofld_nms;
40     int verbose = INTEGER_POINTER(getAttrib(VECTOR_ELT(inp, 5),
41         install("verbose")))[0];
42     int pc=0, i, j, k, is_shpfile, //shp_edge_case_fix,
43         dumpSRS;
44     OGRSpatialReference *poSRS = new OGRSpatialReference;//FIXME VG
45 
46 
47     PROTECT(ans = NEW_CHARACTER(1)); pc++;
48     PROTECT(wkbtype_attr = NEW_INTEGER(1)); pc++;
49     is_shpfile = INTEGER_POINTER(getAttrib(VECTOR_ELT(inp, 3),
50       install("is_shpfile")))[0];
51 //    shp_edge_case_fix = INTEGER_POINTER(getAttrib(VECTOR_ELT(inp, 3),
52 //      install("shp_edge_case_fix")))[0];
53     dumpSRS = INTEGER_POINTER(getAttrib(VECTOR_ELT(inp, 11),
54       install("dumpSRS")))[0];
55 
56     installErrorHandler();
57 #ifdef GDALV2
58     poDriver = GetGDALDriverManager()->GetDriverByName(CHAR(STRING_ELT(VECTOR_ELT(inp, 3), 0)) );
59 #else
60     poDriver = OGRSFDriverRegistrar::GetRegistrar()->GetDriverByName(
61                 CHAR(STRING_ELT(VECTOR_ELT(inp, 3), 0)) );
62 #endif
63     uninstallErrorHandlerAndTriggerError();
64     if( poDriver == NULL )
65     {
66         error("Driver not available");
67     }
68 //  retrieve and set options:
69 //  papszCreateOptions: a StringList of name=value options.
70 //  Options are driver specific.
71 
72     SEXP sOpts = VECTOR_ELT(inp, 9);
73 
74     installErrorHandler();
75     for (i=0; i < length(sOpts); i++) papszCreateOptions = CSLAddString(
76         papszCreateOptions, CHAR(STRING_ELT(sOpts, i)) );
77     uninstallErrorHandlerAndTriggerError();
78 #ifdef RGDALDEBUG
79     installErrorHandler();
80     for (i=0; i < CSLCount(papszCreateOptions); i++)
81         Rprintf("option %d: %s\n", i, CSLGetField(papszCreateOptions, i));
82     uninstallErrorHandlerAndTriggerError();
83 #endif
84 
85     installErrorHandler();
86 #ifdef GDALV2
87     poDS = poDriver->Create( CHAR(STRING_ELT(VECTOR_ELT(inp,
88         1), 0)), 0, 0, 0, GDT_Unknown, papszCreateOptions );
89 #else
90     poDS = poDriver->CreateDataSource( CHAR(STRING_ELT(VECTOR_ELT(inp,
91         1), 0)), papszCreateOptions );
92 #endif
93     uninstallErrorHandlerAndTriggerError();
94     if( poDS == NULL )
95     {
96         installErrorHandler();
97 #ifdef GDALV2
98         GDALClose( poDS );
99 #else
100         OGRDataSource::DestroyDataSource( poDS );
101 #endif
102         CSLDestroy(papszCreateOptions);
103         uninstallErrorHandlerAndTriggerError();
104         error( "Creation of output file failed" );
105     }
106     installErrorHandler();
107     CSLDestroy(papszCreateOptions);
108     uninstallErrorHandlerAndTriggerError();
109 
110 //  define layer characteristics
111 
112     SEXP obj = VECTOR_ELT(inp, 0);
113     int nobs = INTEGER_POINTER(VECTOR_ELT(inp, 4))[0];
114 // const added 070604 RSB
115     const char *cl = CHAR(asChar(getAttrib(obj, R_ClassSymbol)));
116     OGRwkbGeometryType wkbtype = wkbUnknown;
117 
118     if (!strcmp(cl, "SpatialPointsDataFrame")) {
119         wkbtype = wkbPoint;
120         if (verbose) Rprintf("Object initially classed as: wkbPoint\n");
121     } else if (!strcmp(cl, "SpatialLinesDataFrame")) {
122         wkbtype = wkbLineString;
123         if (verbose) Rprintf("Object initially classed as: wkbLineString\n");
124     } else if (!strcmp(cl, "SpatialPolygonsDataFrame")) {
125         wkbtype = wkbPolygon;
126         if (verbose) Rprintf("Object initially classed as: wkbPolygon\n");
127     }
128 
129     SET_STRING_ELT(ans, 0, COPY_TO_USER_STRING(cl));
130 
131 //  check and if necessary set multiple geometries per data frame row
132 //  for line and polygon objects; multi-points not admitted
133 
134     if (wkbtype == wkbLineString) {
135 
136         SEXP lns = GET_SLOT(obj, install("lines"));
137         if (length(lns) != nobs) {
138             installErrorHandler();
139 #ifdef GDALV2
140             GDALClose( poDS );
141 #else
142             OGRDataSource::DestroyDataSource( poDS );
143 #endif
144             uninstallErrorHandlerAndTriggerError();
145             error("number of objects mismatch");
146         }
147         int multi=0, Lns_l;
148 	for (i=0; i<nobs; i++) {
149             Lns_l = length(GET_SLOT(VECTOR_ELT(lns, i), install("Lines")));
150             if (Lns_l > 1) multi=1;
151 	}
152         if (multi > 0) {
153             wkbtype = wkbMultiLineString;
154             if (verbose) Rprintf("Object reclassed as: wkbMultiLineString\n");
155         }
156     }
157 
158     if (wkbtype == wkbPolygon) {
159         SEXP pls = GET_SLOT(obj, install("polygons"));
160         if (length(pls) != nobs) {
161             installErrorHandler();
162 #ifdef GDALV2
163             GDALClose( poDS );
164 #else
165             OGRDataSource::DestroyDataSource( poDS );
166 #endif
167             uninstallErrorHandlerAndTriggerError();
168             error("number of objects mismatch");
169         }
170         int multi=0, Pls_l, icomms=0;
171 	for (i=0; i<nobs; i++) {
172             comms = SP_PREFIX(comment2comm)(VECTOR_ELT(pls, i));
173             if (comms == R_NilValue) {
174                 Pls_l = length(GET_SLOT(VECTOR_ELT(pls, i),
175                     install("Polygons")));
176                 if (Pls_l > 1) {
177                     multi=1;
178                     break;
179                 }
180             } else {
181                 icomms = 1;
182                 if (length(comms) > 1) {
183                     multi=1;
184                     break;
185                 }
186             }
187 	}
188         if (verbose) {
189             if (icomms == 0) Rprintf("No SFS comments in Polygons objects\n");
190             else Rprintf("SFS comments in Polygons objects\n");
191         }
192         if (multi > 0) {
193             wkbtype = wkbMultiPolygon;
194             if (verbose) Rprintf("Object reclassed as: wkbMultiPolygon\n");
195         }
196     }
197 
198 //  retrieve and set spatial reference system
199 //  retrieve and set options:
200 //  papszCreateOptions: a StringList of name=value options.
201 //  Options are driver specific.
202 
203     SEXP sxpOpts = VECTOR_ELT(inp, 10);
204 
205     installErrorHandler();
206     for (i=0; i < length(sxpOpts); i++) papszCreateOptionsLayer = CSLAddString(
207         papszCreateOptionsLayer, CHAR(STRING_ELT(sxpOpts, i)) );
208     uninstallErrorHandlerAndTriggerError();
209 #ifdef RGDALDEBUG
210     installErrorHandler();
211     for (i=0; i < CSLCount(papszCreateOptionsLayer); i++)
212         Rprintf("option %d: %s\n", i, CSLGetField(papszCreateOptionsLayer, i));
213     uninstallErrorHandlerAndTriggerError();
214 #endif
215 
216     SEXP p4s = GET_SLOT(obj, install("proj4string"));
217     SEXP comment = getAttrib(p4s, install("comment"));
218 
219     if (!isNull(comment)) {
220 //        Rprintf("CRS comment: %s\n", CHAR(STRING_ELT(comment, 0)));
221         installErrorHandler();
222 #if GDAL_VERSION_MAJOR == 1 || ( GDAL_VERSION_MAJOR == 2 && GDAL_VERSION_MINOR <= 2 ) // thanks to Even Roualt https://github.com/OSGeo/gdal/issues/681
223 //#if GDAL_VERSION_MAJOR <= 2 && GDAL_VERSION_MINOR <= 2
224         if (poSRS->importFromWkt((char **) CHAR(STRING_ELT(comment, 0))) != OGRERR_NONE) {
225 #else
226         if (poSRS->importFromWkt((const char *) CHAR(STRING_ELT(comment, 0))) != OGRERR_NONE) {
227 #endif
228 #ifdef GDALV2
229                 GDALClose( poDS );
230 #else
231                 OGRDataSource::DestroyDataSource( poDS );
232 #endif
233             poSRS->Release();
234             uninstallErrorHandlerAndTriggerError();
235             error("Can't parse WKT2-style parameter string");
236         }
237         uninstallErrorHandlerAndTriggerError();
238 
239 #if GDAL_VERSION_MAJOR >= 3
240         installErrorHandler();
241 //Rprintf("OGR_write input AxisMappingStrategy %d\n", poSRS->GetAxisMappingStrategy());
242         poSRS->SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
243         uninstallErrorHandlerAndTriggerError();
244 //Rprintf("OGR_write output AxisMappingStrategy %d\n", poSRS->GetAxisMappingStrategy());
245 #endif
246 
247         installErrorHandler();
248         if (dumpSRS) poSRS->dumpReadable();
249         uninstallErrorHandlerAndTriggerError();
250 
251         installErrorHandler();
252         poLayer = poDS->CreateLayer( CHAR(STRING_ELT(VECTOR_ELT(inp, 2),
253             0)), poSRS, wkbtype, papszCreateOptionsLayer );
254         uninstallErrorHandlerAndTriggerError();
255 
256     } else {
257 // const added 070604 RSB
258         const char *PROJ4 = CHAR(STRING_ELT(GET_SLOT(p4s, install("projargs")), 0));
259 
260         if (strcmp(PROJ4, "NA")) {
261             OGRSpatialReference* poSRS = new OGRSpatialReference;//FIXME VG
262             installErrorHandler();
263 
264             if (poSRS->importFromProj4(PROJ4) != OGRERR_NONE) { //FIXME VG
265 #ifdef GDALV2
266                 GDALClose( poDS );
267 #else
268                 OGRDataSource::DestroyDataSource( poDS );
269 #endif
270                 poSRS->Release();
271                 uninstallErrorHandlerAndTriggerError();
272 	        error("Can't parse PROJ.4-style parameter string");
273             }
274             uninstallErrorHandlerAndTriggerError();
275 /*            installErrorHandler();
276             if (LOGICAL_POINTER(VECTOR_ELT(inp, 11))[0]) {
277                     poSRS->morphToESRI();
278             }
279             uninstallErrorHandlerAndTriggerError();*/
280 
281             installErrorHandler();
282             if (dumpSRS) poSRS->dumpReadable();
283             uninstallErrorHandlerAndTriggerError();
284 
285             installErrorHandler();
286             poLayer = poDS->CreateLayer( CHAR(STRING_ELT(VECTOR_ELT(inp, 2),
287                 0)), poSRS, wkbtype, papszCreateOptionsLayer );
288             uninstallErrorHandlerAndTriggerError();
289 
290         } else {
291             installErrorHandler();
292             poLayer = poDS->CreateLayer( CHAR(STRING_ELT(VECTOR_ELT(inp, 2),
293                 0)), NULL, wkbtype, papszCreateOptionsLayer );
294             uninstallErrorHandlerAndTriggerError();
295         }
296     }
297     if( poLayer == NULL )
298     {
299         installErrorHandler();
300         delete poSRS;
301 #ifdef GDALV2
302         GDALClose( poDS );
303 #else
304         OGRDataSource::DestroyDataSource( poDS );
305 #endif
306         uninstallErrorHandlerAndTriggerError();
307         error( "Layer creation failed" );
308     }
309 
310     installErrorHandler();
311     CSLDestroy(papszCreateOptionsLayer);
312     uninstallErrorHandlerAndTriggerError();
313 
314     INTEGER_POINTER(wkbtype_attr)[0] = wkbtype;
315     setAttrib(ans, install("wkbtype_attr"), wkbtype_attr);
316 // create fields in layer
317 
318     int nf = INTEGER_POINTER(VECTOR_ELT(inp, 5))[0];
319     SEXP fld_names = VECTOR_ELT(inp, 6);
320     SEXP ogr_ftype = VECTOR_ELT(inp, 7);
321     int OGR_type;
322 
323 
324     for (i=0; i<nf; i++) {
325         OGR_type = INTEGER_POINTER(ogr_ftype)[i];
326         if (OGR_type != 0 && OGR_type != 2 && OGR_type != 4) {
327             Rprintf("%s %d\n", CHAR(STRING_ELT(fld_names, i)),
328                 (OGRFieldType) OGR_type);
329             installErrorHandler();
330 #ifdef GDALV2
331             GDALClose( poDS );
332 #else
333             OGRDataSource::DestroyDataSource( poDS );
334 #endif
335             uninstallErrorHandlerAndTriggerError();
336             error( "Unknown field type" );
337         }
338         installErrorHandler();
339         OGRFieldDefn oField( CHAR(STRING_ELT(fld_names, i)),
340             (OGRFieldType)  OGR_type);
341 // RSB 081009 FIXME - not working yet, integer flips to real in shapefile
342         if (OGR_type == 0) oField.SetPrecision(0);
343         if( poLayer->CreateField( &oField ) != OGRERR_NONE ) {
344             delete poSRS;
345 #ifdef GDALV2
346             GDALClose( poDS );
347 #else
348             OGRDataSource::DestroyDataSource( poDS );
349 #endif
350             uninstallErrorHandlerAndTriggerError();
351             error( "Creating Name field failed" );
352         }
353         uninstallErrorHandlerAndTriggerError();
354     }
355 
356     installErrorHandler();
357     poDefn = poLayer->GetLayerDefn();
358     int nFields =  poDefn->GetFieldCount();
359     uninstallErrorHandlerAndTriggerError();
360 
361     PROTECT(ofld_nms = NEW_CHARACTER(nFields)); pc++;
362 
363     installErrorHandler();
364     for (i=0; i<nFields; i++) {
365       OGRFieldDefn *poField = poDefn->GetFieldDefn(i);
366       SET_STRING_ELT(ofld_nms, i, COPY_TO_USER_STRING(poField->GetNameRef()));
367     }
368     uninstallErrorHandlerAndTriggerError();
369 
370     setAttrib(ans, install("ofld_nms"), ofld_nms);
371 
372     SEXP ldata = VECTOR_ELT(inp, 8);
373 
374 // Point data
375 
376 #ifdef GDALV2
377     int transaction_ds = poDS->TestCapability(ODsCTransactions);
378 //Rprintf("transaction_ds: %s, %d\n", CHAR(STRING_ELT(VECTOR_ELT(inp, 3), 0)) , transaction_ds);
379 #endif
380 
381     if (wkbtype == wkbPoint) {
382         SEXP crds, dim;
383         crds = GET_SLOT(obj, install("coords"));
384         dim = getAttrib(crds, install("dim"));
385         int z=INTEGER_POINTER(dim)[1];
386         if (INTEGER_POINTER(dim)[0] != nobs) {
387             installErrorHandler();
388             delete poSRS;
389 #ifdef GDALV2
390             GDALClose( poDS );
391 #else
392             OGRDataSource::DestroyDataSource( poDS );
393 #endif
394             uninstallErrorHandlerAndTriggerError();
395             error("number of objects mismatch");
396         }
397 
398         if (verbose) Rprintf("Writing %d wkbPoint objects\n", nobs);
399         installErrorHandler();
400 
401 #ifdef GDALV2
402         int start_trans = transaction_ds &&
403           (poDS->StartTransaction() == OGRERR_NONE);
404 #endif
405 
406         for (i=0; i<nobs; i++) {
407             OGRFeature *poFeature;
408             poFeature = new OGRFeature( poLayer->GetLayerDefn() );
409 
410 // RSB 081009
411             wrtDF(i, nf, fld_names, ldata, ogr_ftype, poFeature, is_shpfile);
412 
413             OGRPoint pt;
414             pt.setX( NUMERIC_POINTER(crds)[i] );
415             pt.setY( NUMERIC_POINTER(crds)[i+nobs] );
416             if (z > 2) pt.setZ( NUMERIC_POINTER(crds)[i+(2*nobs)] );
417 
418             poFeature->SetGeometry( &pt );
419 #ifdef GDALV2
420             if(poFeature->SetFID((GIntBig) INTEGER_POINTER(VECTOR_ELT(inp, 12))[i])  != OGRERR_NONE ) {
421                installErrorHandler();
422                 GDALClose( poDS );
423 #else
424             if(poFeature->SetFID((long) INTEGER_POINTER(VECTOR_ELT(inp, 12))[i])  != OGRERR_NONE ) {
425                 OGRDataSource::DestroyDataSource( poDS );
426 #endif
427                 delete poSRS;
428                 uninstallErrorHandlerAndTriggerError();
429                 error( "Failed to set FID" );
430             }
431 
432             if( poLayer->CreateFeature( poFeature ) != OGRERR_NONE ) {
433 #ifdef GDALV2
434                 GDALClose( poDS );
435 #else
436                 OGRDataSource::DestroyDataSource( poDS );
437 #endif
438                 delete poSRS;
439                 uninstallErrorHandlerAndTriggerError();
440                 error( "Failed to create feature" );
441             }
442 
443              OGRFeature::DestroyFeature( poFeature );
444         }
445 
446 #ifdef GDALV2
447         if (start_trans && poDS->CommitTransaction() != OGRERR_NONE) {
448           GDALClose( poDS );
449           delete poSRS;
450           uninstallErrorHandlerAndTriggerError();
451           error("commit transaction failure");
452         }
453 #endif
454 
455         uninstallErrorHandlerAndTriggerError();
456 
457 // Line data
458 
459     } else if (wkbtype == wkbLineString) {
460 
461         SEXP lns = GET_SLOT(obj, install("lines"));
462         if (length(lns) != nobs) {
463             installErrorHandler();
464 #ifdef GDALV2
465             GDALClose( poDS );
466 #else
467             OGRDataSource::DestroyDataSource( poDS );
468 #endif
469             delete poSRS;
470             uninstallErrorHandlerAndTriggerError();
471             error("number of objects mismatch");
472         }
473 
474         if (verbose) Rprintf("Writing %d wkbLineString objects\n", nobs);
475         installErrorHandler();
476 #ifdef GDALV2
477         int start_trans = transaction_ds &&
478           (poDS->StartTransaction() == OGRERR_NONE);
479 #endif
480 
481 	for (i=0; i<nobs; i++) {
482 
483             OGRFeature *poFeature;
484             poFeature = new OGRFeature( poLayer->GetLayerDefn() );
485 // RSB 081009
486             wrtDF(i, nf, fld_names, ldata, ogr_ftype, poFeature, is_shpfile);
487 
488             SEXP crds, dim;
489             crds = GET_SLOT(VECTOR_ELT(GET_SLOT(VECTOR_ELT(lns, i),
490                 install("Lines")), 0), install("coords"));
491             dim = getAttrib(crds, install("dim"));
492             int ncrds = INTEGER_POINTER(dim)[0];
493             OGRLineString OGRln;
494             for (j=0; j<ncrds; j++)
495                 OGRln.setPoint( j, NUMERIC_POINTER(crds)[j],
496                                    NUMERIC_POINTER(crds)[j+ncrds] );
497 
498             if( poFeature->SetGeometry( &OGRln ) != OGRERR_NONE ) {
499                installErrorHandler();
500 #ifdef GDALV2
501                 GDALClose( poDS );
502 #else
503                 OGRDataSource::DestroyDataSource( poDS );
504 #endif
505                 delete poSRS;
506                 uninstallErrorHandlerAndTriggerError();
507                 error( "Failed to set geometry" );
508             }
509 
510 #ifdef GDALV2
511             if(poFeature->SetFID((GIntBig) INTEGER_POINTER(VECTOR_ELT(inp, 12))[i])  != OGRERR_NONE ) {
512                installErrorHandler();
513                 GDALClose( poDS );
514 #else
515             if(poFeature->SetFID((long) INTEGER_POINTER(VECTOR_ELT(inp, 12))[i])  != OGRERR_NONE ) {
516                installErrorHandler();
517                 OGRDataSource::DestroyDataSource( poDS );
518 #endif
519                 delete poSRS;
520                 uninstallErrorHandlerAndTriggerError();
521                 error( "Failed to set FID" );
522             }
523 
524             if( poLayer->CreateFeature( poFeature ) != OGRERR_NONE ) {
525                installErrorHandler();
526 #ifdef GDALV2
527                 GDALClose( poDS );
528 #else
529                 OGRDataSource::DestroyDataSource( poDS );
530 #endif
531                 delete poSRS;
532                 uninstallErrorHandlerAndTriggerError();
533                 error( "Failed to create feature" );
534             }
535 
536              OGRFeature::DestroyFeature( poFeature );
537         }
538 #ifdef GDALV2
539         if (start_trans && poDS->CommitTransaction() != OGRERR_NONE) {
540           GDALClose( poDS );
541           delete poSRS;
542           uninstallErrorHandlerAndTriggerError();
543           error("commit transaction failure");
544         }
545 #endif
546 
547         uninstallErrorHandlerAndTriggerError();
548 
549 // Multi line data
550 
551     } else if (wkbtype == wkbMultiLineString) {
552 
553         SEXP lns = GET_SLOT(obj, install("lines"));
554         if (length(lns) != nobs) {
555             installErrorHandler();
556 #ifdef GDALV2
557             GDALClose( poDS );
558 #else
559             OGRDataSource::DestroyDataSource( poDS );
560 #endif
561             delete poSRS;
562             uninstallErrorHandlerAndTriggerError();
563             error("number of objects mismatch");
564         }
565         SEXP Lns;
566         int Lns_l;
567         if (verbose) Rprintf("Writing %d wkbMultiLineString objects\n", nobs);
568         installErrorHandler();
569 #ifdef GDALV2
570         int start_trans = transaction_ds &&
571           (poDS->StartTransaction() == OGRERR_NONE);
572 #endif
573 
574 	for (i=0; i<nobs; i++) {
575 
576             OGRFeature *poFeature;
577             poFeature = new OGRFeature( poLayer->GetLayerDefn() );
578 // RSB 081009
579             wrtDF(i, nf, fld_names, ldata, ogr_ftype, poFeature, is_shpfile);
580 
581             Lns = GET_SLOT(VECTOR_ELT(lns, i), install("Lines"));
582             Lns_l = length(Lns);
583 
584             OGRMultiLineString OGRlns;
585 
586             for (k=0; k<Lns_l; k++) {
587                 SEXP crds, dim;
588                 crds = GET_SLOT(VECTOR_ELT(GET_SLOT(VECTOR_ELT(lns, i),
589                     install("Lines")), k), install("coords"));
590                 dim = getAttrib(crds, install("dim"));
591                 int ncrds = INTEGER_POINTER(dim)[0];
592 
593                 OGRLineString OGRln;
594 
595                 for (j=0; j<ncrds; j++)
596                     OGRln.setPoint( j, NUMERIC_POINTER(crds)[j],
597                                        NUMERIC_POINTER(crds)[j+ncrds] );
598 
599                 if( OGRlns.addGeometry( &OGRln ) != OGRERR_NONE ) {
600                    installErrorHandler();
601 #ifdef GDALV2
602                     GDALClose( poDS );
603 #else
604                     OGRDataSource::DestroyDataSource( poDS );
605 #endif
606                     delete poSRS;
607                     uninstallErrorHandlerAndTriggerError();
608                     error( "Failed to add line" );
609                 }
610             }
611 
612             if( poFeature->SetGeometry( &OGRlns ) != OGRERR_NONE ) {
613                installErrorHandler();
614 #ifdef GDALV2
615                 GDALClose( poDS );
616 #else
617                 OGRDataSource::DestroyDataSource( poDS );
618 #endif
619                 uninstallErrorHandlerAndTriggerError();
620                 error( "Failed to set geometry" );
621             }
622 
623 #ifdef GDALV2
624             if(poFeature->SetFID((GIntBig) INTEGER_POINTER(VECTOR_ELT(inp, 12))[i])  != OGRERR_NONE ) {
625                installErrorHandler();
626                 GDALClose( poDS );
627 #else
628             if(poFeature->SetFID((long) INTEGER_POINTER(VECTOR_ELT(inp, 12))[i])  != OGRERR_NONE ) {
629                installErrorHandler();
630                 OGRDataSource::DestroyDataSource( poDS );
631 #endif
632                 uninstallErrorHandlerAndTriggerError();
633                 delete poSRS;
634                 error( "Failed to set FID" );
635             }
636 
637             if( poLayer->CreateFeature( poFeature ) != OGRERR_NONE ) {
638                installErrorHandler();
639 #ifdef GDALV2
640                 GDALClose( poDS );
641 #else
642                 OGRDataSource::DestroyDataSource( poDS );
643 #endif
644                 delete poSRS;
645                 uninstallErrorHandlerAndTriggerError();
646                 error( "Failed to create feature" );
647             }
648 
649              OGRFeature::DestroyFeature( poFeature );
650         }
651 #ifdef GDALV2
652         if (start_trans && poDS->CommitTransaction() != OGRERR_NONE) {
653           GDALClose( poDS );
654           delete poSRS;
655           uninstallErrorHandlerAndTriggerError();
656           error("commit transaction failure");
657         }
658 #endif
659 
660         uninstallErrorHandlerAndTriggerError();
661 
662 // Polygon data
663 
664     } else if (wkbtype == wkbPolygon) {
665 
666         SEXP pls = GET_SLOT(obj, install("polygons"));
667         if (length(pls) != nobs) {
668             installErrorHandler();
669 #ifdef GDALV2
670             GDALClose( poDS );
671 #else
672             OGRDataSource::DestroyDataSource( poDS );
673 #endif
674             delete poSRS;
675             uninstallErrorHandlerAndTriggerError();
676             error("number of objects mismatch");
677         }
678         if (verbose) Rprintf("Writing %d wkbPolygon objects\n", nobs);
679         installErrorHandler();
680 #ifdef GDALV2
681         int start_trans = transaction_ds &&
682           (poDS->StartTransaction() == OGRERR_NONE);
683 #endif
684 
685 	for (i=0; i<nobs; i++) {
686 
687 
688             OGRFeature *poFeature;
689             poFeature = new OGRFeature( poLayer->GetLayerDefn() );
690 // RSB 081009
691             wrtDF(i, nf, fld_names, ldata, ogr_ftype, poFeature, is_shpfile);
692 
693             SEXP crds, dim;
694 
695             OGRPolygon OGRply;
696 
697             PROTECT(comms = SP_PREFIX(comment2comm)(VECTOR_ELT(pls, i)));
698             if (comms == R_NilValue) {
699                 crds = GET_SLOT(VECTOR_ELT(GET_SLOT(VECTOR_ELT(pls, i),
700                     install("Polygons")), 0), install("coords"));
701                 dim = getAttrib(crds, install("dim"));
702                 int ncrds = INTEGER_POINTER(dim)[0];
703                 OGRLinearRing OGRlr;
704                 for (j=0; j<ncrds; j++)
705                     OGRlr.setPoint( j, NUMERIC_POINTER(crds)[j],
706                                    NUMERIC_POINTER(crds)[j+ncrds] );
707                 OGRply.addRing( &OGRlr );
708             } else {
709                 for (k=0; k<length(VECTOR_ELT(comms, 0)); k++) {
710                     crds = GET_SLOT(VECTOR_ELT(GET_SLOT(VECTOR_ELT(pls, i),
711                         install("Polygons")), INTEGER_POINTER(VECTOR_ELT(comms,
712                         0))[k]-R_OFFSET), install("coords"));
713                     dim = getAttrib(crds, install("dim"));
714                     int ncrds = INTEGER_POINTER(dim)[0];
715                     OGRLinearRing OGRlr;
716                     for (j=0; j<ncrds; j++)
717                         OGRlr.setPoint( j, NUMERIC_POINTER(crds)[j],
718                                    NUMERIC_POINTER(crds)[j+ncrds] );
719                     OGRply.addRing( &OGRlr ); // first is Ering, others Iring(s)
720                 }
721             }
722             UNPROTECT(1);
723             if( poFeature->SetGeometry( &OGRply ) != OGRERR_NONE ) {
724                installErrorHandler();
725 #ifdef GDALV2
726                 GDALClose( poDS );
727 #else
728                 OGRDataSource::DestroyDataSource( poDS );
729 #endif
730                 delete poSRS;
731                 uninstallErrorHandlerAndTriggerError();
732                 error( "Failed to set geometry" );
733             }
734 
735 #ifdef GDALV2
736             if(poFeature->SetFID((GIntBig) INTEGER_POINTER(VECTOR_ELT(inp, 12))[i])  != OGRERR_NONE ) {
737                installErrorHandler();
738                 GDALClose( poDS );
739 #else
740             if(poFeature->SetFID((long) INTEGER_POINTER(VECTOR_ELT(inp, 12))[i])  != OGRERR_NONE ) {
741                installErrorHandler();
742                 OGRDataSource::DestroyDataSource( poDS );
743 #endif
744                 uninstallErrorHandlerAndTriggerError();
745                 delete poSRS;
746                 error( "Failed to set FID" );
747             }
748 
749             if( poLayer->CreateFeature( poFeature ) != OGRERR_NONE ) {
750                installErrorHandler();
751 #ifdef GDALV2
752                 GDALClose( poDS );
753 #else
754                 OGRDataSource::DestroyDataSource( poDS );
755 #endif
756                delete poSRS;
757                uninstallErrorHandlerAndTriggerError();
758                error( "Failed to create feature" );
759             }
760 
761              OGRFeature::DestroyFeature( poFeature );
762         }
763 #ifdef GDALV2
764         if (start_trans && poDS->CommitTransaction() != OGRERR_NONE) {
765           GDALClose( poDS );
766           delete poSRS;
767           uninstallErrorHandlerAndTriggerError();
768           error("commit transaction failure");
769         }
770 #endif
771 
772         uninstallErrorHandlerAndTriggerError();
773 
774 // Multi polygon data
775 
776     } else if (wkbtype == wkbMultiPolygon) {
777 
778         SEXP pls = GET_SLOT(obj, install("polygons"));
779         if (length(pls) != nobs) {
780             installErrorHandler();
781 #ifdef GDALV2
782             GDALClose( poDS );
783 #else
784             OGRDataSource::DestroyDataSource( poDS );
785 #endif
786             delete poSRS;
787             uninstallErrorHandlerAndTriggerError();
788             error("number of objects mismatch");
789         }
790 
791         if (verbose) Rprintf("Writing %d wkbMultiPolygon objects\n", nobs);
792         installErrorHandler();
793 #ifdef GDALV2
794         int start_trans = transaction_ds &&
795           (poDS->StartTransaction() == OGRERR_NONE);
796 #endif
797 
798 	for (i=0; i<nobs; i++) {
799             OGRFeature *poFeature =  new OGRFeature( poLayer->GetLayerDefn() );
800             PROTECT(comms = SP_PREFIX(comment2comm)(VECTOR_ELT(pls, i)));
801             if (comms == R_NilValue) {
802                 error("Polygons object without ordering comment");
803             } else {
804                 int nExtRings = length(comms);
805                 OGRMultiPolygon poRet = OGRMultiPolygon();
806                 SEXP PLSi = GET_SLOT(VECTOR_ELT(pls, i), install("Polygons"));
807                 for (int iER=0; iER<nExtRings; iER++) {
808                     OGRPolygon OGRply = OGRPolygon();
809                     int nthisiER = length(VECTOR_ELT(comms, iER));
810 
811                     for (k=0; k<nthisiER; k++) {
812                         int thisk = INTEGER_POINTER(VECTOR_ELT(comms,
813                             iER))[k]-R_OFFSET;
814                         SEXP crds, dim;
815                         crds = GET_SLOT(VECTOR_ELT(PLSi, thisk),
816                             install("coords"));
817                         dim = getAttrib(crds, install("dim"));
818                         int ncrds = INTEGER_POINTER(dim)[0];
819                         OGRLinearRing OGRlr;
820                         for (j=0; j<ncrds; j++)
821                            OGRlr.setPoint( j, NUMERIC_POINTER(crds)[j],
822                                    NUMERIC_POINTER(crds)[j+ncrds] );
823                         OGRply.addRing( &OGRlr );
824                     }
825                     poRet.addGeometry(&OGRply);
826                     OGRply.empty();
827                 }
828                 if( poFeature->SetGeometry( &poRet ) != OGRERR_NONE ) {
829 #ifdef GDALV2
830                     GDALClose( poDS );
831 #else
832                     OGRDataSource::DestroyDataSource( poDS );
833 #endif
834                     delete poSRS;
835                     uninstallErrorHandlerAndTriggerError();
836                     error( "Failed to set geometry" );
837                 }
838             }
839             UNPROTECT(1);
840             uninstallErrorHandlerAndTriggerError();
841             wrtDF(i, nf, fld_names, ldata, ogr_ftype, poFeature, is_shpfile);
842 
843 // FIXME
844 #ifdef GDALV2
845             if(poFeature->SetFID((GIntBig) INTEGER_POINTER(VECTOR_ELT(inp, 12))[i])  != OGRERR_NONE ) {
846                 GDALClose( poDS );
847 #else
848             if(poFeature->SetFID((long) INTEGER_POINTER(VECTOR_ELT(inp, 12))[i])  != OGRERR_NONE ) {
849                 OGRDataSource::DestroyDataSource( poDS );
850 #endif
851                delete poSRS;
852                uninstallErrorHandlerAndTriggerError();
853                error( "Failed to set FID" );
854             }
855 
856             if( poLayer->CreateFeature( poFeature ) != OGRERR_NONE ) {
857 
858 #ifdef GDALV2
859                 GDALClose( poDS );
860 #else
861                 OGRDataSource::DestroyDataSource( poDS );
862 #endif
863                delete poSRS;
864                uninstallErrorHandlerAndTriggerError();
865                error( "Failed to create feature" );
866             }
867 
868              OGRFeature::DestroyFeature( poFeature );
869         } // i
870 
871 #ifdef GDALV2
872         if (start_trans && poDS->CommitTransaction() != OGRERR_NONE) {
873           GDALClose( poDS );
874           delete poSRS;
875           uninstallErrorHandlerAndTriggerError();
876           error("commit transaction failure");
877         }
878 #endif
879     uninstallErrorHandlerAndTriggerError();
880     } // multiPolygon
881 
882 
883     installErrorHandler();
884 #ifdef GDALV2
885         GDALClose( poDS );
886 #else
887     OGRDataSource::DestroyDataSource( poDS );
888 #endif
889     delete poSRS;
890     uninstallErrorHandlerAndTriggerError();
891 
892     UNPROTECT(pc);
893     return(ans);
894 }
895 
896 // RSB 081009
897 void wrtDF(int i, int nf, SEXP fld_names, SEXP ldata,
898      SEXP ogr_ftype, OGRFeature* poFeature, int is_shpfile) {
899      int j, OGR_type;
900      for (j=0; j<nf; j++) {
901          installErrorHandler();
902          OGR_type = INTEGER_POINTER(ogr_ftype)[j];
903          if (OGR_type == 2) {
904              if (!ISNA(NUMERIC_POINTER(VECTOR_ELT(ldata, j))[i]))
905                if (is_shpfile) {
906                  poFeature->SetField(j,
907                      NUMERIC_POINTER(VECTOR_ELT(ldata, j))[i] );
908                } else {
909                  poFeature->SetField( CHAR(STRING_ELT(fld_names, j)),
910                      NUMERIC_POINTER(VECTOR_ELT(ldata, j))[i] );
911                }
912 #if ((GDAL_VERSION_MAJOR == 2 && GDAL_VERSION_MINOR >= 2) || GDAL_VERSION_MAJOR > 2)
913              else poFeature->SetFieldNull(j);
914 #else
915              else poFeature->UnsetField(j);
916 #endif
917          } else if (OGR_type == 4) {
918              if (STRING_ELT(VECTOR_ELT(ldata, j), i) != NA_STRING)
919                if (is_shpfile) {
920                  poFeature->SetField(j,
921                      CHAR(STRING_ELT(VECTOR_ELT(ldata, j), i)) );
922                } else {
923                  poFeature->SetField( CHAR(STRING_ELT(fld_names, j)),
924                      CHAR(STRING_ELT(VECTOR_ELT(ldata, j), i)) );
925                }
926 #if ((GDAL_VERSION_MAJOR == 2 && GDAL_VERSION_MINOR >= 2) || GDAL_VERSION_MAJOR > 2)
927              else poFeature->SetFieldNull(j);
928 #else
929              else poFeature->UnsetField(j);
930 #endif
931          } else if (OGR_type == 0) {
932               if (INTEGER_POINTER(VECTOR_ELT(ldata, j))[i] != NA_INTEGER)
933                if (is_shpfile) {
934                   poFeature->SetField(j,
935                       INTEGER_POINTER(VECTOR_ELT(ldata, j))[i] );
936                } else {
937                   poFeature->SetField( CHAR(STRING_ELT(fld_names, j)),
938                       INTEGER_POINTER(VECTOR_ELT(ldata, j))[i] );
939                }
940 #if ((GDAL_VERSION_MAJOR == 2 && GDAL_VERSION_MINOR >= 2) || GDAL_VERSION_MAJOR > 2)
941              else poFeature->SetFieldNull(j);
942 #else
943              else poFeature->UnsetField(j);
944 #endif
945          }
946          uninstallErrorHandlerAndTriggerError();
947      }
948 }
949 
950 /* FIXME return to this later
951  * SEXP ogrOrganizeSpatialPolygons(SEXP obj) {
952         SEXP ans;
953         SEXP pls = GET_SLOT(obj, install("polygons"));
954         SEXP Lns, Pls;
955         int Lns_l, pls_l, i, k, pc=0;
956         pls_l <- length(pls);
957         PROTECT(ans = NEW_LIST(pls_l)); pc++;
958         installErrorHandler();
959 	for (i=0; i<; i++) {
960             Pls = VECTOR_ELT(pls, i);
961             Lns = GET_SLOT(Pls, install("Polygons"));
962             Lns_l = (int) length(Lns);
963             if (Lns_l == 1) {
964 
965                 SET_VECTOR_ELT(ans, i, Pls);
966             } else {
967                 OGRGeometry* poRet = NULL;
968                 OGRPolygon** papoPolygons = new OGRPolygon*[ Lns_l ];
969                 for (k=0; k<Lns_l; k++) {
970                     papoPolygons[k] = new OGRPolygon();
971                     SEXP crds, dim;
972                     crds = GET_SLOT(VECTOR_ELT(GET_SLOT(VECTOR_ELT(pls, i),
973                         install("Polygons")), k), install("coords"));
974                     dim = getAttrib(crds, install("dim"));
975                     int ncrds = INTEGER_POINTER(dim)[0];
976 
977                     OGRLinearRing *OGRlr = new OGRLinearRing;
978 
979                     for (j=0; j<ncrds; j++)
980                         OGRlr->setPoint( j, NUMERIC_POINTER(crds)[j],
981                                     NUMERIC_POINTER(crds)[j+ncrds] );
982 
983                     papoPolygons[k]->addRingDirectly(OGRlr);
984 
985                 } // k
986                 int isValidGeometry;
987                 poRet = OGRGeometryFactory::organizePolygons(
988                     (OGRGeometry**)papoPolygons, Lns_l, &isValidGeometry );
989                 if (!isValidGeometry) {
990                     warning("OGR_write: uncommented multiring Polygons object %d conversion to SFS invalid", i+R_OFFSET);
991                 }
992 
993                 delete[] papoPolygons;
994             }
995          }
996          uninstallErrorHandlerAndTriggerError();
997 
998 
999 
1000 }*/
1001 
1002 #ifdef __cplusplus
1003 }
1004 #endif
1005 
1006