Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: R Format Driver
4 : * Purpose: Read/write R stats package object format.
5 : * Author: Frank Warmerdam, warmerdam@pobox.com
6 : *
7 : ******************************************************************************
8 : * Copyright (c) 2009, Frank Warmerdam <warmerdam@pobox.com>
9 : * Copyright (c) 2009-2010, Even Rouault <even dot rouault at spatialys.com>
10 : *
11 : * SPDX-License-Identifier: MIT
12 : ****************************************************************************/
13 :
14 : #include "cpl_port.h"
15 : #include "rdataset.h"
16 :
17 : #include <cstddef>
18 : #include <cstdlib>
19 : #include <cstring>
20 : #if HAVE_FCNTL_H
21 : #include <fcntl.h>
22 : #endif
23 :
24 : #include <algorithm>
25 : #include <limits>
26 : #include <string>
27 : #include <utility>
28 :
29 : #include "cpl_conv.h"
30 : #include "cpl_error.h"
31 : #include "cpl_progress.h"
32 : #include "cpl_string.h"
33 : #include "cpl_vsi.h"
34 : #include "gdal.h"
35 : #include "gdal_frmts.h"
36 : #include "gdal_pam.h"
37 : #include "gdal_priv.h"
38 :
39 : // constexpr int R_NILSXP = 0;
40 : constexpr int R_LISTSXP = 2;
41 : constexpr int R_CHARSXP = 9;
42 : constexpr int R_INTSXP = 13;
43 : constexpr int R_REALSXP = 14;
44 : constexpr int R_STRSXP = 16;
45 :
46 : namespace
47 : {
48 :
49 : // TODO(schwehr): Move this to port/? for general use.
50 20 : bool SafeMult(GIntBig a, GIntBig b, GIntBig *result)
51 : {
52 20 : if (a == 0 || b == 0)
53 : {
54 0 : *result = 0;
55 0 : return true;
56 : }
57 :
58 20 : bool result_positive = (a >= 0 && b >= 0) || (a < 0 && b < 0);
59 20 : if (result_positive)
60 : {
61 : // Cannot convert min() to positive.
62 40 : if (a == std::numeric_limits<GIntBig>::min() ||
63 20 : b == std::numeric_limits<GIntBig>::min())
64 : {
65 0 : *result = 0;
66 0 : return false;
67 : }
68 20 : if (a < 0)
69 : {
70 0 : a = -a;
71 0 : b = -b;
72 : }
73 20 : if (a > std::numeric_limits<GIntBig>::max() / b)
74 : {
75 0 : *result = 0;
76 0 : return false;
77 : }
78 20 : *result = a * b;
79 20 : return true;
80 : }
81 :
82 0 : if (b < a)
83 0 : std::swap(a, b);
84 0 : if (a < (std::numeric_limits<GIntBig>::min() + 1) / b)
85 : {
86 0 : *result = 0;
87 0 : return false;
88 : }
89 :
90 0 : *result = a * b;
91 0 : return true;
92 : }
93 :
94 : } // namespace
95 :
96 : /************************************************************************/
97 : /* RRasterBand() */
98 : /************************************************************************/
99 :
100 7 : RRasterBand::RRasterBand(RDataset *poDSIn, int nBandIn,
101 7 : const double *padfMatrixValuesIn)
102 7 : : padfMatrixValues(padfMatrixValuesIn)
103 : {
104 7 : poDS = poDSIn;
105 7 : nBand = nBandIn;
106 :
107 7 : eDataType = GDT_Float64;
108 :
109 7 : nBlockXSize = poDSIn->nRasterXSize;
110 7 : nBlockYSize = 1;
111 7 : }
112 :
113 : /************************************************************************/
114 : /* IReadBlock() */
115 : /************************************************************************/
116 :
117 45 : CPLErr RRasterBand::IReadBlock(int /* nBlockXOff */, int nBlockYOff,
118 : void *pImage)
119 : {
120 45 : memcpy(pImage, padfMatrixValues + nBlockYOff * nBlockXSize,
121 45 : nBlockXSize * 8);
122 45 : return CE_None;
123 : }
124 :
125 : /************************************************************************/
126 : /* ==================================================================== */
127 : /* RDataset() */
128 : /* ==================================================================== */
129 : /************************************************************************/
130 :
131 : /************************************************************************/
132 : /* RDataset() */
133 : /************************************************************************/
134 :
135 10 : RDataset::RDataset()
136 10 : : fp(nullptr), bASCII(FALSE), nStartOfData(0), padfMatrixValues(nullptr)
137 : {
138 10 : }
139 :
140 : /************************************************************************/
141 : /* ~RDataset() */
142 : /************************************************************************/
143 :
144 20 : RDataset::~RDataset()
145 : {
146 10 : FlushCache(true);
147 10 : CPLFree(padfMatrixValues);
148 :
149 10 : if (fp)
150 10 : VSIFCloseL(fp);
151 20 : }
152 :
153 : /************************************************************************/
154 : /* ASCIIFGets() */
155 : /* */
156 : /* Fetch one line from an ASCII source into osLastStringRead. */
157 : /************************************************************************/
158 :
159 1385 : const char *RDataset::ASCIIFGets()
160 :
161 : {
162 1385 : char chNextChar = '\0';
163 :
164 1385 : osLastStringRead.resize(0);
165 :
166 3952 : do
167 : {
168 5337 : chNextChar = '\n';
169 5337 : VSIFReadL(&chNextChar, 1, 1, fp);
170 5337 : if (chNextChar != '\n')
171 3952 : osLastStringRead += chNextChar;
172 5337 : } while (chNextChar != '\n' && chNextChar != '\0');
173 :
174 2770 : return osLastStringRead;
175 : }
176 :
177 : /************************************************************************/
178 : /* ReadInteger() */
179 : /************************************************************************/
180 :
181 190 : int RDataset::ReadInteger()
182 :
183 : {
184 190 : if (bASCII)
185 : {
186 95 : return atoi(ASCIIFGets());
187 : }
188 :
189 95 : GInt32 nValue = 0;
190 :
191 95 : if (VSIFReadL(&nValue, 4, 1, fp) != 1)
192 0 : return -1;
193 95 : CPL_MSBPTR32(&nValue);
194 :
195 95 : return nValue;
196 : }
197 :
198 : /************************************************************************/
199 : /* ReadFloat() */
200 : /************************************************************************/
201 :
202 1280 : double RDataset::ReadFloat()
203 :
204 : {
205 1280 : if (bASCII)
206 : {
207 1280 : return CPLAtof(ASCIIFGets());
208 : }
209 :
210 0 : double dfValue = 0.0;
211 :
212 0 : if (VSIFReadL(&dfValue, 8, 1, fp) != 1)
213 0 : return -1;
214 0 : CPL_MSBPTR64(&dfValue);
215 :
216 0 : return dfValue;
217 : }
218 :
219 : /************************************************************************/
220 : /* ReadString() */
221 : /************************************************************************/
222 :
223 20 : const char *RDataset::ReadString()
224 :
225 : {
226 20 : if (ReadInteger() % 256 != R_CHARSXP)
227 : {
228 0 : osLastStringRead = "";
229 0 : return "";
230 : }
231 :
232 20 : const int nLenSigned = ReadInteger();
233 20 : if (nLenSigned < 0)
234 : {
235 0 : osLastStringRead = "";
236 0 : return "";
237 : }
238 20 : const size_t nLen = static_cast<size_t>(nLenSigned);
239 :
240 20 : char *pachWrkBuf = static_cast<char *>(VSIMalloc(nLen));
241 20 : if (pachWrkBuf == nullptr)
242 : {
243 0 : osLastStringRead = "";
244 0 : return "";
245 : }
246 20 : if (VSIFReadL(pachWrkBuf, 1, nLen, fp) != nLen)
247 : {
248 0 : osLastStringRead = "";
249 0 : CPLFree(pachWrkBuf);
250 0 : return "";
251 : }
252 :
253 20 : if (bASCII)
254 : {
255 : // Suck up newline and any extra junk.
256 10 : ASCIIFGets();
257 : }
258 :
259 20 : osLastStringRead.assign(pachWrkBuf, nLen);
260 20 : CPLFree(pachWrkBuf);
261 :
262 20 : return osLastStringRead;
263 : }
264 :
265 : /************************************************************************/
266 : /* ReadPair() */
267 : /************************************************************************/
268 :
269 30 : bool RDataset::ReadPair(CPLString &osObjName, int &nObjCode)
270 :
271 : {
272 30 : nObjCode = ReadInteger();
273 30 : if (nObjCode == 254)
274 10 : return true;
275 :
276 20 : if ((nObjCode % 256) != R_LISTSXP)
277 : {
278 0 : CPLError(CE_Failure, CPLE_OpenFailed,
279 : "Did not find expected object pair object.");
280 0 : return false;
281 : }
282 :
283 20 : int nPairCount = ReadInteger();
284 20 : if (nPairCount != 1)
285 : {
286 0 : CPLError(CE_Failure, CPLE_OpenFailed,
287 : "Did not find expected pair count of 1.");
288 0 : return false;
289 : }
290 :
291 : // Read the object name.
292 20 : const char *pszName = ReadString();
293 20 : if (pszName == nullptr || pszName[0] == '\0')
294 0 : return false;
295 :
296 20 : osObjName = pszName;
297 :
298 : // Confirm that we have a numeric matrix object.
299 20 : nObjCode = ReadInteger();
300 :
301 20 : return true;
302 : }
303 :
304 : /************************************************************************/
305 : /* Identify() */
306 : /************************************************************************/
307 :
308 52846 : int RDataset::Identify(GDALOpenInfo *poOpenInfo)
309 : {
310 52846 : if (poOpenInfo->nHeaderBytes < 50)
311 48572 : return FALSE;
312 :
313 : // If the extension is .rda and the file type is gzip
314 : // compressed we assume it is a gzipped R binary file.
315 4308 : if (memcmp(poOpenInfo->pabyHeader, "\037\213\b", 3) == 0 &&
316 34 : EQUAL(CPLGetExtension(poOpenInfo->pszFilename), "rda"))
317 6 : return TRUE;
318 :
319 : // Is this an ASCII or XDR binary R file?
320 4268 : if (!STARTS_WITH_CI((const char *)poOpenInfo->pabyHeader, "RDA2\nA\n") &&
321 4258 : !STARTS_WITH_CI((const char *)poOpenInfo->pabyHeader, "RDX2\nX\n"))
322 4254 : return FALSE;
323 :
324 14 : return TRUE;
325 : }
326 :
327 : /************************************************************************/
328 : /* Open() */
329 : /************************************************************************/
330 :
331 10 : GDALDataset *RDataset::Open(GDALOpenInfo *poOpenInfo)
332 : {
333 : #ifdef FUZZING_BUILD_MODE_UNSAFE_FOR_PRODUCTION
334 : if (poOpenInfo->pabyHeader == nullptr)
335 : return nullptr;
336 : #else
337 : // During fuzzing, do not use Identify to reject crazy content.
338 10 : if (!Identify(poOpenInfo))
339 0 : return nullptr;
340 : #endif
341 :
342 : // Confirm the requested access is supported.
343 10 : if (poOpenInfo->eAccess == GA_Update)
344 : {
345 0 : CPLError(CE_Failure, CPLE_NotSupported,
346 : "The R driver does not support update access to existing"
347 : " datasets.");
348 0 : return nullptr;
349 : }
350 :
351 : // Do we need to route the file through the decompression machinery?
352 10 : const bool bCompressed =
353 10 : memcmp(poOpenInfo->pabyHeader, "\037\213\b", 3) == 0;
354 : const CPLString osAdjustedFilename =
355 30 : std::string(bCompressed ? "/vsigzip/" : "") + poOpenInfo->pszFilename;
356 :
357 : // Establish this as a dataset and open the file using VSI*L.
358 20 : auto poDS = std::make_unique<RDataset>();
359 :
360 10 : poDS->fp = VSIFOpenL(osAdjustedFilename, "r");
361 10 : if (poDS->fp == nullptr)
362 : {
363 0 : return nullptr;
364 : }
365 :
366 10 : poDS->bASCII = STARTS_WITH_CI(
367 : reinterpret_cast<char *>(poOpenInfo->pabyHeader), "RDA2\nA\n");
368 :
369 : // Confirm this is a version 2 file.
370 10 : VSIFSeekL(poDS->fp, 7, SEEK_SET);
371 10 : if (poDS->ReadInteger() != R_LISTSXP)
372 : {
373 0 : CPLError(CE_Failure, CPLE_OpenFailed,
374 : "It appears %s is not a version 2 R object file after all!",
375 : poOpenInfo->pszFilename);
376 0 : return nullptr;
377 : }
378 :
379 : // Skip the version values.
380 10 : poDS->ReadInteger();
381 10 : poDS->ReadInteger();
382 :
383 : // Confirm we have a numeric vector object in a pairlist.
384 20 : CPLString osObjName;
385 10 : int nObjCode = 0;
386 :
387 10 : if (!poDS->ReadPair(osObjName, nObjCode))
388 : {
389 0 : return nullptr;
390 : }
391 :
392 10 : if (nObjCode % 256 != R_REALSXP)
393 : {
394 0 : CPLError(CE_Failure, CPLE_OpenFailed,
395 : "Failed to find expected numeric vector object.");
396 0 : return nullptr;
397 : }
398 :
399 10 : poDS->SetMetadataItem("R_OBJECT_NAME", osObjName);
400 :
401 : // Read the count.
402 10 : const int nValueCount = poDS->ReadInteger();
403 10 : if (nValueCount < 0)
404 : {
405 0 : CPLError(CE_Failure, CPLE_AppDefined, "nValueCount < 0: %d",
406 : nValueCount);
407 0 : return nullptr;
408 : }
409 :
410 10 : poDS->nStartOfData = VSIFTellL(poDS->fp);
411 :
412 : // TODO(schwehr): Factor in the size of doubles.
413 : VSIStatBufL stat;
414 : const int dStatSuccess =
415 10 : VSIStatExL(osAdjustedFilename, &stat, VSI_STAT_SIZE_FLAG);
416 20 : if (dStatSuccess != 0 || static_cast<vsi_l_offset>(nValueCount) >
417 10 : stat.st_size - poDS->nStartOfData)
418 : {
419 0 : CPLError(CE_Failure, CPLE_AppDefined,
420 : "Corrupt file. "
421 : "Object claims to be larger than available bytes. "
422 : "%d > " CPL_FRMT_GUIB,
423 0 : nValueCount, stat.st_size - poDS->nStartOfData);
424 0 : return nullptr;
425 : }
426 :
427 : // Read/Skip ahead to attributes.
428 10 : if (poDS->bASCII)
429 : {
430 10 : poDS->padfMatrixValues =
431 5 : static_cast<double *>(VSIMalloc2(nValueCount, sizeof(double)));
432 5 : if (poDS->padfMatrixValues == nullptr)
433 : {
434 0 : CPLError(CE_Failure, CPLE_AppDefined, "Cannot allocate %d doubles",
435 : nValueCount);
436 0 : return nullptr;
437 : }
438 1285 : for (int iValue = 0; iValue < nValueCount; iValue++)
439 1280 : poDS->padfMatrixValues[iValue] = poDS->ReadFloat();
440 : }
441 : else
442 : {
443 5 : VSIFSeekL(poDS->fp, 8 * nValueCount, SEEK_CUR);
444 : }
445 :
446 : // Read pairs till we run out, trying to find a few items that
447 : // have special meaning to us.
448 10 : poDS->nRasterXSize = 0;
449 10 : poDS->nRasterYSize = 0;
450 10 : int nBandCount = 0;
451 :
452 20 : while (poDS->ReadPair(osObjName, nObjCode) && nObjCode != 254)
453 : {
454 10 : if (osObjName == "dim" && nObjCode % 256 == R_INTSXP)
455 : {
456 10 : const int nCount = poDS->ReadInteger();
457 10 : if (nCount == 2)
458 : {
459 0 : poDS->nRasterXSize = poDS->ReadInteger();
460 0 : poDS->nRasterYSize = poDS->ReadInteger();
461 0 : nBandCount = 1;
462 : }
463 10 : else if (nCount == 3)
464 : {
465 10 : poDS->nRasterXSize = poDS->ReadInteger();
466 10 : poDS->nRasterYSize = poDS->ReadInteger();
467 10 : nBandCount = poDS->ReadInteger();
468 : }
469 : else
470 : {
471 0 : CPLError(CE_Failure, CPLE_AppDefined,
472 : "R 'dim' dimension wrong.");
473 0 : return nullptr;
474 : }
475 : }
476 0 : else if (nObjCode % 256 == R_REALSXP)
477 : {
478 0 : int nCount = poDS->ReadInteger();
479 0 : while (nCount > 0 && !VSIFEofL(poDS->fp))
480 : {
481 0 : nCount--;
482 0 : poDS->ReadFloat();
483 : }
484 : }
485 0 : else if (nObjCode % 256 == R_INTSXP)
486 : {
487 0 : int nCount = poDS->ReadInteger();
488 0 : while (nCount > 0 && !VSIFEofL(poDS->fp))
489 : {
490 0 : nCount--;
491 0 : poDS->ReadInteger();
492 : }
493 : }
494 0 : else if (nObjCode % 256 == R_STRSXP)
495 : {
496 0 : int nCount = poDS->ReadInteger();
497 0 : while (nCount > 0 && !VSIFEofL(poDS->fp))
498 : {
499 0 : nCount--;
500 0 : poDS->ReadString();
501 : }
502 : }
503 0 : else if (nObjCode % 256 == R_CHARSXP)
504 : {
505 0 : poDS->ReadString();
506 : }
507 : }
508 :
509 10 : if (poDS->nRasterXSize == 0)
510 : {
511 0 : CPLError(CE_Failure, CPLE_AppDefined,
512 : "Failed to find dim dimension information for R dataset.");
513 0 : return nullptr;
514 : }
515 :
516 20 : if (!GDALCheckDatasetDimensions(poDS->nRasterXSize, poDS->nRasterYSize) ||
517 10 : !GDALCheckBandCount(nBandCount, TRUE))
518 : {
519 0 : return nullptr;
520 : }
521 :
522 10 : GIntBig result = 0;
523 10 : bool ok = SafeMult(nBandCount, poDS->nRasterXSize, &result);
524 10 : ok &= SafeMult(result, poDS->nRasterYSize, &result);
525 :
526 10 : if (!ok || nValueCount < result)
527 : {
528 0 : CPLError(CE_Failure, CPLE_AppDefined, "Not enough pixel data.");
529 0 : return nullptr;
530 : }
531 :
532 : // Create the raster band object(s).
533 24 : for (int iBand = 0; iBand < nBandCount; iBand++)
534 : {
535 0 : std::unique_ptr<GDALRasterBand> poBand;
536 :
537 14 : if (poDS->bASCII)
538 7 : poBand = std::make_unique<RRasterBand>(
539 7 : poDS.get(), iBand + 1,
540 14 : poDS->padfMatrixValues + static_cast<size_t>(iBand) *
541 7 : poDS->nRasterXSize *
542 14 : poDS->nRasterYSize);
543 : else
544 : {
545 14 : poBand = RawRasterBand::Create(
546 7 : poDS.get(), iBand + 1, poDS->fp,
547 14 : poDS->nStartOfData +
548 7 : static_cast<vsi_l_offset>(poDS->nRasterXSize) *
549 7 : poDS->nRasterYSize * 8 * iBand,
550 7 : 8, poDS->nRasterXSize * 8, GDT_Float64,
551 : RawRasterBand::ByteOrder::ORDER_BIG_ENDIAN,
552 7 : RawRasterBand::OwnFP::NO);
553 7 : if (!poBand)
554 0 : return nullptr;
555 : }
556 14 : poDS->SetBand(iBand + 1, std::move(poBand));
557 : }
558 :
559 : // Initialize any PAM information.
560 10 : poDS->SetDescription(poOpenInfo->pszFilename);
561 10 : poDS->TryLoadXML();
562 :
563 : // Check for overviews.
564 10 : poDS->oOvManager.Initialize(poDS.get(), poOpenInfo->pszFilename);
565 :
566 10 : return poDS.release();
567 : }
568 :
569 : /************************************************************************/
570 : /* GDALRegister_R() */
571 : /************************************************************************/
572 :
573 1595 : void GDALRegister_R()
574 :
575 : {
576 1595 : if (GDALGetDriverByName("R") != nullptr)
577 302 : return;
578 :
579 1293 : GDALDriver *poDriver = new GDALDriver();
580 :
581 1293 : poDriver->SetDescription("R");
582 1293 : poDriver->SetMetadataItem(GDAL_DCAP_RASTER, "YES");
583 1293 : poDriver->SetMetadataItem(GDAL_DMD_LONGNAME, "R Object Data Store");
584 1293 : poDriver->SetMetadataItem(GDAL_DMD_HELPTOPIC, "drivers/raster/r.html");
585 1293 : poDriver->SetMetadataItem(GDAL_DMD_EXTENSION, "rda");
586 1293 : poDriver->SetMetadataItem(GDAL_DMD_CREATIONDATATYPES, "Float32");
587 1293 : poDriver->SetMetadataItem(
588 : GDAL_DMD_CREATIONOPTIONLIST,
589 : "<CreationOptionList>"
590 : " <Option name='ASCII' type='boolean' description='For ASCII output, "
591 : "default NO'/>"
592 : " <Option name='COMPRESS' type='boolean' description='Produced "
593 : "Compressed output, default YES'/>"
594 1293 : "</CreationOptionList>");
595 :
596 1293 : poDriver->SetMetadataItem(GDAL_DCAP_VIRTUALIO, "YES");
597 :
598 1293 : poDriver->pfnOpen = RDataset::Open;
599 1293 : poDriver->pfnIdentify = RDataset::Identify;
600 1293 : poDriver->pfnCreateCopy = RCreateCopy;
601 :
602 1293 : GetGDALDriverManager()->RegisterDriver(poDriver);
603 : }
|
