SAGA API v9.10
Loading...
Searching...
No Matches
grid_operation.cpp
Go to the documentation of this file.
1
3// //
4// SAGA //
5// //
6// System for Automated Geoscientific Analyses //
7// //
8// Application Programming Interface //
9// //
10// Library: SAGA_API //
11// //
12//-------------------------------------------------------//
13// //
14// grid_operation.cpp //
15// //
16// Copyright (C) 2005 by Olaf Conrad //
17// //
18//-------------------------------------------------------//
19// //
20// This file is part of 'SAGA - System for Automated //
21// Geoscientific Analyses'. //
22// //
23// This library is free software; you can redistribute //
24// it and/or modify it under the terms of the GNU Lesser //
25// General Public License as published by the Free //
26// Software Foundation, either version 2.1 of the //
27// License, or (at your option) any later version. //
28// //
29// This library is distributed in the hope that it will //
30// be useful, but WITHOUT ANY WARRANTY; without even the //
31// implied warranty of MERCHANTABILITY or FITNESS FOR A //
32// PARTICULAR PURPOSE. See the GNU Lesser General Public //
33// License for more details. //
34// //
35// You should have received a copy of the GNU Lesser //
36// General Public License along with this program; if //
37// not, see <http://www.gnu.org/licenses/>. //
38// //
39//-------------------------------------------------------//
40// //
41// contact: Olaf Conrad //
42// Institute of Geography //
43// University of Goettingen //
44// Goldschmidtstr. 5 //
45// 37077 Goettingen //
46// Germany //
47// //
48// e-mail: oconrad@saga-gis.org //
49// //
51
52//---------------------------------------------------------
53#include <memory.h>
54
55#include "grid.h"
56
57
59// //
60// Data Assignments - Value //
61// //
63
64//---------------------------------------------------------
66{
67 double Value = Get_NoData_Value();
68
69 #pragma omp parallel for
70 for(int y=0; y<Get_NY(); y++) for(int x=0; x<Get_NX(); x++)
71 {
72 Set_Value(x, y, Value, false);
73 }
74
75 return( true );
76}
77
78//---------------------------------------------------------
79bool CSG_Grid::Assign(double Value)
80{
81 if( !is_Valid() )
82 {
83 return( false );
84 }
85
86 if( Value == 0. && !is_Cached() )
87 {
88 #pragma omp parallel for
89 for(int y=0; y<Get_NY(); y++)
90 {
91 memset(m_Values[y], 0, Get_nLineBytes());
92 }
93 }
94 else
95 {
96 #pragma omp parallel for
97 for(int y=0; y<Get_NY(); y++) for(int x=0; x<Get_NX(); x++)
98 {
99 Set_Value(x, y, Value);
100 }
101 }
102
103 //-----------------------------------------------------
105
106 m_Statistics.Invalidate();
107
108 Set_Update_Flag(false);
109
110 return( true );
111}
112
113
115// //
116// Data Assignments - Grid //
117// //
119
120//---------------------------------------------------------
121bool CSG_Grid::Assign(CSG_Data_Object *pObject, bool bProgress)
122{
123 return( pObject && Assign(pObject->asGrid(), CSG_Grid_Resampling::Undefined) );
124}
125
126bool CSG_Grid::Assign(CSG_Grid *pGrid, CSG_Grid_Resampling Method, bool bProgress)
127{
128 if( !is_Valid() || !pGrid || !pGrid->is_Valid() || is_Intersecting(pGrid->Get_Extent()) == INTERSECTION_None )
129 {
130 return( false );
131 }
132
133 bool bResult = false;
134
135 //---------------------------------------------------------
136 if( m_System == pGrid->m_System )
137 {
138 #pragma omp parallel for
139 for(int y=0; y<Get_NY(); y++) for(int x=0; x<Get_NX(); x++)
140 {
141 if( pGrid->is_NoData(x, y) )
142 {
143 Set_NoData(x, y);
144 }
145 else
146 {
147 Set_Value(x, y, pGrid->asDouble(x, y));
148 }
149 }
150
151 bResult = true; bProgress = false;
152 }
153
154 //---------------------------------------------------------
155 else if( Get_Cellsize() == pGrid->Get_Cellsize() // No-Scaling...
156 && fmod(Get_XMin() - pGrid->Get_XMin(), Get_Cellsize()) == 0.
157 && fmod(Get_YMin() - pGrid->Get_YMin(), Get_Cellsize()) == 0. )
158 {
159 bResult = _Assign_Interpolated(pGrid, CSG_Grid_Resampling::NearestNeighbour, bProgress);
160 }
161
162 //---------------------------------------------------------
163 else switch( Method )
164 {
168 case CSG_Grid_Resampling::Bicubic_2 : bResult = _Assign_Interpolated(pGrid, Method , bProgress); break;
169
171 case CSG_Grid_Resampling::Mean_Cells : bResult = _Assign_MeanValue (pGrid, Method != CSG_Grid_Resampling::Mean_Nodes, bProgress); break;
172
174 case CSG_Grid_Resampling::Maximum : bResult = _Assign_ExtremeValue(pGrid, Method == CSG_Grid_Resampling::Maximum , bProgress); break;
175
176 case CSG_Grid_Resampling::Majority : bResult = _Assign_Majority (pGrid , bProgress); break;
177
178 default : bResult = Get_Cellsize() < pGrid->Get_Cellsize()
179 ? _Assign_Interpolated(pGrid, CSG_Grid_Resampling::Bicubic_2, bProgress) // Down-Scaling...
180 : _Assign_MeanValue (pGrid, true , bProgress); // Up-Scaling...
181 break;
182 }
183
184 //---------------------------------------------------------
185 if( bResult )
186 {
187 Set_Unit(pGrid->Get_Unit());
188
189 if( pGrid->Get_Projection().is_Okay() )
190 {
191 Get_Projection() = pGrid->Get_Projection();
192 }
193
194 Get_History() = pGrid->Get_History();
195 }
196
197 //---------------------------------------------------------
198 if( bProgress )
199 {
201 }
202
203 return( bResult );
204}
205
206//---------------------------------------------------------
207bool CSG_Grid::_Assign_Interpolated(CSG_Grid *pGrid, CSG_Grid_Resampling Interpolation, bool bProgress)
208{
209 for(int y=0; y<Get_NY() && (!bProgress || SG_UI_Process_Set_Progress(y, Get_NY())); y++)
210 {
211 double py = Get_YMin() + y * Get_Cellsize();
212
213 #pragma omp parallel for
214 for(int x=0; x<Get_NX(); x++)
215 {
216 double px = Get_XMin() + x * Get_Cellsize(), z;
217
218 if( pGrid->Get_Value(px, py, z, Interpolation) )
219 {
220 Set_Value(x, y, z);
221 }
222 else
223 {
224 Set_NoData(x, y);
225 }
226 }
227 }
228
229 return( true );
230}
231
232//---------------------------------------------------------
233bool CSG_Grid::_Assign_ExtremeValue(CSG_Grid *pGrid, bool bMaximum, bool bProgress)
234{
235 if( Get_Cellsize() < pGrid->Get_Cellsize() )
236 {
237 return( false );
238 }
239
240 //-----------------------------------------------------
242
243 double d = pGrid->Get_Cellsize() / Get_Cellsize();
244
245 double ax = 0.5 + (pGrid->Get_XMin() - Get_XMin()) / Get_Cellsize();
246 double py = 0.5 + (pGrid->Get_YMin() - Get_YMin()) / Get_Cellsize();
247
248 for(int y=0; y<pGrid->Get_NY() && (!bProgress || SG_UI_Process_Set_Progress(y, pGrid->Get_NY())); y++, py+=d)
249 {
250 int iy = (int)floor(py);
251
252 if( iy >= 0 && iy < Get_NY() )
253 {
254 #pragma omp parallel for
255 for(int x=0; x<pGrid->Get_NX(); x++)
256 {
257 if( !pGrid->is_NoData(x, y) )
258 {
259 int ix = (int)floor(ax + x * d);
260
261 if( ix >= 0 && ix < Get_NX() )
262 {
263 double z = pGrid->asDouble(x, y);
264
265 if( is_NoData(ix, iy)
266 || (bMaximum == true && z > asDouble(ix, iy))
267 || (bMaximum == false && z < asDouble(ix, iy)) )
268 {
269 Set_Value(ix, iy, z);
270 }
271 }
272 }
273 }
274 }
275 }
276
277 //-----------------------------------------------------
278 return( true );
279}
280
281//---------------------------------------------------------
282bool CSG_Grid::_Assign_MeanValue(CSG_Grid *pGrid, bool bAreaProportional, bool bProgress)
283{
284 if( Get_Cellsize() < pGrid->Get_Cellsize() )
285 {
286 return( false );
287 }
288
289 //-----------------------------------------------------
290 double d = Get_Cellsize() / pGrid->Get_Cellsize();
291
292 double ox = (Get_XMin(true) - pGrid->Get_XMin()) / pGrid->Get_Cellsize();
293 double py = (Get_YMin(true) - pGrid->Get_YMin()) / pGrid->Get_Cellsize();
294
295 for(int y=0; y<Get_NY() && (!bProgress || SG_UI_Process_Set_Progress(y, Get_NY())); y++, py+=d)
296 {
297 int ay = (int)(bAreaProportional ? floor(py ) : ceil (py ));
298 int by = (int)(bAreaProportional ? ceil (py + d) : floor(py + d));
299
300 //-------------------------------------------------
301 #pragma omp parallel for
302 for(int x=0; x<Get_NX(); x++)
303 {
304 double px = ox + x * d;
305
306 int ax = (int)(bAreaProportional ? floor(px ) : ceil (px ));
307 int bx = (int)(bAreaProportional ? ceil (px + d) : floor(px + d));
308
309 CSG_Rect rMean(px, py, px + d, py + d); CSG_Simple_Statistics s;
310
311 for(int iy=ay; iy<=by; iy++)
312 {
313 if( iy >= 0 && iy < pGrid->Get_NY() )
314 {
315 for(int ix=ax; ix<=bx; ix++)
316 {
317 if( ix >= 0 && ix < pGrid->Get_NX() && !pGrid->is_NoData(ix, iy) )
318 {
319 if( bAreaProportional )
320 {
321 CSG_Rect r(ix - 0.5, iy - 0.5, ix + 0.5, iy + 0.5);
322
323 if( r.Intersect(rMean) )
324 {
325 s.Add_Value(pGrid->asDouble(ix, iy), r.Get_Area());
326 }
327 }
328 else
329 {
330 s.Add_Value(pGrid->asDouble(ix, iy));
331 }
332 }
333 }
334 }
335 }
336
337 //---------------------------------------------
338 if( s.Get_Count() > 0 )
339 {
340 Set_Value(x, y, s.Get_Mean());
341 }
342 else
343 {
344 Set_NoData(x, y);
345 }
346 }
347 }
348
349 //-----------------------------------------------------
350 return( true );
351}
352
353
355// //
357
358//---------------------------------------------------------
359bool CSG_Grid::_Assign_Majority(CSG_Grid *pGrid, bool bProgress)
360{
361 if( Get_Cellsize() < pGrid->Get_Cellsize() )
362 {
363 return( false );
364 }
365
366 //-----------------------------------------------------
368
370
371 //-----------------------------------------------------
372 int ay, by = (int)(1. + (((0 - 0.5) * Get_Cellsize() + Get_YMin()) - pGrid->Get_YMin()) / pGrid->Get_Cellsize());
373
374 for(int y=0; y<Get_NY() && (!bProgress || SG_UI_Process_Set_Progress(y, Get_NY())); y++)
375 {
376 ay = by; by = (int)(1. + (((y + 0.5) * Get_Cellsize() + Get_YMin()) - pGrid->Get_YMin()) / pGrid->Get_Cellsize());
377
378 if( ay < pGrid->Get_NY() && by > 0 )
379 {
380 if( ay < 0 ) { ay = 0; } if( by > pGrid->Get_NY() ) { by = pGrid->Get_NY(); }
381
382 int ax, bx = (int)(1. + (((0 - 0.5) * Get_Cellsize() + Get_XMin()) - pGrid->Get_XMin()) / pGrid->Get_Cellsize());
383
384 for(int x=0; x<Get_NX(); x++)
385 {
386 ax = bx; bx = (int)(1. + (((x + 0.5) * Get_Cellsize() + Get_XMin()) - pGrid->Get_XMin()) / pGrid->Get_Cellsize());
387
388 if( ax < pGrid->Get_NX() && bx > 0 )
389 {
390 CSG_Unique_Number_Statistics s;
391
392 if( ax < 0 ) { ax = 0; } if( bx > pGrid->Get_NX() ) { bx = pGrid->Get_NX(); }
393
394 for(int iy=ay; iy<by; iy++) for(int ix=ax; ix<bx; ix++)
395 {
396 if( !pGrid->is_NoData(ix, iy) )
397 {
398 s += pGrid->asDouble(ix, iy);
399 }
400 }
401
402 int n; double z;
403
404 if( s.Get_Majority(z, n) )//&& n > 1 )
405 {
406 Set_Value(x, y, z);
407 }
408 }
409 }
410 }
411 }
412
413 //-----------------------------------------------------
414 return( true );
415}
416
417
419// //
420// Operatoren //
421// //
423
424//---------------------------------------------------------
426{
428
429 return( *this );
430}
431
433{
434 Assign(Value);
435
436 return( *this );
437}
438
439//---------------------------------------------------------
441{
442 CSG_Grid g(*this);
443
444 return( g._Operation_Arithmetic(Grid, GRID_OPERATION_Addition) );
445}
446
448{
449 CSG_Grid g(*this);
450
451 return( g._Operation_Arithmetic(Value, GRID_OPERATION_Addition) );
452}
453
455{
456 return( _Operation_Arithmetic(Grid, GRID_OPERATION_Addition) );
457}
458
460{
461 return( _Operation_Arithmetic(Value, GRID_OPERATION_Addition) );
462}
463
465{
466 return( _Operation_Arithmetic(Grid, GRID_OPERATION_Addition) );
467}
468
469CSG_Grid & CSG_Grid::Add (double Value)
470{
471 return( _Operation_Arithmetic(Value, GRID_OPERATION_Addition) );
472}
473
474//---------------------------------------------------------
476{
477 CSG_Grid g(*this);
478
479 return( g._Operation_Arithmetic(Grid, GRID_OPERATION_Subtraction) );
480}
481
483{
484 CSG_Grid g(*this);
485
486 return( g._Operation_Arithmetic(Value, GRID_OPERATION_Subtraction) );
487}
488
490{
491 return( _Operation_Arithmetic(Grid, GRID_OPERATION_Subtraction) );
492}
493
495{
496 return( _Operation_Arithmetic(Value, GRID_OPERATION_Subtraction) );
497}
498
500{
501 return( _Operation_Arithmetic(Grid, GRID_OPERATION_Subtraction) );
502}
503
505{
506 return( _Operation_Arithmetic(Value, GRID_OPERATION_Subtraction) );
507}
508
509//---------------------------------------------------------
511{
512 CSG_Grid g(*this);
513
514 return( g._Operation_Arithmetic(Grid, GRID_OPERATION_Multiplication) );
515}
516
518{
519 CSG_Grid g(*this);
520
521 return( g._Operation_Arithmetic(Value, GRID_OPERATION_Multiplication) );
522}
523
525{
526 return( _Operation_Arithmetic(Grid, GRID_OPERATION_Multiplication) );
527}
528
530{
531 return( _Operation_Arithmetic(Value, GRID_OPERATION_Multiplication) );
532}
533
535{
536 return( _Operation_Arithmetic(Grid, GRID_OPERATION_Multiplication) );
537}
538
540{
541 return( _Operation_Arithmetic(Value, GRID_OPERATION_Multiplication) );
542}
543
544//---------------------------------------------------------
546{
547 CSG_Grid g(*this);
548
549 return( g._Operation_Arithmetic(Grid, GRID_OPERATION_Division) );
550}
551
553{
554 CSG_Grid g(*this);
555
556 return( g._Operation_Arithmetic(Value, GRID_OPERATION_Division) );
557}
558
560{
561 return( _Operation_Arithmetic(Grid, GRID_OPERATION_Division) );
562}
563
565{
566 return( _Operation_Arithmetic(Value, GRID_OPERATION_Division) );
567}
568
570{
571 return( _Operation_Arithmetic(Grid, GRID_OPERATION_Division) );
572}
573
575{
576 return( _Operation_Arithmetic(Value, GRID_OPERATION_Division) );
577}
578
579
581// //
582// Operatoren //
583// //
585
586//---------------------------------------------------------
587CSG_Grid & CSG_Grid::_Operation_Arithmetic(const CSG_Grid &Grid, TSG_Grid_Operation Operation)
588{
589 if( is_Intersecting(Grid.Get_Extent()) )
590 {
591 CSG_Grid_Resampling Interpolation =
592 Get_Cellsize() == Grid.Get_Cellsize() && fmod(Get_XMin() - Grid.Get_XMin(), Get_Cellsize()) == 0.
593 && Get_Cellsize() == Grid.Get_Cellsize() && fmod(Get_YMin() - Grid.Get_YMin(), Get_Cellsize()) == 0.
596
597 #pragma omp parallel for
598 for(int y=0; y<Get_NY(); y++)
599 {
600 double yWorld = Get_YMin() + y * Get_Cellsize();
601
602 for(int x=0; x<Get_NX(); x++)
603 {
604 if( !is_NoData(x, y) )
605 {
606 double xWorld = Get_XMin() + x * Get_Cellsize(), Value;
607
608 if( Grid.Get_Value(xWorld, yWorld, Value, Interpolation) )
609 {
610 switch( Operation )
611 {
612 case GRID_OPERATION_Addition : Add_Value(x, y, Value); break;
613 case GRID_OPERATION_Subtraction : Add_Value(x, y, -Value); break;
614 case GRID_OPERATION_Multiplication: Mul_Value(x, y, Value); break;
616 if( Value != 0. )
617 {
618 Mul_Value(x, y, 1. / Value);
619 }
620 else
621 {
622 Set_NoData(x, y);
623 }
624 break;
625 }
626 }
627 }
628 }
629 }
630 }
631
632 return( *this );
633}
634
635//---------------------------------------------------------
636CSG_Grid & CSG_Grid::_Operation_Arithmetic(double Value, TSG_Grid_Operation Operation)
637{
638 switch( Operation )
639 {
641 if( Value == 0. )
642 return( *this );
643 break;
644
646 if( Value == 0. )
647 return( *this );
648 Value = -Value;
649 break;
650
652 if( Value == 1. )
653 return( *this );
654 break;
655
657 if( Value == 0. )
658 return( *this );
659 Value = 1. / Value;
660 break;
661 }
662
663 //-----------------------------------------------------
664 #pragma omp parallel for
665 for(int y=0; y<Get_NY(); y++)
666 {
667 for(int x=0; x<Get_NX(); x++)
668 {
669 if( !is_NoData(x, y) )
670 {
671 switch( Operation )
672 {
674 case GRID_OPERATION_Subtraction : Add_Value(x, y, Value); break;
675
677 case GRID_OPERATION_Division : Mul_Value(x, y, Value); break;
678 }
679 }
680 }
681 }
682
683 //-----------------------------------------------------
684 return( *this );
685}
686
687
689// //
690// Grid-Operations - A //
691// //
693
694//---------------------------------------------------------
696{
697 if( is_Valid() && Get_Range() > 0. )
698 {
699 double Min = Get_Min();
700 double Max = Get_Max();
701
702 #pragma omp parallel for
703 for(int y=0; y<Get_NY(); y++)
704 {
705 for(int x=0; x<Get_NX(); x++)
706 {
707 if( !is_NoData(x, y) )
708 {
709 Set_Value(x, y, Max - (asDouble(x, y) - Min));
710 }
711 }
712 }
713 }
714}
715
716//---------------------------------------------------------
718{
719 if( is_Valid() )
720 {
721 #pragma omp parallel for
722 for(int x=0; x<Get_NX(); x++)
723 {
724 for(int yA=0, yB=Get_NY()-1; yA<yB; yA++, yB--)
725 {
726 double d = asDouble(x, yA);
727 Set_Value(x, yA, asDouble(x, yB));
728 Set_Value(x, yB, d);
729 }
730 }
731 }
732}
733
734//---------------------------------------------------------
736{
737 if( is_Valid() )
738 {
739 #pragma omp parallel for
740 for(int y=0; y<Get_NY(); y++)
741 {
742 for(int xA=0, xB=Get_NX()-1; xA<xB; xA++, xB--)
743 {
744 double d = asDouble(xA, y);
745 Set_Value(xA, y, asDouble(xB, y));
746 Set_Value(xB, y, d);
747 }
748 }
749 }
750}
751
752
754// //
755// Grid-Operations - B //
756// //
758
759//---------------------------------------------------------
761{
762 if( is_Valid() && Get_Range() > 0. )
763 {
764 double Min = Get_Min ();
765 double Range = Get_Range();
766
767 #pragma omp parallel for
768 for(int y=0; y<Get_NY(); y++)
769 {
770 for(int x=0; x<Get_NX(); x++)
771 {
772 if( !is_NoData(x, y) )
773 {
774 Set_Value(x, y, (asDouble(x, y) - Min) / Range);
775 }
776 }
777 }
778
779 return( true );
780 }
781
782 return( false );
783}
784
785//---------------------------------------------------------
787{
788 if( is_Valid() && Minimum < Maximum )
789 {
790 #pragma omp parallel for
791 for(int y=0; y<Get_NY(); y++)
792 {
793 for(int x=0; x<Get_NX(); x++)
794 {
795 if( !is_NoData(x, y) )
796 {
797 Set_Value(x, y, Minimum + asDouble(x, y) * (Maximum - Minimum));
798 }
799 }
800 }
801
802 return( true );
803 }
804
805 return( false );
806}
807
808//---------------------------------------------------------
810{
811 if( is_Valid() && Get_StdDev() > 0. )
812 {
813 double Mean = Get_Mean ();
814 double StdDev = Get_StdDev();
815
816 #pragma omp parallel for
817 for(int y=0; y<Get_NY(); y++)
818 {
819 for(int x=0; x<Get_NX(); x++)
820 {
821 if( !is_NoData(x, y) )
822 {
823 Set_Value(x, y, (asDouble(x, y) - Mean) / StdDev);
824 }
825 }
826 }
827
828 return( true );
829 }
830
831 return( false );
832}
833
834//---------------------------------------------------------
835bool CSG_Grid::DeStandardise(double Mean, double StdDev)
836{
837 if( is_Valid() && StdDev > 0. )
838 {
839 #pragma omp parallel for
840 for(int y=0; y<Get_NY(); y++)
841 {
842 for(int x=0; x<Get_NX(); x++)
843 {
844 if( !is_NoData(x, y) )
845 {
846 Set_Value(x, y, Mean + asDouble(x, y) * StdDev);
847 }
848 }
849 }
850
851 return( true );
852 }
853
854 return( false );
855}
856
857
859// //
860// //
861// //
863
864//---------------------------------------------------------
876//---------------------------------------------------------
877int CSG_Grid::Get_Gradient_NeighborDir(int x, int y, bool bDown, bool bNoEdges) const
878{
879 int Direction = -1;
880
881 if( is_InGrid(x, y) )
882 {
883 double z = asDouble(x, y), dzMax = 0.;
884
885 for(int i=0; i<8; i++)
886 {
887 int ix = m_System.Get_xTo(i, x);
888 int iy = m_System.Get_yTo(i, y);
889
890 if( is_InGrid(ix, iy) )
891 {
892 double dz = (z - asDouble(ix, iy)) / m_System.Get_Length(i);
893
894 if( (!bDown || dz > 0.) && (Direction < 0 || dzMax < dz) )
895 {
896 dzMax = dz;
897 Direction = i;
898 }
899 }
900 else if( bNoEdges )
901 {
902 return( -1 );
903 }
904 }
905 }
906
907 return( Direction );
908}
909
910//---------------------------------------------------------
918//---------------------------------------------------------
919bool CSG_Grid::Get_Gradient(int x, int y, double &Slope, double &Aspect) const
920{
921 if( is_InGrid(x, y) )
922 {
923 double z = asDouble(x, y), dz[4];
924
925 for(int i=0, iDir=0, ix, iy; i<4; i++, iDir+=2)
926 {
927 if( is_InGrid(
928 ix = m_System.Get_xTo (iDir, x),
929 iy = m_System.Get_yTo (iDir, y)) )
930 {
931 dz[i] = asDouble(ix, iy) - z;
932 }
933 else if( is_InGrid(
934 ix = m_System.Get_xFrom(iDir, x),
935 iy = m_System.Get_yFrom(iDir, y)) )
936 {
937 dz[i] = z - asDouble(ix, iy);
938 }
939 else
940 {
941 dz[i] = 0.;
942 }
943 }
944
945 double G = (dz[0] - dz[2]) / (2. * Get_Cellsize());
946 double H = (dz[1] - dz[3]) / (2. * Get_Cellsize());
947
948 Slope = atan(sqrt(G*G + H*H));
949 Aspect = G != 0. ? M_PI_180 + atan2(H, G) : H > 0. ? M_PI_270 : H < 0. ? M_PI_090 : -1.;
950
951 return( true );
952 }
953
954 Slope = 0.;
955 Aspect = -1.;
956
957 return( false );
958}
959
960//---------------------------------------------------------
967//---------------------------------------------------------
968bool CSG_Grid::Get_Gradient(double x, double y, double &Slope, double &Aspect, CSG_Grid_Resampling Interpolation) const
969{
970 double z, iz, dz[4];
971
972 if( Get_Value(x, y, z, Interpolation) )
973 {
974 for(int i=0, iDir=0; i<4; i++, iDir+=2)
975 {
976 if( Get_Value(
977 x + Get_Cellsize() * m_System.Get_xTo (iDir),
978 y + Get_Cellsize() * m_System.Get_yTo (iDir), iz, Interpolation) )
979 {
980 dz[i] = iz - z;
981 }
982 else if( Get_Value(
983 x + Get_Cellsize() * m_System.Get_xFrom(iDir),
984 y + Get_Cellsize() * m_System.Get_yFrom(iDir), iz, Interpolation) )
985 {
986 dz[i] = z - iz;
987 }
988 else
989 {
990 dz[i] = 0.;
991 }
992 }
993
994 double G = (dz[0] - dz[2]) / (2. * Get_Cellsize());
995 double H = (dz[1] - dz[3]) / (2. * Get_Cellsize());
996
997 Slope = atan(sqrt(G*G + H*H));
998 Aspect = G != 0. ? M_PI_180 + atan2(H, G) : H > 0. ? M_PI_270 : H < 0. ? M_PI_090 : -1.;
999
1000 return( true );
1001 }
1002
1003 Slope = 0.;
1004 Aspect = -1.;
1005
1006 return( false );
1007}
1008
1009//---------------------------------------------------------
1016//---------------------------------------------------------
1017bool CSG_Grid::Get_Gradient(const TSG_Point &p, double &Incline, double &Azimuth, CSG_Grid_Resampling Interpolation) const
1018{
1019 return( Get_Gradient(p.x, p.y, Incline, Azimuth, Interpolation) );
1020}
1021
1022
1024// //
1025// //
1026// //
1028
1029//---------------------------------------------------------
bool SG_UI_Process_Set_Ready(void)
bool SG_UI_Process_Set_Progress(int Position, int Range)
void Set_Update_Flag(bool bOn=true)
Definition dataobject.h:285
double Get_NoData_Value(bool bUpper=false) const
Definition dataobject.h:253
CSG_MetaData & Get_History(void)
Definition dataobject.h:236
virtual bool Set_NoData_Value(double Value)
CSG_Projection & Get_Projection(void)
class CSG_Grid * asGrid(bool bPolymorph=false) const
int Get_NY(void) const
Definition grid.h:564
virtual CSG_Grid & operator*=(const CSG_Grid &Grid)
double Get_Max(void)
Definition grid.cpp:1027
bool DeNormalise(double Minimum, double Maximum)
double Get_YMin(bool bCells=false) const
Definition grid.h:576
bool is_Cached(void) const
Definition grid.h:633
void Invert(void)
virtual CSG_Grid & operator+=(const CSG_Grid &Grid)
virtual CSG_Grid & Subtract(const CSG_Grid &Grid)
double Get_Mean(void)
Definition grid.cpp:1017
TSG_Intersection is_Intersecting(const CSG_Rect &Extent) const
Definition grid.cpp:452
bool Standardise(void)
double Get_XMin(bool bCells=false) const
Definition grid.h:572
void Mirror(void)
virtual void Set_NoData(int x, int y)
Definition grid.h:730
virtual CSG_Grid operator-(const CSG_Grid &Grid) const
bool DeStandardise(double Mean, double StdDev)
virtual CSG_Grid & operator=(const CSG_Grid &Grid)
virtual CSG_Grid & Add(const CSG_Grid &Grid)
double Get_Cellsize(void) const
Definition grid.h:567
virtual void Mul_Value(int x, int y, double Value)
Definition grid.h:839
virtual bool is_NoData(int x, int y) const
Definition grid.h:727
int Get_Gradient_NeighborDir(int x, int y, bool bDown=true, bool bNoEdges=true) const
virtual CSG_Grid operator*(const CSG_Grid &Grid) const
bool Normalise(void)
CSG_Grid(void)
Definition grid.cpp:136
virtual bool is_Valid(void) const
Definition grid.cpp:441
int Get_nLineBytes(void) const
Definition grid.h:550
const CSG_String & Get_Unit(void) const
Definition grid.h:553
virtual CSG_Grid & Divide(const CSG_Grid &Grid)
virtual double asDouble(sLong i, bool bScaled=true) const
Definition grid.h:793
virtual void Add_Value(int x, int y, double Value)
Definition grid.h:836
void Flip(void)
virtual CSG_Grid operator/(const CSG_Grid &Grid) const
double Get_Min(void)
Definition grid.cpp:1022
virtual CSG_Grid operator+(const CSG_Grid &Grid) const
virtual bool Assign(double Value=0.)
virtual CSG_Grid & Multiply(const CSG_Grid &Grid)
double Get_Value(double x, double y, CSG_Grid_Resampling Resampling=CSG_Grid_Resampling::Bicubic_2, bool bByteWise=false) const
Definition grid.cpp:539
void Set_Unit(const CSG_String &Unit)
Definition grid.cpp:400
virtual const CSG_Rect & Get_Extent(void)
Definition grid.h:561
int Get_NX(void) const
Definition grid.h:563
bool is_InGrid(int x, int y, bool bCheckNoData=true) const
Definition grid.h:621
double Get_StdDev(void)
Definition grid.cpp:1037
virtual void Set_Value(sLong i, double Value, bool bScaled=true)
Definition grid.h:843
virtual CSG_Grid & operator/=(const CSG_Grid &Grid)
bool Assign_NoData(void)
virtual CSG_Grid & operator-=(const CSG_Grid &Grid)
double Get_Range(void)
Definition grid.cpp:1032
bool Get_Gradient(int x, int y, double &Slope, double &Aspect) const
void Destroy(void)
Definition metadata.cpp:140
bool is_Okay(void) const
Definition geo_tools.h:863
double Get_Mean(void)
Definition mat_tools.h:753
sLong Get_Count(void) const
Definition mat_tools.h:745
void Add_Value(double Value, double Weight=1.)
virtual bool Get_Majority(double &Value) const
Definition mat_tools.h:870
struct SSG_Point TSG_Point
@ INTERSECTION_None
Definition geo_tools.h:102
TSG_Grid_Operation
Definition grid.h:194
@ GRID_OPERATION_Multiplication
Definition grid.h:197
@ GRID_OPERATION_Division
Definition grid.h:198
@ GRID_OPERATION_Addition
Definition grid.h:195
@ GRID_OPERATION_Subtraction
Definition grid.h:196
CSG_Grid_Resampling
Definition grid.h:156
#define M_PI_090
Definition mat_tools.h:100
#define M_PI_180
Definition mat_tools.h:102
#define M_PI_270
Definition mat_tools.h:104
double x
Definition geo_tools.h:129
double y
Definition geo_tools.h:129