1 /******************************************************************************
2 * Copyright (c) 2016, Hobu Inc. (info@hobu.co)
3 *
4 * All rights reserved.
5 *
6 * Redistribution and use in source and binary forms, with or without
7 * modification, are permitted provided that the following
8 * conditions are met:
9 *
10 *     * Redistributions of source code must retain the above copyright
11 *       notice, this list of conditions and the following disclaimer.
12 *     * Redistributions in binary form must reproduce the above copyright
13 *       notice, this list of conditions and the following disclaimer in
14 *       the documentation and/or other materials provided
15 *       with the distribution.
16 *     * Neither the name of Hobu, Inc. nor the
17 *       names of its contributors may be used to endorse or promote
18 *       products derived from this software without specific prior
19 *       written permission.
20 *
21 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
22 * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
23 * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
24 * FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
25 * COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
26 * INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
27 * BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS
28 * OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED
29 * AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
30 * OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT
31 * OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY
32 * OF SUCH DAMAGE.
33 ****************************************************************************/
34 
35 #include <pdal/pdal_test_main.hpp>
36 #include <pdal/util/FileUtils.hpp>
37 #include <pdal/private/gdal/Raster.hpp>
38 #include <filters/RangeFilter.hpp>
39 #include <io/BufferReader.hpp>
40 #include <io/FauxReader.hpp>
41 #include <io/GDALWriter.hpp>
42 #include <io/LasReader.hpp>
43 #include <io/TextReader.hpp>
44 #include <io/private/GDALGrid.hpp>
45 #include "Support.hpp"
46 
47 #include <iostream>
48 #include <sstream>
49 
50 namespace pdal
51 {
52 
53 namespace
54 {
55 
runGdalWriter(const Options & wo,const std::string & infile,const std::string & outfile,const std::string & values)56 void runGdalWriter(const Options& wo, const std::string& infile,
57     const std::string& outfile, const std::string& values)
58 {
59     auto run = [=](bool streamMode)
60     {
61         FileUtils::deleteFile(outfile);
62 
63         Options ro;
64         ro.add("filename", infile);
65 
66         TextReader r;
67         r.setOptions(ro);
68 
69         GDALWriter w;
70         w.setOptions(wo);
71         w.setInput(r);
72 
73         if (streamMode)
74         {
75             FixedPointTable t(100);
76 
77             w.prepare(t);
78             w.execute(t);
79         }
80         else
81         {
82             PointTable t;
83 
84             w.prepare(t);
85             w.execute(t);
86         }
87 
88         std::istringstream iss(values);
89 
90         std::vector<double> arr;
91         while (true)
92         {
93             double d;
94             iss >> d;
95             if (!iss)
96                 break;
97             arr.push_back(d);
98         }
99 
100         gdal::Raster raster(outfile, "GTiff");
101         if (raster.open() != gdal::GDALError::None)
102         {
103             throw pdal_error(raster.errorMsg());
104         }
105         std::vector<double> data;
106         raster.readBand(data, 1);
107         int row = 0;
108         int col = 0;
109 
110         for (size_t i = 0; i < arr.size(); ++i)
111         {
112             EXPECT_NEAR(arr[i], data[i], .001) << "Error row/col = " <<
113                 row << "/" << col << std::endl;
114             if (++col == raster.width())
115             {
116                 col = 0;
117                 row++;
118             }
119         }
120     };
121 
122     run(false);
123     run(true);
124 }
125 
runGdalWriter2(const Options & wo,const std::string & outfile,const std::string & values,bool stream)126 void runGdalWriter2(const Options& wo, const std::string& outfile,
127     const std::string& values, bool stream)
128 {
129     FileUtils::deleteFile(outfile);
130 
131     Options ro;
132     ro.add("filename", Support::datapath("gdal/grid.txt"));
133 
134     Options ro2;
135     ro2.add("filename", Support::datapath("gdal/grid2.txt"));
136 
137     TextReader r;
138     r.setOptions(ro);
139 
140     TextReader r2;
141     r2.setOptions(ro2);
142 
143     GDALWriter w;
144     w.setOptions(wo);
145     w.setInput(r);
146     w.setInput(r2);
147 
148     if (!stream)
149     {
150         PointTable t;
151 
152         w.prepare(t);
153         w.execute(t);
154     }
155     else
156     {
157         FixedPointTable t(10);
158 
159         w.prepare(t);
160         w.execute(t);
161     }
162 
163     using namespace gdal;
164 
165     std::istringstream iss(values);
166 
167     std::vector<double> arr;
168     while (true)
169     {
170         double d;
171         iss >> d;
172         if (!iss)
173             break;
174         arr.push_back(d);
175     }
176 
177     gdal::Raster raster(outfile, "GTiff");
178     if (raster.open() != gdal::GDALError::None)
179     {
180         throw pdal_error(raster.errorMsg());
181     }
182     std::vector<double> data;
183     raster.readBand(data, 1);
184     for (size_t i = 0; i < arr.size(); ++i)
185         EXPECT_NEAR(arr[i], data[i], .001);
186 }
187 
188 }
189 
TEST(GDALWriterTest,min)190 TEST(GDALWriterTest, min)
191 {
192     std::string infile = Support::datapath("gdal/grid.txt");
193     std::string outfile = Support::temppath("tmp.tif");
194 
195     Options wo;
196     wo.add("gdaldriver", "GTiff");
197     wo.add("output_type", "min");
198     wo.add("resolution", 1);
199     wo.add("radius", .7071);
200     wo.add("filename", outfile);
201 
202     const std::string output =
203         "5.000 -9999.000     7.000     8.000     8.900 "
204         "4.000 -9999.000     6.000     7.000     8.000 "
205         "3.000     4.000     5.000     5.400     6.400 "
206         "2.000     3.000     4.000     4.400     5.400 "
207         "1.000     2.000     3.000     4.000     5.000 ";
208 
209     runGdalWriter(wo, infile, outfile, output);
210 }
211 
TEST(GDALWriterTest,min2)212 TEST(GDALWriterTest, min2)
213 {
214     std::string outfile = Support::temppath("tmp.tif");
215 
216     Options wo;
217     wo.add("gdaldriver", "GTiff");
218     wo.add("output_type", "min");
219     wo.add("resolution", 1);
220     wo.add("radius", .7071);
221     wo.add("filename", outfile);
222 
223     Options wo2 = wo;
224 //    wo2.add("bounds", "([-2, 4.7],[-2, 6.5])");
225     wo2.add("origin_x", -2);
226     wo2.add("origin_y", -2);
227     wo2.add("width", 7);
228     wo2.add("height", 9);
229 
230     const std::string output =
231         "-9999.000 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00    -1.00"
232         "-9999.000 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 "
233         "-9999.000 -9999.00     5.00 -9999.00     7.00     8.00     8.90 "
234         "-9999.000 -9999.00     4.00 -9999.00     6.00     7.00     8.00 "
235         "-9999.000 -9999.00     3.00     4.00     5.00     5.40     6.40 "
236         "-9999.000 -9999.00     2.00     3.00     4.00     4.40     5.40 "
237         "-9999.000 -9999.00     1.00     2.00     3.00     4.00     5.00 "
238         "   -1.000    -1.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 "
239         "   -1.000    -1.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00";
240 
241     runGdalWriter2(wo, outfile, output, false);
242     runGdalWriter2(wo2, outfile, output, false);
243     runGdalWriter2(wo2, outfile, output, true);
244 }
245 
TEST(GDALWriterTest,minWindow)246 TEST(GDALWriterTest, minWindow)
247 {
248     std::string infile = Support::datapath("gdal/grid.txt");
249     std::string outfile = Support::temppath("tmp.tif");
250 
251     Options wo;
252     wo.add("gdaldriver", "GTiff");
253     wo.add("output_type", "min");
254     wo.add("resolution", 1);
255     wo.add("radius", .7071);
256     wo.add("filename", outfile);
257     wo.add("window_size", 2);
258 
259     const std::string output =
260         "5.000     5.457     7.000     8.000     8.900 "
261         "4.000     4.848     6.000     7.000     8.000 "
262         "3.000     4.000     5.000     5.400     6.400 "
263         "2.000     3.000     4.000     4.400     5.400 "
264         "1.000     2.000     3.000     4.000     5.000 ";
265 
266     runGdalWriter(wo, infile, outfile, output);
267 }
268 
TEST(GDALWriterTest,max)269 TEST(GDALWriterTest, max)
270 {
271     std::string infile = Support::datapath("gdal/grid.txt");
272     std::string outfile = Support::temppath("tmp.tif");
273 
274     Options wo;
275     wo.add("gdaldriver", "GTiff");
276     wo.add("output_type", "max");
277     wo.add("resolution", 1);
278     wo.add("radius", .7071);
279     wo.add("filename", outfile);
280 
281     const std::string output =
282         "5.000 -9999.000     7.000     8.000     9.100 "
283         "4.000 -9999.000     6.000     7.000     8.000 "
284         "3.000     4.000     5.000     6.000     7.000 "
285         "2.000     3.000     4.400     5.400     6.400 "
286         "1.000     2.000     3.000     4.400     5.400 ";
287 
288     runGdalWriter(wo, infile, outfile, output);
289 }
290 
TEST(GDALWriterTest,maxWindow)291 TEST(GDALWriterTest, maxWindow)
292 {
293     std::string infile = Support::datapath("gdal/grid.txt");
294     std::string outfile = Support::temppath("tmp.tif");
295 
296     Options wo;
297     wo.add("gdaldriver", "GTiff");
298     wo.add("output_type", "max");
299     wo.add("resolution", 1);
300     wo.add("radius", .7071);
301     wo.add("filename", outfile);
302     wo.add("window_size", 2);
303 
304     const std::string output =
305         "5.000     5.500     7.000     8.000     9.100 "
306         "4.000     4.942     6.000     7.000     8.000 "
307         "3.000     4.000     5.000     6.000     7.000 "
308         "2.000     3.000     4.400     5.400     6.400 "
309         "1.000     2.000     3.000     4.400     5.400 ";
310 
311     runGdalWriter(wo, infile, outfile, output);
312 }
313 
TEST(GDALWriterTest,mean)314 TEST(GDALWriterTest, mean)
315 {
316     std::string infile = Support::datapath("gdal/grid.txt");
317     std::string outfile = Support::temppath("tmp.tif");
318 
319     Options wo;
320     wo.add("gdaldriver", "GTiff");
321     wo.add("output_type", "mean");
322     wo.add("resolution", 1);
323     wo.add("radius", .7071);
324     wo.add("filename", outfile);
325 
326     const std::string output =
327         "5.000 -9999.000     7.000     8.000     8.967 "
328         "4.000 -9999.000     6.000     7.000     8.000 "
329         "3.000     4.000     5.000     5.700     6.700 "
330         "2.000     3.000     4.200     4.920     5.800 "
331         "1.000     2.000     3.000     4.200     5.200 ";
332 
333     runGdalWriter(wo, infile, outfile, output);
334 }
335 
TEST(GDALWriterTest,meanWindow)336 TEST(GDALWriterTest, meanWindow)
337 {
338     std::string infile = Support::datapath("gdal/grid.txt");
339     std::string outfile = Support::temppath("tmp.tif");
340 
341     Options wo;
342     wo.add("gdaldriver", "GTiff");
343     wo.add("output_type", "mean");
344     wo.add("resolution", 1);
345     wo.add("radius", .7071);
346     wo.add("filename", outfile);
347     wo.add("window_size", 2);
348 
349     const std::string output =
350         "5.000     5.478     7.000     8.000     8.967 "
351         "4.000     4.896     6.000     7.000     8.000 "
352         "3.000     4.000     5.000     5.700     6.700 "
353         "2.000     3.000     4.200     4.920     5.800 "
354         "1.000     2.000     3.000     4.200     5.200 ";
355 
356     runGdalWriter(wo, infile, outfile, output);
357 }
358 
TEST(GDALWriterTest,idw)359 TEST(GDALWriterTest, idw)
360 {
361     std::string infile = Support::datapath("gdal/grid.txt");
362     std::string outfile = Support::temppath("tmp.tif");
363 
364     Options wo;
365     wo.add("gdaldriver", "GTiff");
366     wo.add("output_type", "idw");
367     wo.add("resolution", 1);
368     wo.add("radius", .7071);
369     wo.add("filename", outfile);
370 
371     const std::string output =
372         "5.000 -9999.000     7.000     8.000     9.000 "
373         "4.000 -9999.000     6.000     7.000     8.000 "
374         "3.000     4.000     5.000     6.000     7.000 "
375         "2.000     3.000     4.000     5.000     6.000 "
376         "1.000     2.000     3.000     4.000     5.000 ";
377 
378     runGdalWriter(wo, infile, outfile, output);
379 }
380 
TEST(GDALWriterTest,idwWindow)381 TEST(GDALWriterTest, idwWindow)
382 {
383     std::string infile = Support::datapath("gdal/grid.txt");
384     std::string outfile = Support::temppath("tmp.tif");
385 
386     Options wo;
387     wo.add("gdaldriver", "GTiff");
388     wo.add("output_type", "idw");
389     wo.add("resolution", 1);
390     wo.add("radius", .7071);
391     wo.add("filename", outfile);
392     wo.add("window_size", 2);
393 
394     const std::string output =
395         "5.000     5.500     7.000     8.000     9.000 "
396         "4.000     4.905     6.000     7.000     8.000 "
397         "3.000     4.000     5.000     6.000     7.000 "
398         "2.000     3.000     4.000     5.000     6.000 "
399         "1.000     2.000     3.000     4.000     5.000 ";
400 
401     runGdalWriter(wo, infile, outfile, output);
402 }
403 
TEST(GDALWriterTest,count)404 TEST(GDALWriterTest, count)
405 {
406     std::string infile = Support::datapath("gdal/grid.txt");
407     std::string outfile = Support::temppath("tmp.tif");
408 
409     Options wo;
410     wo.add("gdaldriver", "GTiff");
411     wo.add("output_type", "count");
412     wo.add("resolution", 1);
413     wo.add("radius", .7071);
414     wo.add("filename", outfile);
415 
416     const std::string output =
417         "1.000     0.000     1.000     1.000     3.000 "
418         "1.000     0.000     1.000     1.000     1.000 "
419         "1.000     1.000     1.000     2.000     2.000 "
420         "1.000     1.000     2.000     5.000     4.000 "
421         "1.000     1.000     1.000     2.000     2.000 ";
422 
423     runGdalWriter(wo, infile, outfile, output);
424 }
425 
TEST(GDALWriterTest,stdev)426 TEST(GDALWriterTest, stdev)
427 {
428     std::string infile = Support::datapath("gdal/grid.txt");
429     std::string outfile = Support::temppath("tmp.tif");
430 
431     Options wo;
432     wo.add("gdaldriver", "GTiff");
433     wo.add("output_type", "stdev");
434     wo.add("resolution", 1);
435     wo.add("radius", .7071);
436     wo.add("filename", outfile);
437 
438     const std::string output =
439         "0.000 -9999.000     0.000     0.000     0.094 "
440         "0.000 -9999.000     0.000     0.000     0.000 "
441         "0.000     0.000     0.000     0.300     0.300 "
442         "0.000     0.000     0.200     0.449     0.424 "
443         "0.000     0.000     0.000     0.200     0.200 ";
444 
445     runGdalWriter(wo, infile, outfile, output);
446 }
447 
TEST(GDALWriterTest,stdevWindow)448 TEST(GDALWriterTest, stdevWindow)
449 {
450     std::string infile = Support::datapath("gdal/grid.txt");
451     std::string outfile = Support::temppath("tmp.tif");
452 
453     Options wo;
454     wo.add("gdaldriver", "GTiff");
455     wo.add("output_type", "stdev");
456     wo.add("resolution", 1);
457     wo.add("radius", .7071);
458     wo.add("filename", outfile);
459     wo.add("window_size", 2);
460 
461     const std::string output =
462         "0.000     0.021     0.000     0.000     0.094 "
463         "0.000     0.045     0.000     0.000     0.000 "
464         "0.000     0.000     0.000     0.300     0.300 "
465         "0.000     0.000     0.200     0.449     0.424 "
466         "0.000     0.000     0.000     0.200     0.200 ";
467 
468     runGdalWriter(wo, infile, outfile, output);
469 }
470 
TEST(GDALWriterTest,additionalDim)471 TEST(GDALWriterTest, additionalDim)
472 {
473     std::string outfile(Support::temppath("out.tif"));
474 
475     FileUtils::deleteFile(outfile);
476 
477     LasReader r;
478     Options ro;
479 
480     ro.add("filename", Support::datapath("las/extrabytes.las"));
481 
482     GDALWriter w;
483     Options wo;
484 
485     wo.add("dimension", "Flags1");
486     wo.add("resolution", 2);
487     wo.add("radius", 2.7);
488     wo.add("filename", outfile);
489 
490     r.setOptions(ro);
491     w.setOptions(wo);
492     w.setInput(r);
493 
494     PointTable t;
495 
496     EXPECT_NO_THROW(w.prepare(t));
497 
498     Options wo2;
499 
500     wo2.add("dimension", "Flag55");
501     wo2.add("resolution", 2);
502     wo2.add("radius", 2.7);
503     wo2.add("filename", outfile);
504 
505     w.setOptions(wo2);
506 
507     PointTable t2;
508     EXPECT_THROW(w.prepare(t2), pdal_error);
509 }
510 
TEST(GDALWriterTest,btbad)511 TEST(GDALWriterTest, btbad)
512 {
513     std::string infile = Support::datapath("gdal/grid.txt");
514     std::string outfile = Support::temppath("tmp.tif");
515 
516     Options wo;
517     wo.add("gdaldriver", "BT");
518     wo.add("output_type", "min");
519     wo.add("resolution", 1);
520     wo.add("radius", .7071);
521     wo.add("data_type", "double");
522     wo.add("filename", outfile);
523 
524     const std::string output =
525         "5.000 -9999.000     7.000     8.000     8.900 "
526         "4.000 -9999.000     6.000     7.000     8.000 "
527         "3.000     4.000     5.000     5.000     6.000 "
528         "2.000     3.000     4.000     4.000     5.000 "
529         "1.000     2.000     3.000     4.000     5.000 ";
530 
531     EXPECT_THROW(runGdalWriter(wo, infile, outfile, output), pdal_error);
532 }
533 
TEST(GDALWriterTest,btint)534 TEST(GDALWriterTest, btint)
535 {
536     std::string outfile = Support::temppath("tmp.bt");
537     FileUtils::deleteFile(outfile);
538 
539     Options ro;
540     ro.add("filename", Support::datapath("gdal/grid.txt"));
541 
542     TextReader r;
543     r.setOptions(ro);
544 
545     GDALWriter w;
546     Options wo;
547     wo.add("gdaldriver", "BT");
548     wo.add("output_type", "min");
549     wo.add("resolution", 1);
550     wo.add("radius", .7071);
551     wo.add("data_type", "int32");
552     wo.add("filename", outfile);
553     w.setOptions(wo);
554     w.setInput(r);
555 
556     PointTable t;
557 
558     w.prepare(t);
559     w.execute(t);
560 
561     using namespace gdal;
562 
563     const std::string values =
564         "5.000 -9999.000     7.000     8.000     9.000 "
565         "4.000 -9999.000     6.000     7.000     8.000 "
566         "3.000     4.000     5.000     5.000     6.000 "
567         "2.000     3.000     4.000     4.000     5.000 "
568         "1.000     2.000     3.000     4.000     5.000 ";
569     std::istringstream iss(values);
570 
571     std::vector<double> arr;
572     while (true)
573     {
574         double d;
575         iss >> d;
576         if (!iss)
577             break;
578         arr.push_back(d);
579     }
580 
581     gdal::Raster raster(outfile, "BT");
582     if (raster.open() != gdal::GDALError::None)
583     {
584         throw pdal_error(raster.errorMsg());
585     }
586     std::vector<int32_t> data;
587     raster.readBand(data, 1);
588     for (size_t i = 0; i < arr.size(); ++i)
589         EXPECT_NEAR(arr[i], data[i], .001);
590 }
591 
TEST(GDALWriterTest,no_points)592 TEST(GDALWriterTest, no_points)
593 {
594     std::string outfile(Support::temppath("out.tif"));
595     FileUtils::deleteFile(outfile);
596 
597     LasReader r;
598     Options ro;
599     ro.add("filename", Support::datapath("las/no-points.las"));
600     r.setOptions(ro);
601 
602     RangeFilter f;
603     f.setInput(r);
604     Options fo;
605     fo.add("limits", "Classification[2:2]");
606     f.setOptions(fo);
607 
608     GDALWriter w;
609     w.setInput(f);
610     Options wo;
611     wo.add("resolution", 2);
612     wo.add("filename", outfile);
613     w.setOptions(wo);
614 
615     PointTable t;
616     w.prepare(t);
617     EXPECT_THROW(w.execute(t), pdal_error);
618 }
619 
620 // Check that we don't crash with bad X/Y's in grid_bounds.txt and a preset
621 // grid bounds.
TEST(GDALWriterTest,bounds)622 TEST(GDALWriterTest, bounds)
623 {
624     std::string infile = Support::datapath("gdal/grid_bounds.txt");
625     std::string outfile = Support::temppath("tmp.tif");
626 
627     Options wo;
628     wo.add("gdaldriver", "GTiff");
629     wo.add("output_type", "mean");
630     wo.add("resolution", 1);
631     wo.add("radius", .7071);
632     wo.add("filename", outfile);
633 //    wo.add("bounds", "([0, 4.5],[0, 4.5])");
634     wo.add("origin_x", 0);
635     wo.add("origin_y", 0);
636     wo.add("height", 5);
637     wo.add("width", 5);
638 
639     const std::string output =
640         "5.000 -9999.000     7.000     8.000     8.967 "
641         "4.000 -9999.000     6.000     7.000     8.000 "
642         "3.000     4.000     5.000     5.700     6.700 "
643         "2.000     3.000     4.200     4.920     5.800 "
644         "1.000     2.000     3.000     4.200     5.200 ";
645 
646     runGdalWriter(wo, infile, outfile, output);
647 }
648 
649 // Make sure we reset bounds when starting a new file.
TEST(GDALWriterTest,issue_2074)650 TEST(GDALWriterTest, issue_2074)
651 {
652     FileUtils::deleteFile(Support::temppath("gdal.tif"));
653     FileUtils::deleteFile(Support::temppath("gdal1.tif"));
654     FileUtils::deleteFile(Support::temppath("gdal2.tif"));
655 
656     PointTable t;
657     t.layout()->registerDim(Dimension::Id::X);
658     t.layout()->registerDim(Dimension::Id::Y);
659     t.layout()->registerDim(Dimension::Id::Z);
660 
661     PointViewPtr v1(new PointView(t));
662     PointViewPtr v2(new PointView(t));
663 
664     v1->setField(Dimension::Id::X, 0, 0);
665     v1->setField(Dimension::Id::Y, 0, 0);
666     v1->setField(Dimension::Id::Z, 0, 0);
667     v1->setField(Dimension::Id::X, 1, 1);
668     v1->setField(Dimension::Id::Y, 1, 1);
669     v1->setField(Dimension::Id::Z, 1, 0);
670 
671     v2->setField(Dimension::Id::X, 0, 2);
672     v2->setField(Dimension::Id::Y, 0, 2);
673     v2->setField(Dimension::Id::Z, 0, 0);
674     v2->setField(Dimension::Id::X, 1, 3);
675     v2->setField(Dimension::Id::Y, 1, 3);
676     v2->setField(Dimension::Id::Z, 1, 0);
677 
678     BufferReader r;
679     r.addView(v1);
680     r.addView(v2);
681 
682     Options wo;
683     wo.add("gdaldriver", "GTiff");
684     wo.add("output_type", "mean");
685     wo.add("filename", Support::temppath("gdal#.tif"));
686     wo.add("resolution", .5);
687 
688     GDALWriter w;
689 
690     w.setOptions(wo);
691     w.setInput(r);
692 
693     w.prepare(t);
694     w.execute(t);
695 
696     gdal::Raster raster(Support::temppath("gdal1.tif"), "GTiff");
697     if (raster.open() != gdal::GDALError::None)
698     {
699         throw pdal_error(raster.errorMsg());
700     }
701 
702     gdal::Raster raster2(Support::temppath("gdal2.tif"), "GTiff");
703     if (raster2.open() != gdal::GDALError::None)
704     {
705         throw pdal_error(raster2.errorMsg());
706     }
707     // Make sure the output files are the same size.
708     EXPECT_EQ(raster.width(), 3);
709     EXPECT_EQ(raster.height(), 3);
710     EXPECT_EQ(raster.width(), raster2.width());
711     EXPECT_EQ(raster.height(), raster2.height());
712 
713     // Make sure that we accumulate data for two views to one output raster.
714     Options wo2;
715     wo2.add("gdaldriver", "GTiff");
716     wo2.add("output_type", "mean");
717     wo2.add("filename", Support::temppath("gdal.tif"));
718     wo2.add("resolution", .5);
719 
720     GDALWriter w2;
721     w2.setOptions(wo2);
722     w2.setInput(r);
723 
724     w2.prepare(t);
725     w2.execute(t);
726 
727     gdal::Raster raster3(Support::temppath("gdal.tif"), "GTiff");
728     if (raster3.open() != gdal::GDALError::None)
729     {
730         throw pdal_error(raster3.errorMsg());
731     }
732     EXPECT_EQ(raster3.width(), 7);
733     EXPECT_EQ(raster3.height(), 7);
734 }
735 
736 // If the radius is sufficiently large, make sure the grid is filled.
TEST(GDALWriterTest,issue_2545)737 TEST(GDALWriterTest, issue_2545)
738 {
739     TextReader r;
740     Options rOpts;
741 
742     rOpts.add("filename", Support::datapath("gdal/issue_2545.txt"));
743     r.setOptions(rOpts);
744 
745     std::string outfile(Support::temppath("gdal.tif"));
746 
747     GDALWriter w;
748     Options wOpts;
749     wOpts.add("resolution", 1);
750     wOpts.add("radius", 10);
751     wOpts.add("output_type", "idw");
752     wOpts.add("gdaldriver", "GTiff");
753     wOpts.add("origin_x", .5);
754     wOpts.add("origin_y", .5);
755     wOpts.add("width", 7);
756     wOpts.add("height", 7);
757 //    wOpts.add("bounds", BOX2D(.5, .5, 6.5, 6.5));
758     wOpts.add("filename", outfile);
759 
760     w.setOptions(wOpts);
761     w.setInput(r);
762 
763     PointTable t;
764     w.prepare(t);
765     w.execute(t);
766 
767     gdal::Raster raster(outfile, "GTiff");
768     if (raster.open() != gdal::GDALError::None)
769         throw pdal_error(raster.errorMsg());
770 
771     EXPECT_EQ(raster.width(), 7);
772     EXPECT_EQ(raster.height(), 7);
773 
774     auto index = [](size_t row, size_t col)
775     {
776         return (row * 7) + col;
777     };
778 
779     std::vector<double> data;
780     raster.readBand(data, 1);
781     EXPECT_EQ(data[index(1, 0)], 3.0);
782     EXPECT_EQ(data[index(2, 4)], 0.);
783     EXPECT_EQ(data[index(5, 5)], 10.0);
784     EXPECT_EQ(data[index(6, 0)], 5.0);
785     for (size_t i = 0; i < data.size(); ++i)
786         EXPECT_TRUE(data[i] <= 10.0 && data[i] >= 0);
787 }
788 
789 // Test that using the alternate grid specification works
TEST(GDALWriterTest,alternate_grid)790 TEST(GDALWriterTest, alternate_grid)
791 {
792     FauxReader r;
793     Options opts;
794     opts.add("count", 1000);
795     opts.add("mode", "constant");
796     r.setOptions(opts);
797 
798     std::string outfile(Support::temppath("test.tif"));
799 
800     {
801         Options opts;
802         opts.add("origin_x", "10.0");
803         opts.add("filename", outfile);
804         opts.add("resolution", 1.1);
805 
806         PointTable t;
807         GDALWriter w;
808         w.setOptions(opts);
809         w.setInput(r);
810         EXPECT_THROW(w.prepare(t), pdal_error);
811     }
812 
813     {
814         Options opts;
815         opts.add("origin_x", "10.0");
816         opts.add("origin_y", "20.0");
817         opts.add("width", 10);
818         opts.add("height", 5);
819         opts.add("bounds", "([1, 10],[5, 25])");
820         opts.add("filename", outfile);
821         opts.add("resolution", 1.1);
822 
823         PointTable t;
824         GDALWriter w;
825         w.setOptions(opts);
826         w.setInput(r);
827         EXPECT_THROW(w.prepare(t), pdal_error);
828     }
829 
830     {
831         Options opts;
832         opts.add("origin_x", "10.0");
833         opts.add("origin_y", "20.0");
834         opts.add("width", 10);
835         opts.add("height", 5);
836         opts.add("filename", outfile);
837         opts.add("resolution", 1);
838 
839         PointTable t;
840         GDALWriter w;
841         w.setOptions(opts);
842         w.setInput(r);
843         w.prepare(t);
844         w.execute(t);
845 
846         gdal::Raster raster(outfile, "GTiff");
847         raster.open();
848         EXPECT_EQ(raster.width(), 10);
849         EXPECT_EQ(raster.height(), 5);
850     }
851 }
852 
TEST(GDALWriterTest,srs)853 TEST(GDALWriterTest, srs)
854 {
855     auto test = [](const std::string& sourceSrs, const std::string& defaultSrs,
856         const std::string& overrideSrs, const std::string& testSrs)
857     {
858         std::string outfile(Support::temppath("out.tif"));
859 
860         TextReader t;
861         Options to;
862         to.add("filename", Support::datapath("text/with_edgeofflightline.txt"));
863         if (sourceSrs.size())
864             to.add("override_srs", sourceSrs);
865         t.setOptions(to);
866 
867         GDALWriter w;
868         Options wo;
869         wo.add("filename", outfile);
870         wo.add("origin_x", 0);
871         wo.add("origin_y", 0);
872         wo.add("width", 1000);
873         wo.add("height", 1000);
874         wo.add("resolution", 10);
875         if (defaultSrs.size())
876             wo.add("default_srs", defaultSrs);
877         if (overrideSrs.size())
878             wo.add("override_srs", overrideSrs);
879         w.setOptions(wo);
880         w.setInput(t);
881 
882         FileUtils::deleteFile(outfile);
883         PointTable table;
884         w.prepare(table);
885         w.execute(table);
886 
887         gdal::Raster raster(outfile);
888         raster.open();
889         EXPECT_EQ(raster.getSpatialRef(), testSrs);
890     };
891 
892     test("", "EPSG:4326", "", "EPSG:4326");
893     test("", "", "EPSG:4326", "EPSG:4326");
894     test("EPSG:4326", "EPSG:2030", "", "EPSG:4326");
895     test("EPSG:4326", "", "EPSG:2030", "EPSG:2030");
896     EXPECT_THROW(test("EPSG:4326", "EPSG:4326", "EPSG:2030", "EPSG:2030"), pdal_error);
897 }
898 
899 } // namespace pdal
900