SAGA API  v9.6
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 void 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 
76 //---------------------------------------------------------
77 bool CSG_Grid::Assign(double Value)
78 {
79  if( !is_Valid() )
80  {
81  return( false );
82  }
83 
84  if( Value == 0. && !is_Cached() )
85  {
86  #pragma omp parallel for
87  for(int y=0; y<Get_NY(); y++)
88  {
89  memset(m_Values[y], 0, Get_nLineBytes());
90  }
91  }
92  else
93  {
94  #pragma omp parallel for
95  for(int y=0; y<Get_NY(); y++)
96  {
97  for(int x=0; x<Get_NX(); x++)
98  {
99  Set_Value(x, y, Value);
100  }
101  }
102  }
103 
104  //-----------------------------------------------------
105  Get_History().Destroy();
106 
107  m_Statistics.Invalidate();
108 
109  Set_Update_Flag(false);
110 
111  return( true );
112 }
113 
114 
116 // //
117 // Data Assignments - Grid //
118 // //
120 
121 //---------------------------------------------------------
122 bool CSG_Grid::Assign(CSG_Data_Object *pObject)
123 {
124  return( pObject && pObject->is_Valid() && pObject->Get_ObjectType() == Get_ObjectType()
125  && Assign((CSG_Grid *)pObject, GRID_RESAMPLING_Undefined) );
126 }
127 
128 bool CSG_Grid::Assign(CSG_Grid *pGrid, TSG_Grid_Resampling Interpolation)
129 {
130  //-----------------------------------------------------
131  if( !is_Valid() || !pGrid || !pGrid->is_Valid() || is_Intersecting(pGrid->Get_Extent()) == INTERSECTION_None )
132  {
133  return( false );
134  }
135 
136  bool bResult = false;
137 
138  //---------------------------------------------------------
139  if( m_System == pGrid->m_System )
140  {
141  for(int y=0; y<Get_NY() && SG_UI_Process_Set_Progress(y, Get_NY()); y++)
142  {
143  #pragma omp parallel for
144  for(int x=0; x<Get_NX(); x++)
145  {
146  if( pGrid->is_NoData(x, y) )
147  {
148  Set_NoData(x, y);
149  }
150  else
151  {
152  Set_Value(x, y, pGrid->asDouble(x, y));
153  }
154  }
155  }
156 
157  bResult = true;
158  }
159 
160  else if( Get_Cellsize() == pGrid->Get_Cellsize() // No-Scaling...
161  && fmod(Get_XMin() - pGrid->Get_XMin(), Get_Cellsize()) == 0.
162  && fmod(Get_YMin() - pGrid->Get_YMin(), Get_Cellsize()) == 0. )
163  {
164  bResult = _Assign_Interpolated(pGrid, GRID_RESAMPLING_NearestNeighbour);
165  }
166 
167  //---------------------------------------------------------
168  else switch( Interpolation )
169  {
170  case GRID_RESAMPLING_NearestNeighbour:
171  case GRID_RESAMPLING_Bilinear:
172  case GRID_RESAMPLING_BicubicSpline:
173  case GRID_RESAMPLING_BSpline:
174  bResult = _Assign_Interpolated(pGrid, Interpolation);
175  break;
176 
177  case GRID_RESAMPLING_Mean_Nodes:
178  case GRID_RESAMPLING_Mean_Cells:
179  bResult = _Assign_MeanValue (pGrid, Interpolation != GRID_RESAMPLING_Mean_Nodes);
180  break;
181 
182  case GRID_RESAMPLING_Minimum:
183  case GRID_RESAMPLING_Maximum:
184  bResult = _Assign_ExtremeValue(pGrid, Interpolation == GRID_RESAMPLING_Maximum);
185  break;
186 
187  case GRID_RESAMPLING_Majority:
188  bResult = _Assign_Majority (pGrid);
189  break;
190 
191  default:
192  if( Get_Cellsize() < pGrid->Get_Cellsize() ) // Down-Scaling...
193  {
194  bResult = _Assign_Interpolated(pGrid, GRID_RESAMPLING_BSpline);
195  }
196  else // Up-Scaling...
197  {
198  bResult = _Assign_MeanValue(pGrid, Interpolation != GRID_RESAMPLING_Mean_Nodes);
199  }
200  break;
201  }
202 
203  //---------------------------------------------------------
204  if( bResult )
205  {
206  Set_Unit(pGrid->Get_Unit());
207 
208  if( pGrid->Get_Projection().is_Okay() )
209  {
210  Get_Projection() = pGrid->Get_Projection();
211  }
212 
213  Get_History() = pGrid->Get_History();
214  }
215 
216  //---------------------------------------------------------
217  SG_UI_Process_Set_Ready();
218 
219  return( bResult );
220 }
221 
222 //---------------------------------------------------------
223 bool CSG_Grid::_Assign_Interpolated(CSG_Grid *pGrid, TSG_Grid_Resampling Interpolation)
224 {
225  double py = Get_YMin();
226 
227  for(int y=0; y<Get_NY() && SG_UI_Process_Set_Progress(y, Get_NY()); y++, py+=Get_Cellsize())
228  {
229  #pragma omp parallel for
230  for(int x=0; x<Get_NX(); x++)
231  {
232  double z;
233 
234  if( pGrid->Get_Value(Get_XMin() + x * Get_Cellsize(), py, z, Interpolation) )
235  {
236  Set_Value(x, y, z);
237  }
238  else
239  {
240  Set_NoData(x, y);
241  }
242  }
243  }
244 
245  return( true );
246 }
247 
248 //---------------------------------------------------------
249 bool CSG_Grid::_Assign_ExtremeValue(CSG_Grid *pGrid, bool bMaximum)
250 {
251  if( Get_Cellsize() < pGrid->Get_Cellsize() )
252  {
253  return( false );
254  }
255 
256  //-----------------------------------------------------
257  Assign_NoData();
258 
259  double ax = 0.5 + (pGrid->Get_XMin() - Get_XMin()) / Get_Cellsize();
260  double py = 0.5 + (pGrid->Get_YMin() - Get_YMin()) / Get_Cellsize();
261 
262  double d = pGrid->Get_Cellsize() / Get_Cellsize();
263 
264  for(int y=0; y<pGrid->Get_NY() && SG_UI_Process_Set_Progress(y, pGrid->Get_NY()); y++, py+=d)
265  {
266  int iy = (int)floor(py);
267 
268  if( iy >= 0 && iy < Get_NY() )
269  {
270  #pragma omp parallel for
271  for(int x=0; x<pGrid->Get_NX(); x++)
272  {
273  if( !pGrid->is_NoData(x, y) )
274  {
275  int ix = (int)floor(ax + x * d);
276 
277  if( ix >= 0 && ix < Get_NX() )
278  {
279  double z = pGrid->asDouble(x, y);
280 
281  if( is_NoData(ix, iy)
282  || (bMaximum == true && z > asDouble(ix, iy))
283  || (bMaximum == false && z < asDouble(ix, iy)) )
284  {
285  Set_Value(ix, iy, z);
286  }
287  }
288  }
289  }
290  }
291  }
292 
293  //-----------------------------------------------------
294  return( true );
295 }
296 
297 //---------------------------------------------------------
298 bool CSG_Grid::_Assign_MeanValue(CSG_Grid *pGrid, bool bAreaProportional)
299 {
300  if( Get_Cellsize() < pGrid->Get_Cellsize() )
301  {
302  return( false );
303  }
304 
305  //-----------------------------------------------------
306  double d = Get_Cellsize() / pGrid->Get_Cellsize();
307 
308  double ox = (Get_XMin(true) - pGrid->Get_XMin()) / pGrid->Get_Cellsize();
309  double py = (Get_YMin(true) - pGrid->Get_YMin()) / pGrid->Get_Cellsize();
310 
311  for(int y=0; y<Get_NY() && SG_UI_Process_Set_Progress(y, Get_NY()); y++, py+=d)
312  {
313  int ay = (int)(bAreaProportional ? floor(py ) : ceil (py ));
314  int by = (int)(bAreaProportional ? ceil (py + d) : floor(py + d));
315 
316  //-------------------------------------------------
317  #pragma omp parallel for
318  for(int x=0; x<Get_NX(); x++)
319  {
320  double px = ox + x * d;
321 
322  int ax = (int)(bAreaProportional ? floor(px ) : ceil (px ));
323  int bx = (int)(bAreaProportional ? ceil (px + d) : floor(px + d));
324 
325  CSG_Rect rMean(px, py, px + d, py + d);
326 
327  CSG_Simple_Statistics s;
328 
329  for(int iy=ay; iy<=by; iy++)
330  {
331  if( iy >= 0 && iy < pGrid->Get_NY() )
332  {
333  for(int ix=ax; ix<=bx; ix++)
334  {
335  if( ix >= 0 && ix < pGrid->Get_NX() && !pGrid->is_NoData(ix, iy) )
336  {
337  if( bAreaProportional )
338  {
339  CSG_Rect r(ix - 0.5, iy - 0.5, ix + 0.5, iy + 0.5);
340 
341  if( r.Intersect(rMean) )
342  {
343  s.Add_Value(pGrid->asDouble(ix, iy), r.Get_Area());
344  }
345  }
346  else
347  {
348  s.Add_Value(pGrid->asDouble(ix, iy));
349  }
350  }
351  }
352  }
353  }
354 
355  //---------------------------------------------
356  if( s.Get_Count() > 0 )
357  {
358  Set_Value(x, y, s.Get_Mean());
359  }
360  else
361  {
362  Set_NoData(x, y);
363  }
364  }
365  }
366 
367  //-----------------------------------------------------
368  return( true );
369 }
370 
371 
373 // //
375 
376 //---------------------------------------------------------
377 bool CSG_Grid::_Assign_Majority(CSG_Grid *pGrid)
378 {
379  if( Get_Cellsize() < pGrid->Get_Cellsize() )
380  {
381  return( false );
382  }
383 
384  //-----------------------------------------------------
385  Set_NoData_Value(pGrid->Get_NoData_Value());
386 
387  Assign_NoData();
388 
389  //-----------------------------------------------------
390  int ay, by = (int)(1. + (((0 - 0.5) * Get_Cellsize() + Get_YMin()) - pGrid->Get_YMin()) / pGrid->Get_Cellsize());
391 
392  for(int y=0; y<Get_NY() && SG_UI_Process_Set_Progress(y, Get_NY()); y++)
393  {
394  ay = by;
395  by = (int)(1. + (((y + 0.5) * Get_Cellsize() + Get_YMin()) - pGrid->Get_YMin()) / pGrid->Get_Cellsize());
396 
397  if( ay < pGrid->Get_NY() && by > 0 )
398  {
399  if( ay < 0 )
400  {
401  ay = 0;
402  }
403 
404  if( by > pGrid->Get_NY() )
405  {
406  by = pGrid->Get_NY();
407  }
408 
409  int ax, bx = (int)(1. + (((0 - 0.5) * Get_Cellsize() + Get_XMin()) - pGrid->Get_XMin()) / pGrid->Get_Cellsize());
410 
411  for(int x=0; x<Get_NX(); x++)
412  {
413  ax = bx;
414  bx = (int)(1. + (((x + 0.5) * Get_Cellsize() + Get_XMin()) - pGrid->Get_XMin()) / pGrid->Get_Cellsize());
415 
416  if( ax < pGrid->Get_NX() && bx > 0 )
417  {
418  CSG_Unique_Number_Statistics s;
419 
420  if( ax < 0 )
421  {
422  ax = 0;
423  }
424 
425  if( bx > pGrid->Get_NX() )
426  {
427  bx = pGrid->Get_NX();
428  }
429 
430  for(int iy=ay; iy<by; iy++)
431  {
432  for(int ix=ax; ix<bx; ix++)
433  {
434  if( !pGrid->is_NoData(ix, iy) )
435  {
436  s += pGrid->asDouble(ix, iy);
437  }
438  }
439  }
440 
441  int n; double z;
442 
443  if( s.Get_Majority(z, n) )//&& n > 1 )
444  {
445  Set_Value(x, y, z);
446  }
447  }
448  }
449  }
450  }
451 
452  //-----------------------------------------------------
453  return( true );
454 }
455 
456 
458 // //
459 // Operatoren //
460 // //
462 
463 //---------------------------------------------------------
464 CSG_Grid & CSG_Grid::operator = (const CSG_Grid &Grid)
465 {
466  Assign((CSG_Grid *)&Grid, GRID_RESAMPLING_Undefined);
467 
468  return( *this );
469 }
470 
471 CSG_Grid & CSG_Grid::operator = (double Value)
472 {
473  Assign(Value);
474 
475  return( *this );
476 }
477 
478 //---------------------------------------------------------
479 CSG_Grid CSG_Grid::operator + (const CSG_Grid &Grid) const
480 {
481  CSG_Grid g(*this);
482 
483  return( g._Operation_Arithmetic(Grid, GRID_OPERATION_Addition) );
484 }
485 
486 CSG_Grid CSG_Grid::operator + (double Value) const
487 {
488  CSG_Grid g(*this);
489 
490  return( g._Operation_Arithmetic(Value, GRID_OPERATION_Addition) );
491 }
492 
493 CSG_Grid & CSG_Grid::operator += (const CSG_Grid &Grid)
494 {
495  return( _Operation_Arithmetic(Grid, GRID_OPERATION_Addition) );
496 }
497 
498 CSG_Grid & CSG_Grid::operator += (double Value)
499 {
500  return( _Operation_Arithmetic(Value, GRID_OPERATION_Addition) );
501 }
502 
503 CSG_Grid & CSG_Grid::Add (const CSG_Grid &Grid)
504 {
505  return( _Operation_Arithmetic(Grid, GRID_OPERATION_Addition) );
506 }
507 
508 CSG_Grid & CSG_Grid::Add (double Value)
509 {
510  return( _Operation_Arithmetic(Value, GRID_OPERATION_Addition) );
511 }
512 
513 //---------------------------------------------------------
514 CSG_Grid CSG_Grid::operator - (const CSG_Grid &Grid) const
515 {
516  CSG_Grid g(*this);
517 
518  return( g._Operation_Arithmetic(Grid, GRID_OPERATION_Subtraction) );
519 }
520 
521 CSG_Grid CSG_Grid::operator - (double Value) const
522 {
523  CSG_Grid g(*this);
524 
525  return( g._Operation_Arithmetic(Value, GRID_OPERATION_Subtraction) );
526 }
527 
528 CSG_Grid & CSG_Grid::operator -= (const CSG_Grid &Grid)
529 {
530  return( _Operation_Arithmetic(Grid, GRID_OPERATION_Subtraction) );
531 }
532 
533 CSG_Grid & CSG_Grid::operator -= (double Value)
534 {
535  return( _Operation_Arithmetic(Value, GRID_OPERATION_Subtraction) );
536 }
537 
538 CSG_Grid & CSG_Grid::Subtract (const CSG_Grid &Grid)
539 {
540  return( _Operation_Arithmetic(Grid, GRID_OPERATION_Subtraction) );
541 }
542 
543 CSG_Grid & CSG_Grid::Subtract (double Value)
544 {
545  return( _Operation_Arithmetic(Value, GRID_OPERATION_Subtraction) );
546 }
547 
548 //---------------------------------------------------------
549 CSG_Grid CSG_Grid::operator * (const CSG_Grid &Grid) const
550 {
551  CSG_Grid g(*this);
552 
553  return( g._Operation_Arithmetic(Grid, GRID_OPERATION_Multiplication) );
554 }
555 
556 CSG_Grid CSG_Grid::operator * (double Value) const
557 {
558  CSG_Grid g(*this);
559 
560  return( g._Operation_Arithmetic(Value, GRID_OPERATION_Multiplication) );
561 }
562 
563 CSG_Grid & CSG_Grid::operator *= (const CSG_Grid &Grid)
564 {
565  return( _Operation_Arithmetic(Grid, GRID_OPERATION_Multiplication) );
566 }
567 
568 CSG_Grid & CSG_Grid::operator *= (double Value)
569 {
570  return( _Operation_Arithmetic(Value, GRID_OPERATION_Multiplication) );
571 }
572 
573 CSG_Grid & CSG_Grid::Multiply (const CSG_Grid &Grid)
574 {
575  return( _Operation_Arithmetic(Grid, GRID_OPERATION_Multiplication) );
576 }
577 
578 CSG_Grid & CSG_Grid::Multiply (double Value)
579 {
580  return( _Operation_Arithmetic(Value, GRID_OPERATION_Multiplication) );
581 }
582 
583 //---------------------------------------------------------
584 CSG_Grid CSG_Grid::operator / (const CSG_Grid &Grid) const
585 {
586  CSG_Grid g(*this);
587 
588  return( g._Operation_Arithmetic(Grid, GRID_OPERATION_Division) );
589 }
590 
591 CSG_Grid CSG_Grid::operator / (double Value) const
592 {
593  CSG_Grid g(*this);
594 
595  return( g._Operation_Arithmetic(Value, GRID_OPERATION_Division) );
596 }
597 
598 CSG_Grid & CSG_Grid::operator /= (const CSG_Grid &Grid)
599 {
600  return( _Operation_Arithmetic(Grid, GRID_OPERATION_Division) );
601 }
602 
603 CSG_Grid & CSG_Grid::operator /= (double Value)
604 {
605  return( _Operation_Arithmetic(Value, GRID_OPERATION_Division) );
606 }
607 
608 CSG_Grid & CSG_Grid::Divide (const CSG_Grid &Grid)
609 {
610  return( _Operation_Arithmetic(Grid, GRID_OPERATION_Division) );
611 }
612 
613 CSG_Grid & CSG_Grid::Divide (double Value)
614 {
615  return( _Operation_Arithmetic(Value, GRID_OPERATION_Division) );
616 }
617 
618 
620 // //
621 // Operatoren //
622 // //
624 
625 //---------------------------------------------------------
626 CSG_Grid & CSG_Grid::_Operation_Arithmetic(const CSG_Grid &Grid, TSG_Grid_Operation Operation)
627 {
628  if( is_Intersecting(Grid.Get_Extent()) )
629  {
630  TSG_Grid_Resampling Interpolation =
631  Get_Cellsize() == Grid.Get_Cellsize() && fmod(Get_XMin() - Grid.Get_XMin(), Get_Cellsize()) == 0.
632  && Get_Cellsize() == Grid.Get_Cellsize() && fmod(Get_YMin() - Grid.Get_YMin(), Get_Cellsize()) == 0.
633  ? GRID_RESAMPLING_NearestNeighbour
634  : GRID_RESAMPLING_BSpline;
635 
636  #pragma omp parallel for
637  for(int y=0; y<Get_NY(); y++)
638  {
639  double yWorld = Get_YMin() + y * Get_Cellsize();
640 
641  for(int x=0; x<Get_NX(); x++)
642  {
643  if( !is_NoData(x, y) )
644  {
645  double xWorld = Get_XMin() + x * Get_Cellsize(), Value;
646 
647  if( Grid.Get_Value(xWorld, yWorld, Value, Interpolation) )
648  {
649  switch( Operation )
650  {
651  case GRID_OPERATION_Addition : Add_Value(x, y, Value); break;
652  case GRID_OPERATION_Subtraction : Add_Value(x, y, -Value); break;
653  case GRID_OPERATION_Multiplication: Mul_Value(x, y, Value); break;
654  case GRID_OPERATION_Division :
655  if( Value != 0. )
656  {
657  Mul_Value(x, y, 1. / Value);
658  }
659  else
660  {
661  Set_NoData(x, y);
662  }
663  break;
664  }
665  }
666  }
667  }
668  }
669  }
670 
671  return( *this );
672 }
673 
674 //---------------------------------------------------------
675 CSG_Grid & CSG_Grid::_Operation_Arithmetic(double Value, TSG_Grid_Operation Operation)
676 {
677  switch( Operation )
678  {
679  case GRID_OPERATION_Addition :
680  if( Value == 0. )
681  return( *this );
682  break;
683 
684  case GRID_OPERATION_Subtraction :
685  if( Value == 0. )
686  return( *this );
687  Value = -Value;
688  break;
689 
690  case GRID_OPERATION_Multiplication:
691  if( Value == 1. )
692  return( *this );
693  break;
694 
695  case GRID_OPERATION_Division :
696  if( Value == 0. )
697  return( *this );
698  Value = 1. / Value;
699  break;
700  }
701 
702  //-----------------------------------------------------
703  #pragma omp parallel for
704  for(int y=0; y<Get_NY(); y++)
705  {
706  for(int x=0; x<Get_NX(); x++)
707  {
708  if( !is_NoData(x, y) )
709  {
710  switch( Operation )
711  {
712  case GRID_OPERATION_Addition :
713  case GRID_OPERATION_Subtraction : Add_Value(x, y, Value); break;
714 
715  case GRID_OPERATION_Multiplication:
716  case GRID_OPERATION_Division : Mul_Value(x, y, Value); break;
717  }
718  }
719  }
720  }
721 
722  //-----------------------------------------------------
723  return( *this );
724 }
725 
726 
728 // //
729 // Grid-Operations - A //
730 // //
732 
733 //---------------------------------------------------------
734 void CSG_Grid::Invert(void)
735 {
736  if( is_Valid() && Get_Range() > 0. )
737  {
738  double Min = Get_Min();
739  double Max = Get_Max();
740 
741  #pragma omp parallel for
742  for(int y=0; y<Get_NY(); y++)
743  {
744  for(int x=0; x<Get_NX(); x++)
745  {
746  if( !is_NoData(x, y) )
747  {
748  Set_Value(x, y, Max - (asDouble(x, y) - Min));
749  }
750  }
751  }
752  }
753 }
754 
755 //---------------------------------------------------------
756 void CSG_Grid::Flip(void)
757 {
758  if( is_Valid() )
759  {
760  #pragma omp parallel for
761  for(int x=0; x<Get_NX(); x++)
762  {
763  for(int yA=0, yB=Get_NY()-1; yA<yB; yA++, yB--)
764  {
765  double d = asDouble(x, yA);
766  Set_Value(x, yA, asDouble(x, yB));
767  Set_Value(x, yB, d);
768  }
769  }
770  }
771 }
772 
773 //---------------------------------------------------------
774 void CSG_Grid::Mirror(void)
775 {
776  if( is_Valid() )
777  {
778  #pragma omp parallel for
779  for(int y=0; y<Get_NY(); y++)
780  {
781  for(int xA=0, xB=Get_NX()-1; xA<xB; xA++, xB--)
782  {
783  double d = asDouble(xA, y);
784  Set_Value(xA, y, asDouble(xB, y));
785  Set_Value(xB, y, d);
786  }
787  }
788  }
789 }
790 
791 
793 // //
794 // Grid-Operations - B //
795 // //
797 
798 //---------------------------------------------------------
799 bool CSG_Grid::Normalise(void)
800 {
801  if( is_Valid() && Get_Range() > 0. )
802  {
803  double Min = Get_Min ();
804  double Range = Get_Range();
805 
806  #pragma omp parallel for
807  for(int y=0; y<Get_NY(); y++)
808  {
809  for(int x=0; x<Get_NX(); x++)
810  {
811  if( !is_NoData(x, y) )
812  {
813  Set_Value(x, y, (asDouble(x, y) - Min) / Range);
814  }
815  }
816  }
817 
818  return( true );
819  }
820 
821  return( false );
822 }
823 
824 //---------------------------------------------------------
825 bool CSG_Grid::DeNormalise(double Minimum, double Maximum)
826 {
827  if( is_Valid() && Minimum < Maximum )
828  {
829  #pragma omp parallel for
830  for(int y=0; y<Get_NY(); y++)
831  {
832  for(int x=0; x<Get_NX(); x++)
833  {
834  if( !is_NoData(x, y) )
835  {
836  Set_Value(x, y, Minimum + asDouble(x, y) * (Maximum - Minimum));
837  }
838  }
839  }
840 
841  return( true );
842  }
843 
844  return( false );
845 }
846 
847 //---------------------------------------------------------
848 bool CSG_Grid::Standardise(void)
849 {
850  if( is_Valid() && Get_StdDev() > 0. )
851  {
852  double Mean = Get_Mean ();
853  double StdDev = Get_StdDev();
854 
855  #pragma omp parallel for
856  for(int y=0; y<Get_NY(); y++)
857  {
858  for(int x=0; x<Get_NX(); x++)
859  {
860  if( !is_NoData(x, y) )
861  {
862  Set_Value(x, y, (asDouble(x, y) - Mean) / StdDev);
863  }
864  }
865  }
866 
867  return( true );
868  }
869 
870  return( false );
871 }
872 
873 //---------------------------------------------------------
874 bool CSG_Grid::DeStandardise(double Mean, double StdDev)
875 {
876  if( is_Valid() && StdDev > 0. )
877  {
878  #pragma omp parallel for
879  for(int y=0; y<Get_NY(); y++)
880  {
881  for(int x=0; x<Get_NX(); x++)
882  {
883  if( !is_NoData(x, y) )
884  {
885  Set_Value(x, y, Mean + asDouble(x, y) * StdDev);
886  }
887  }
888  }
889 
890  return( true );
891  }
892 
893  return( false );
894 }
895 
896 
898 // //
899 // //
900 // //
902 
903 //---------------------------------------------------------
915 //---------------------------------------------------------
916 int CSG_Grid::Get_Gradient_NeighborDir(int x, int y, bool bDown, bool bNoEdges) const
917 {
918  int Direction = -1;
919 
920  if( is_InGrid(x, y) )
921  {
922  double z = asDouble(x, y), dzMax = 0.;
923 
924  for(int i=0; i<8; i++)
925  {
926  int ix = m_System.Get_xTo(i, x);
927  int iy = m_System.Get_yTo(i, y);
928 
929  if( is_InGrid(ix, iy) )
930  {
931  double dz = (z - asDouble(ix, iy)) / m_System.Get_Length(i);
932 
933  if( (!bDown || dz > 0.) && (Direction < 0 || dzMax < dz) )
934  {
935  dzMax = dz;
936  Direction = i;
937  }
938  }
939  else if( bNoEdges )
940  {
941  return( -1 );
942  }
943  }
944  }
945 
946  return( Direction );
947 }
948 
949 //---------------------------------------------------------
957 //---------------------------------------------------------
958 bool CSG_Grid::Get_Gradient(int x, int y, double &Slope, double &Aspect) const
959 {
960  if( is_InGrid(x, y) )
961  {
962  double z = asDouble(x, y), dz[4];
963 
964  for(int i=0, iDir=0, ix, iy; i<4; i++, iDir+=2)
965  {
966  if( is_InGrid(
967  ix = m_System.Get_xTo (iDir, x),
968  iy = m_System.Get_yTo (iDir, y)) )
969  {
970  dz[i] = asDouble(ix, iy) - z;
971  }
972  else if( is_InGrid(
973  ix = m_System.Get_xFrom(iDir, x),
974  iy = m_System.Get_yFrom(iDir, y)) )
975  {
976  dz[i] = z - asDouble(ix, iy);
977  }
978  else
979  {
980  dz[i] = 0.;
981  }
982  }
983 
984  double G = (dz[0] - dz[2]) / (2. * Get_Cellsize());
985  double H = (dz[1] - dz[3]) / (2. * Get_Cellsize());
986 
987  Slope = atan(sqrt(G*G + H*H));
988  Aspect = G != 0. ? M_PI_180 + atan2(H, G) : H > 0. ? M_PI_270 : H < 0. ? M_PI_090 : -1.;
989 
990  return( true );
991  }
992 
993  Slope = 0.;
994  Aspect = -1.;
995 
996  return( false );
997 }
998 
999 //---------------------------------------------------------
1006 //---------------------------------------------------------
1007 bool CSG_Grid::Get_Gradient(double x, double y, double &Slope, double &Aspect, TSG_Grid_Resampling Interpolation) const
1008 {
1009  double z, iz, dz[4];
1010 
1011  if( Get_Value(x, y, z, Interpolation) )
1012  {
1013  for(int i=0, iDir=0; i<4; i++, iDir+=2)
1014  {
1015  if( Get_Value(
1016  x + Get_Cellsize() * m_System.Get_xTo (iDir),
1017  y + Get_Cellsize() * m_System.Get_yTo (iDir), iz, Interpolation) )
1018  {
1019  dz[i] = iz - z;
1020  }
1021  else if( Get_Value(
1022  x + Get_Cellsize() * m_System.Get_xFrom(iDir),
1023  y + Get_Cellsize() * m_System.Get_yFrom(iDir), iz, Interpolation) )
1024  {
1025  dz[i] = z - iz;
1026  }
1027  else
1028  {
1029  dz[i] = 0.;
1030  }
1031  }
1032 
1033  double G = (dz[0] - dz[2]) / (2. * Get_Cellsize());
1034  double H = (dz[1] - dz[3]) / (2. * Get_Cellsize());
1035 
1036  Slope = atan(sqrt(G*G + H*H));
1037  Aspect = G != 0. ? M_PI_180 + atan2(H, G) : H > 0. ? M_PI_270 : H < 0. ? M_PI_090 : -1.;
1038 
1039  return( true );
1040  }
1041 
1042  Slope = 0.;
1043  Aspect = -1.;
1044 
1045  return( false );
1046 }
1047 
1048 //---------------------------------------------------------
1055 //---------------------------------------------------------
1056 bool CSG_Grid::Get_Gradient(const TSG_Point &p, double &Incline, double &Azimuth, TSG_Grid_Resampling Interpolation) const
1057 {
1058  return( Get_Gradient(p.x, p.y, Incline, Azimuth, Interpolation) );
1059 }
1060 
1061 
1063 // //
1064 // //
1065 // //
1067 
1068 //---------------------------------------------------------
CSG_MetaData::Destroy
void Destroy(void)
Definition: metadata.cpp:140
CSG_Rect
Definition: geo_tools.h:471
CSG_Grid::DeStandardise
bool DeStandardise(double Mean, double StdDev)
Definition: grid_operation.cpp:874
CSG_Grid::is_InGrid
bool is_InGrid(int x, int y, bool bCheckNoData=true) const
Definition: grid.h:595
CSG_Grid::Add_Value
virtual void Add_Value(int x, int y, double Value)
Definition: grid.h:810
CSG_Grid_System::Get_xTo
static int Get_xTo(int Direction, int x=0)
Definition: grid.h:312
GRID_RESAMPLING_NearestNeighbour
@ GRID_RESAMPLING_NearestNeighbour
Definition: grid.h:157
CSG_Grid::Get_ObjectType
virtual TSG_Data_Object_Type Get_ObjectType(void) const
Definition: grid.h:516
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:598
CSG_Grid::Get_nLineBytes
int Get_nLineBytes(void) const
Definition: grid.h:524
CSG_Grid::is_Intersecting
TSG_Intersection is_Intersecting(const CSG_Rect &Extent) const
Definition: grid.cpp:447
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:799
CSG_Unique_Number_Statistics
Definition: mat_tools.h:837
TSG_Grid_Resampling
TSG_Grid_Resampling
Definition: grid.h:156
CSG_Grid::Flip
void Flip(void)
Definition: grid_operation.cpp:756
CSG_Grid::operator+=
virtual CSG_Grid & operator+=(const CSG_Grid &Grid)
Definition: grid_operation.cpp:493
GRID_RESAMPLING_Undefined
@ GRID_RESAMPLING_Undefined
Definition: grid.h:168
CSG_Projection::is_Okay
bool is_Okay(void) const
Definition: geo_tools.h:857
GRID_OPERATION_Multiplication
@ GRID_OPERATION_Multiplication
Definition: grid.h:177
grid.h
CSG_Grid_System::Get_yTo
static int Get_yTo(int Direction, int y=0)
Definition: grid.h:326
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:584
CSG_Grid::Mul_Value
virtual void Mul_Value(int x, int y, double Value)
Definition: grid.h:813
CSG_Data_Object::Get_ObjectType
virtual TSG_Data_Object_Type Get_ObjectType(void) const =0
Returns the object type as defined by TSG_Data_Object_Type. Used for run time type checking.
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:573
CSG_Grid::Add
virtual CSG_Grid & Add(const CSG_Grid &Grid)
Definition: grid_operation.cpp:503
CSG_Grid::Get_YMin
double Get_YMin(bool bCells=false) const
Definition: grid.h:550
CSG_Grid::Get_Cellsize
double Get_Cellsize(void) const
Definition: grid.h:541
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:1006
CSG_Grid::Subtract
virtual CSG_Grid & Subtract(const CSG_Grid &Grid)
Definition: grid_operation.cpp:538
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:436
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:1016
CSG_Grid::Set_Value
virtual void Set_Value(sLong i, double Value, bool bScaled=true)
Definition: grid.h:817
CSG_Unique_Number_Statistics::Get_Majority
bool Get_Majority(double &Value) const
Definition: mat_tools.h:858
CSG_Grid::Get_Gradient
bool Get_Gradient(int x, int y, double &Slope, double &Aspect) const
Definition: grid_operation.cpp:958
CSG_Grid::Get_Unit
const CSG_String & Get_Unit(void) const
Definition: grid.h:527
CSG_Grid::Get_NY
int Get_NY(void) const
Definition: grid.h:538
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::Assign
virtual bool Assign(double Value=0.0)
Definition: grid_operation.cpp:77
CSG_Grid::Get_XMin
double Get_XMin(bool bCells=false) const
Definition: grid.h:546
CSG_Grid::Standardise
bool Standardise(void)
Definition: grid_operation.cpp:848
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:534
CSG_Grid::operator-=
virtual CSG_Grid & operator-=(const CSG_Grid &Grid)
Definition: grid_operation.cpp:528
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:514
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:351
CSG_Grid::Assign_NoData
void Assign_NoData(void)
Definition: grid_operation.cpp:65
CSG_Grid::DeNormalise
bool DeNormalise(double Minimum, double Maximum)
Definition: grid_operation.cpp:825
CSG_Grid::operator+
virtual CSG_Grid operator+(const CSG_Grid &Grid) const
Definition: grid_operation.cpp:479
CSG_Grid_System::Get_yFrom
static int Get_yFrom(int Direction, int y=0)
Definition: grid.h:341
CSG_Grid::operator=
virtual CSG_Grid & operator=(const CSG_Grid &Grid)
Definition: grid_operation.cpp:464
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:704
CSG_Grid::Get_Gradient_NeighborDir
int Get_Gradient_NeighborDir(int x, int y, bool bDown=true, bool bNoEdges=true) const
Definition: grid_operation.cpp:916
GRID_OPERATION_Division
@ GRID_OPERATION_Division
Definition: grid.h:178
CSG_Grid::Get_Range
double Get_Range(void)
Definition: grid.cpp:1021
GRID_RESAMPLING_Majority
@ GRID_RESAMPLING_Majority
Definition: grid.h:166
M_PI_090
#define M_PI_090
Definition: mat_tools.h:100
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:767
CSG_Grid::Get_Min
double Get_Min(void)
Definition: grid.cpp:1011
CSG_Grid::Get_NX
int Get_NX(void) const
Definition: grid.h:537
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:340
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:608
CSG_Grid::is_Cached
bool is_Cached(void) const
Definition: grid.h:607
CSG_Grid
Definition: grid.h:475
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:549
CSG_Grid::Mirror
void Mirror(void)
Definition: grid_operation.cpp:774
CSG_Grid::Invert
void Invert(void)
Definition: grid_operation.cpp:734
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:395
CSG_Grid::operator*=
virtual CSG_Grid & operator*=(const CSG_Grid &Grid)
Definition: grid_operation.cpp:563
CSG_Grid::is_NoData
virtual bool is_NoData(int x, int y) const
Definition: grid.h:701
CSG_Grid::Get_Extent
virtual const CSG_Rect & Get_Extent(void)
Definition: grid.h:535
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:1026
CSG_Data_Object::is_Valid
virtual bool is_Valid(void) const =0
TSG_Grid_Operation
TSG_Grid_Operation
Definition: grid.h:174