SAGA API  v9.8
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 //---------------------------------------------------------
65 bool CSG_Grid::Assign_NoData(void)
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 //---------------------------------------------------------
79 bool 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  //-----------------------------------------------------
104  Get_History().Destroy();
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 //---------------------------------------------------------
121 bool CSG_Grid::Assign(CSG_Data_Object *pObject, bool bProgress)
122 {
123  return( pObject && Assign(pObject->asGrid(), GRID_RESAMPLING_Undefined) );
124 }
125 
126 bool CSG_Grid::Assign(CSG_Grid *pGrid, TSG_Grid_Resampling Interpolation, 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, GRID_RESAMPLING_NearestNeighbour, bProgress);
160  }
161 
162  //---------------------------------------------------------
163  else switch( Interpolation )
164  {
165  case GRID_RESAMPLING_NearestNeighbour:
166  case GRID_RESAMPLING_Bilinear :
167  case GRID_RESAMPLING_BicubicSpline :
168  case GRID_RESAMPLING_BSpline : bResult = _Assign_Interpolated(pGrid, Interpolation , bProgress); break;
169 
170  case GRID_RESAMPLING_Mean_Nodes :
171  case GRID_RESAMPLING_Mean_Cells : bResult = _Assign_MeanValue (pGrid, Interpolation != GRID_RESAMPLING_Mean_Nodes, bProgress); break;
172 
173  case GRID_RESAMPLING_Minimum :
174  case GRID_RESAMPLING_Maximum : bResult = _Assign_ExtremeValue(pGrid, Interpolation == GRID_RESAMPLING_Maximum , bProgress); break;
175 
176  case GRID_RESAMPLING_Majority : bResult = _Assign_Majority (pGrid , bProgress); break;
177 
178  default : bResult = Get_Cellsize() < pGrid->Get_Cellsize()
179  ? _Assign_Interpolated(pGrid, GRID_RESAMPLING_BSpline, 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  {
200  SG_UI_Process_Set_Ready();
201  }
202 
203  return( bResult );
204 }
205 
206 //---------------------------------------------------------
207 bool CSG_Grid::_Assign_Interpolated(CSG_Grid *pGrid, TSG_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 //---------------------------------------------------------
233 bool 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  //-----------------------------------------------------
241  Assign_NoData();
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 //---------------------------------------------------------
282 bool 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 //---------------------------------------------------------
359 bool CSG_Grid::_Assign_Majority(CSG_Grid *pGrid, bool bProgress)
360 {
361  if( Get_Cellsize() < pGrid->Get_Cellsize() )
362  {
363  return( false );
364  }
365 
366  //-----------------------------------------------------
367  Set_NoData_Value(pGrid->Get_NoData_Value());
368 
369  Assign_NoData();
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 //---------------------------------------------------------
425 CSG_Grid & CSG_Grid::operator = (const CSG_Grid &Grid)
426 {
427  Assign(Grid.asGrid(), GRID_RESAMPLING_Undefined, false);
428 
429  return( *this );
430 }
431 
432 CSG_Grid & CSG_Grid::operator = (double Value)
433 {
434  Assign(Value);
435 
436  return( *this );
437 }
438 
439 //---------------------------------------------------------
440 CSG_Grid CSG_Grid::operator + (const CSG_Grid &Grid) const
441 {
442  CSG_Grid g(*this);
443 
444  return( g._Operation_Arithmetic(Grid, GRID_OPERATION_Addition) );
445 }
446 
447 CSG_Grid CSG_Grid::operator + (double Value) const
448 {
449  CSG_Grid g(*this);
450 
451  return( g._Operation_Arithmetic(Value, GRID_OPERATION_Addition) );
452 }
453 
454 CSG_Grid & CSG_Grid::operator += (const CSG_Grid &Grid)
455 {
456  return( _Operation_Arithmetic(Grid, GRID_OPERATION_Addition) );
457 }
458 
459 CSG_Grid & CSG_Grid::operator += (double Value)
460 {
461  return( _Operation_Arithmetic(Value, GRID_OPERATION_Addition) );
462 }
463 
464 CSG_Grid & CSG_Grid::Add (const CSG_Grid &Grid)
465 {
466  return( _Operation_Arithmetic(Grid, GRID_OPERATION_Addition) );
467 }
468 
469 CSG_Grid & CSG_Grid::Add (double Value)
470 {
471  return( _Operation_Arithmetic(Value, GRID_OPERATION_Addition) );
472 }
473 
474 //---------------------------------------------------------
475 CSG_Grid CSG_Grid::operator - (const CSG_Grid &Grid) const
476 {
477  CSG_Grid g(*this);
478 
479  return( g._Operation_Arithmetic(Grid, GRID_OPERATION_Subtraction) );
480 }
481 
482 CSG_Grid CSG_Grid::operator - (double Value) const
483 {
484  CSG_Grid g(*this);
485 
486  return( g._Operation_Arithmetic(Value, GRID_OPERATION_Subtraction) );
487 }
488 
489 CSG_Grid & CSG_Grid::operator -= (const CSG_Grid &Grid)
490 {
491  return( _Operation_Arithmetic(Grid, GRID_OPERATION_Subtraction) );
492 }
493 
494 CSG_Grid & CSG_Grid::operator -= (double Value)
495 {
496  return( _Operation_Arithmetic(Value, GRID_OPERATION_Subtraction) );
497 }
498 
499 CSG_Grid & CSG_Grid::Subtract (const CSG_Grid &Grid)
500 {
501  return( _Operation_Arithmetic(Grid, GRID_OPERATION_Subtraction) );
502 }
503 
504 CSG_Grid & CSG_Grid::Subtract (double Value)
505 {
506  return( _Operation_Arithmetic(Value, GRID_OPERATION_Subtraction) );
507 }
508 
509 //---------------------------------------------------------
510 CSG_Grid CSG_Grid::operator * (const CSG_Grid &Grid) const
511 {
512  CSG_Grid g(*this);
513 
514  return( g._Operation_Arithmetic(Grid, GRID_OPERATION_Multiplication) );
515 }
516 
517 CSG_Grid CSG_Grid::operator * (double Value) const
518 {
519  CSG_Grid g(*this);
520 
521  return( g._Operation_Arithmetic(Value, GRID_OPERATION_Multiplication) );
522 }
523 
524 CSG_Grid & CSG_Grid::operator *= (const CSG_Grid &Grid)
525 {
526  return( _Operation_Arithmetic(Grid, GRID_OPERATION_Multiplication) );
527 }
528 
529 CSG_Grid & CSG_Grid::operator *= (double Value)
530 {
531  return( _Operation_Arithmetic(Value, GRID_OPERATION_Multiplication) );
532 }
533 
534 CSG_Grid & CSG_Grid::Multiply (const CSG_Grid &Grid)
535 {
536  return( _Operation_Arithmetic(Grid, GRID_OPERATION_Multiplication) );
537 }
538 
539 CSG_Grid & CSG_Grid::Multiply (double Value)
540 {
541  return( _Operation_Arithmetic(Value, GRID_OPERATION_Multiplication) );
542 }
543 
544 //---------------------------------------------------------
545 CSG_Grid CSG_Grid::operator / (const CSG_Grid &Grid) const
546 {
547  CSG_Grid g(*this);
548 
549  return( g._Operation_Arithmetic(Grid, GRID_OPERATION_Division) );
550 }
551 
552 CSG_Grid CSG_Grid::operator / (double Value) const
553 {
554  CSG_Grid g(*this);
555 
556  return( g._Operation_Arithmetic(Value, GRID_OPERATION_Division) );
557 }
558 
559 CSG_Grid & CSG_Grid::operator /= (const CSG_Grid &Grid)
560 {
561  return( _Operation_Arithmetic(Grid, GRID_OPERATION_Division) );
562 }
563 
564 CSG_Grid & CSG_Grid::operator /= (double Value)
565 {
566  return( _Operation_Arithmetic(Value, GRID_OPERATION_Division) );
567 }
568 
569 CSG_Grid & CSG_Grid::Divide (const CSG_Grid &Grid)
570 {
571  return( _Operation_Arithmetic(Grid, GRID_OPERATION_Division) );
572 }
573 
574 CSG_Grid & CSG_Grid::Divide (double Value)
575 {
576  return( _Operation_Arithmetic(Value, GRID_OPERATION_Division) );
577 }
578 
579 
581 // //
582 // Operatoren //
583 // //
585 
586 //---------------------------------------------------------
587 CSG_Grid & CSG_Grid::_Operation_Arithmetic(const CSG_Grid &Grid, TSG_Grid_Operation Operation)
588 {
589  if( is_Intersecting(Grid.Get_Extent()) )
590  {
591  TSG_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.
594  ? GRID_RESAMPLING_NearestNeighbour
595  : GRID_RESAMPLING_BSpline;
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;
615  case GRID_OPERATION_Division :
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 //---------------------------------------------------------
636 CSG_Grid & CSG_Grid::_Operation_Arithmetic(double Value, TSG_Grid_Operation Operation)
637 {
638  switch( Operation )
639  {
640  case GRID_OPERATION_Addition :
641  if( Value == 0. )
642  return( *this );
643  break;
644 
645  case GRID_OPERATION_Subtraction :
646  if( Value == 0. )
647  return( *this );
648  Value = -Value;
649  break;
650 
651  case GRID_OPERATION_Multiplication:
652  if( Value == 1. )
653  return( *this );
654  break;
655 
656  case GRID_OPERATION_Division :
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  {
673  case GRID_OPERATION_Addition :
674  case GRID_OPERATION_Subtraction : Add_Value(x, y, Value); break;
675 
676  case GRID_OPERATION_Multiplication:
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 //---------------------------------------------------------
695 void CSG_Grid::Invert(void)
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 //---------------------------------------------------------
717 void CSG_Grid::Flip(void)
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 //---------------------------------------------------------
735 void CSG_Grid::Mirror(void)
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 //---------------------------------------------------------
760 bool CSG_Grid::Normalise(void)
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 //---------------------------------------------------------
786 bool CSG_Grid::DeNormalise(double Minimum, double Maximum)
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 //---------------------------------------------------------
809 bool CSG_Grid::Standardise(void)
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 //---------------------------------------------------------
835 bool 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 //---------------------------------------------------------
877 int 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 //---------------------------------------------------------
919 bool 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 //---------------------------------------------------------
968 bool CSG_Grid::Get_Gradient(double x, double y, double &Slope, double &Aspect, TSG_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 //---------------------------------------------------------
1017 bool CSG_Grid::Get_Gradient(const TSG_Point &p, double &Incline, double &Azimuth, TSG_Grid_Resampling Interpolation) const
1018 {
1019  return( Get_Gradient(p.x, p.y, Incline, Azimuth, Interpolation) );
1020 }
1021 
1022 
1024 // //
1025 // //
1026 // //
1028 
1029 //---------------------------------------------------------
CSG_MetaData::Destroy
void Destroy(void)
Definition: metadata.cpp:140
CSG_Rect
Definition: geo_tools.h:474
CSG_Grid::DeStandardise
bool DeStandardise(double Mean, double StdDev)
Definition: grid_operation.cpp:835
CSG_Grid::is_InGrid
bool is_InGrid(int x, int y, bool bCheckNoData=true) const
Definition: grid.h:601
CSG_Grid::Add_Value
virtual void Add_Value(int x, int y, double Value)
Definition: grid.h:816
CSG_Grid_System::Get_xTo
static int Get_xTo(int Direction, int x=0)
Definition: grid.h:318
GRID_RESAMPLING_NearestNeighbour
@ GRID_RESAMPLING_NearestNeighbour
Definition: grid.h:157
CSG_Data_Object::Get_History
CSG_MetaData & Get_History(void)
Definition: dataobject.h:236
GRID_RESAMPLING_BicubicSpline
@ GRID_RESAMPLING_BicubicSpline
Definition: grid.h:159
CSG_Grid::operator/=
virtual CSG_Grid & operator/=(const CSG_Grid &Grid)
Definition: grid_operation.cpp:559
CSG_Grid::Get_nLineBytes
int Get_nLineBytes(void) const
Definition: grid.h:530
CSG_Grid::is_Intersecting
TSG_Intersection is_Intersecting(const CSG_Rect &Extent) const
Definition: grid.cpp:452
CSG_Data_Object::Get_NoData_Value
double Get_NoData_Value(bool bUpper=false) const
Definition: dataobject.h:253
CSG_Grid::Normalise
bool Normalise(void)
Definition: grid_operation.cpp:760
CSG_Unique_Number_Statistics
Definition: mat_tools.h:847
TSG_Grid_Resampling
TSG_Grid_Resampling
Definition: grid.h:156
CSG_Grid::Flip
void Flip(void)
Definition: grid_operation.cpp:717
CSG_Grid::operator+=
virtual CSG_Grid & operator+=(const CSG_Grid &Grid)
Definition: grid_operation.cpp:454
GRID_RESAMPLING_Undefined
@ GRID_RESAMPLING_Undefined
Definition: grid.h:168
CSG_Projection::is_Okay
bool is_Okay(void) const
Definition: geo_tools.h:863
GRID_OPERATION_Multiplication
@ GRID_OPERATION_Multiplication
Definition: grid.h:177
CSG_Grid::Assign
virtual bool Assign(double Value=0.)
Definition: grid_operation.cpp:79
grid.h
CSG_Grid_System::Get_yTo
static int Get_yTo(int Direction, int y=0)
Definition: grid.h:332
GRID_RESAMPLING_Maximum
@ GRID_RESAMPLING_Maximum
Definition: grid.h:165
CSG_Grid::operator/
virtual CSG_Grid operator/(const CSG_Grid &Grid) const
Definition: grid_operation.cpp:545
CSG_Grid::Mul_Value
virtual void Mul_Value(int x, int y, double Value)
Definition: grid.h:819
CSG_Unique_Number_Statistics::Get_Majority
virtual bool Get_Majority(double &Value) const
Definition: mat_tools.h:868
GRID_OPERATION_Subtraction
@ GRID_OPERATION_Subtraction
Definition: grid.h:176
SSG_Point
Definition: geo_tools.h:128
CSG_Data_Object::Set_Update_Flag
void Set_Update_Flag(bool bOn=true)
Definition: dataobject.h:285
GRID_OPERATION_Addition
@ GRID_OPERATION_Addition
Definition: grid.h:175
CSG_Grid::Multiply
virtual CSG_Grid & Multiply(const CSG_Grid &Grid)
Definition: grid_operation.cpp:534
CSG_Grid::Add
virtual CSG_Grid & Add(const CSG_Grid &Grid)
Definition: grid_operation.cpp:464
CSG_Grid::Get_YMin
double Get_YMin(bool bCells=false) const
Definition: grid.h:556
CSG_Grid::Get_Cellsize
double Get_Cellsize(void) const
Definition: grid.h:547
CSG_Data_Object
Definition: dataobject.h:180
GRID_RESAMPLING_Bilinear
@ GRID_RESAMPLING_Bilinear
Definition: grid.h:158
CSG_Grid::Get_Mean
double Get_Mean(void)
Definition: grid.cpp:1011
CSG_Grid::Subtract
virtual CSG_Grid & Subtract(const CSG_Grid &Grid)
Definition: grid_operation.cpp:499
CSG_Simple_Statistics::Get_Count
sLong Get_Count(void) const
Definition: mat_tools.h:743
CSG_Grid::is_Valid
virtual bool is_Valid(void) const
Definition: grid.cpp:441
GRID_RESAMPLING_Mean_Cells
@ GRID_RESAMPLING_Mean_Cells
Definition: grid.h:163
INTERSECTION_None
@ INTERSECTION_None
Definition: geo_tools.h:102
CSG_Grid::Get_Max
double Get_Max(void)
Definition: grid.cpp:1021
CSG_Grid::Set_Value
virtual void Set_Value(sLong i, double Value, bool bScaled=true)
Definition: grid.h:823
CSG_Grid::Get_Gradient
bool Get_Gradient(int x, int y, double &Slope, double &Aspect) const
Definition: grid_operation.cpp:919
CSG_Grid::Get_Unit
const CSG_String & Get_Unit(void) const
Definition: grid.h:533
CSG_Grid::Get_NY
int Get_NY(void) const
Definition: grid.h:544
CSG_Simple_Statistics::Invalidate
void Invalidate(void)
Definition: mat_tools.cpp:447
M_PI_180
#define M_PI_180
Definition: mat_tools.h:102
CSG_Grid::Get_XMin
double Get_XMin(bool bCells=false) const
Definition: grid.h:552
CSG_Grid::Standardise
bool Standardise(void)
Definition: grid_operation.cpp:809
CSG_Grid::Get_Value
double Get_Value(double x, double y, TSG_Grid_Resampling Resampling=GRID_RESAMPLING_BSpline, bool bByteWise=false) const
Definition: grid.cpp:539
CSG_Grid::operator-=
virtual CSG_Grid & operator-=(const CSG_Grid &Grid)
Definition: grid_operation.cpp:489
GRID_RESAMPLING_BSpline
@ GRID_RESAMPLING_BSpline
Definition: grid.h:160
CSG_Grid::operator-
virtual CSG_Grid operator-(const CSG_Grid &Grid) const
Definition: grid_operation.cpp:475
CSG_Simple_Statistics::Get_Mean
double Get_Mean(void)
Definition: mat_tools.h:751
CSG_Grid_System::Get_Length
double Get_Length(int Direction) const
Definition: grid.h:357
CSG_Grid::DeNormalise
bool DeNormalise(double Minimum, double Maximum)
Definition: grid_operation.cpp:786
CSG_Grid::operator+
virtual CSG_Grid operator+(const CSG_Grid &Grid) const
Definition: grid_operation.cpp:440
CSG_Grid_System::Get_yFrom
static int Get_yFrom(int Direction, int y=0)
Definition: grid.h:347
CSG_Grid::operator=
virtual CSG_Grid & operator=(const CSG_Grid &Grid)
Definition: grid_operation.cpp:425
M_PI_270
#define M_PI_270
Definition: mat_tools.h:104
CSG_Grid::Set_NoData
virtual void Set_NoData(int x, int y)
Definition: grid.h:710
CSG_Grid::Get_Gradient_NeighborDir
int Get_Gradient_NeighborDir(int x, int y, bool bDown=true, bool bNoEdges=true) const
Definition: grid_operation.cpp:877
GRID_OPERATION_Division
@ GRID_OPERATION_Division
Definition: grid.h:178
CSG_Grid::Get_Range
double Get_Range(void)
Definition: grid.cpp:1026
GRID_RESAMPLING_Majority
@ GRID_RESAMPLING_Majority
Definition: grid.h:166
M_PI_090
#define M_PI_090
Definition: mat_tools.h:100
CSG_Data_Object::asGrid
class CSG_Grid * asGrid(bool bPolymorph=false) const
Definition: dataobject.cpp:550
GRID_RESAMPLING_Mean_Nodes
@ GRID_RESAMPLING_Mean_Nodes
Definition: grid.h:162
CSG_Simple_Statistics::Add_Value
void Add_Value(double Value, double Weight=1.)
Definition: mat_tools.cpp:527
CSG_Grid::asDouble
virtual double asDouble(sLong i, bool bScaled=true) const
Definition: grid.h:773
CSG_Grid::Get_Min
double Get_Min(void)
Definition: grid.cpp:1016
CSG_Grid::Get_NX
int Get_NX(void) const
Definition: grid.h:543
SG_UI_Process_Set_Progress
bool SG_UI_Process_Set_Progress(int Position, int Range)
Definition: api_callback.cpp:255
SSG_Point::x
double x
Definition: geo_tools.h:129
CSG_Grid_System::Get_xFrom
static int Get_xFrom(int Direction, int x=0)
Definition: grid.h:346
SSG_Point::y
double y
Definition: geo_tools.h:129
CSG_Grid::Divide
virtual CSG_Grid & Divide(const CSG_Grid &Grid)
Definition: grid_operation.cpp:569
CSG_Grid::is_Cached
bool is_Cached(void) const
Definition: grid.h:613
CSG_Grid
Definition: grid.h:481
GRID_RESAMPLING_Minimum
@ GRID_RESAMPLING_Minimum
Definition: grid.h:164
CSG_Data_Object::Set_NoData_Value
virtual bool Set_NoData_Value(double Value)
Definition: dataobject.cpp:572
CSG_Grid::operator*
virtual CSG_Grid operator*(const CSG_Grid &Grid) const
Definition: grid_operation.cpp:510
CSG_Grid::Mirror
void Mirror(void)
Definition: grid_operation.cpp:735
CSG_Grid::Invert
void Invert(void)
Definition: grid_operation.cpp:695
CSG_Grid::Assign_NoData
bool Assign_NoData(void)
Definition: grid_operation.cpp:65
SG_UI_Process_Set_Ready
bool SG_UI_Process_Set_Ready(void)
Definition: api_callback.cpp:305
CSG_Grid::Set_Unit
void Set_Unit(const CSG_String &Unit)
Definition: grid.cpp:400
CSG_Grid::operator*=
virtual CSG_Grid & operator*=(const CSG_Grid &Grid)
Definition: grid_operation.cpp:524
CSG_Grid::is_NoData
virtual bool is_NoData(int x, int y) const
Definition: grid.h:707
CSG_Grid::Get_Extent
virtual const CSG_Rect & Get_Extent(void)
Definition: grid.h:541
CSG_Simple_Statistics
Definition: mat_tools.h:723
CSG_Data_Object::Get_Projection
CSG_Projection & Get_Projection(void)
Definition: dataobject.cpp:637
CSG_Grid::Get_StdDev
double Get_StdDev(void)
Definition: grid.cpp:1031
TSG_Grid_Operation
TSG_Grid_Operation
Definition: grid.h:174