1
2 ///////////////////////////////////////////////////////////
3 // //
4 // SAGA //
5 // //
6 // System for Automated Geoscientific Analyses //
7 // //
8 // Tool Library //
9 // climate_tools //
10 // //
11 //-------------------------------------------------------//
12 // //
13 // Grid_Level_Interpolation.cpp //
14 // //
15 // Copyright (C) 2012 by //
16 // Olaf Conrad //
17 // //
18 //-------------------------------------------------------//
19 // //
20 // This file is part of 'SAGA - System for Automated //
21 // Geoscientific Analyses'. SAGA is free software; you //
22 // can redistribute it and/or modify it under the terms //
23 // of the GNU General Public License as published by the //
24 // Free Software Foundation, either version 2 of the //
25 // License, or (at your option) any later version. //
26 // //
27 // SAGA is distributed in the hope that it will be //
28 // useful, but WITHOUT ANY WARRANTY; without even the //
29 // implied warranty of MERCHANTABILITY or FITNESS FOR A //
30 // PARTICULAR PURPOSE. See the GNU General Public //
31 // License for more details. //
32 // //
33 // You should have received a copy of the GNU General //
34 // Public License along with this program; if not, see //
35 // <http://www.gnu.org/licenses/>. //
36 // //
37 //-------------------------------------------------------//
38 // //
39 // e-mail: oconrad@saga-gis.org //
40 // //
41 // contact: Olaf Conrad //
42 // Institute of Geography //
43 // University of Hamburg //
44 // Germany //
45 // //
46 ///////////////////////////////////////////////////////////
47
48 //---------------------------------------------------------
49 #include "grid_levels_interpolation.h"
50
51
52 ///////////////////////////////////////////////////////////
53 // //
54 // //
55 // //
56 ///////////////////////////////////////////////////////////
57
58 //---------------------------------------------------------
CGrid_Levels_Interpolation(void)59 CGrid_Levels_Interpolation::CGrid_Levels_Interpolation(void)
60 {
61 //-----------------------------------------------------
62 Set_Author("O.Conrad (c) 2012");
63
64 //-----------------------------------------------------
65 Parameters.Add_Grid_List("",
66 "VARIABLE" , _TL("Variable"),
67 _TL(""),
68 PARAMETER_INPUT
69 );
70
71 Parameters.Add_Choice("",
72 "X_SOURCE" , _TL("Get Heights from ..."),
73 _TL(""),
74 CSG_String::Format("%s|%s",
75 _TL("table"),
76 _TL("grid list")
77 ), 1
78 );
79
80 Parameters.Add_Grid_List("",
81 "X_GRIDS" , _TL("Level Heights"),
82 _TL(""),
83 PARAMETER_INPUT_OPTIONAL
84 );
85
86 Parameters.Add_Grid("",
87 "X_GRIDS_CHECK" , _TL("Minimum Height"),
88 _TL("if set, only values with level heights above DEM will be used"),
89 PARAMETER_INPUT_OPTIONAL, true
90 );
91
92 Parameters.Add_FixedTable("",
93 "X_TABLE" , _TL("Level Heights"),
94 _TL("")
95 )->asTable()->Add_Field(_TL("Height"), SG_DATATYPE_Double);
96
97 Parameters.Add_Choice("",
98 "H_METHOD" , _TL("Horizontal Interpolation Method"),
99 _TL(""),
100 CSG_String::Format("%s|%s|%s|%s",
101 _TL("Nearest Neighbour"),
102 _TL("Bilinear Interpolation"),
103 _TL("Bicubic Spline Interpolation"),
104 _TL("B-Spline Interpolation")
105 ), 3
106 );
107
108 Parameters.Add_Choice("",
109 "V_METHOD" , _TL("Vertical Interpolation Method"),
110 _TL(""),
111 CSG_String::Format("%s|%s|%s",
112 _TL("linear"),
113 _TL("spline"),
114 _TL("polynomial trend")
115 ), 0
116 );
117
118 Parameters.Add_Bool("V_METHOD",
119 "COEFFICIENTS" , _TL("Coefficient Interpolation"),
120 _TL(""),
121 false
122 );
123
124 Parameters.Add_Bool("V_METHOD",
125 "LINEAR_SORTED" , _TL("Sorted Levels"),
126 _TL(""),
127 false
128 );
129
130 Parameters.Add_Bool("V_METHOD",
131 "SPLINE_ALL" , _TL("Pre-analyze"),
132 _TL(""),
133 false
134 );
135
136 Parameters.Add_Int("V_METHOD",
137 "TREND_ORDER" , _TL("Polynomial Order"),
138 _TL(""),
139 3, 1, true
140 );
141
142 //-----------------------------------------------------
143 for(int i=0; i<10; i++)
144 {
145 Parameters("X_TABLE")->asTable()->Add_Record()->Set_Value(0, i + 1);
146 }
147
148 Add_Parameters("INTERNAL", "", "");
149
150 m_Coeff = NULL;
151 }
152
153
154 ///////////////////////////////////////////////////////////
155 // //
156 ///////////////////////////////////////////////////////////
157
158 //---------------------------------------------------------
On_Parameters_Enable(CSG_Parameters * pParameters,CSG_Parameter * pParameter)159 int CGrid_Levels_Interpolation::On_Parameters_Enable(CSG_Parameters *pParameters, CSG_Parameter *pParameter)
160 {
161 if( pParameter->Cmp_Identifier("V_METHOD") )
162 {
163 pParameters->Set_Enabled("COEFFICIENTS" , pParameter->asInt() != 1 && Parameters("SURFACE")); // not available for splines, needs reference surface
164 pParameters->Set_Enabled("LINEAR_SORTED", pParameter->asInt() == 0);
165 pParameters->Set_Enabled("SPLINE_ALL" , pParameter->asInt() == 1);
166 pParameters->Set_Enabled("TREND_ORDER" , pParameter->asInt() >= 2);
167 }
168
169 if( pParameter->Cmp_Identifier("X_SOURCE") )
170 {
171 pParameters->Set_Enabled("X_TABLE" , pParameter->asInt() == 0);
172 pParameters->Set_Enabled("X_GRIDS" , pParameter->asInt() == 1);
173 pParameters->Set_Enabled("X_GRIDS_CHECK", pParameter->asInt() == 1);
174 }
175
176 return( CSG_Tool_Grid::On_Parameters_Enable(pParameters, pParameter) );
177 }
178
179
180 ///////////////////////////////////////////////////////////
181 // //
182 ///////////////////////////////////////////////////////////
183
184 //---------------------------------------------------------
Initialize(const CSG_Rect & Extent)185 bool CGrid_Levels_Interpolation::Initialize(const CSG_Rect &Extent)
186 {
187 //-----------------------------------------------------
188 m_pVariables = Parameters("VARIABLE" )->asGridList();
189 m_pXGrids = Parameters("X_GRIDS" )->asGridList();
190 m_pXTable = Parameters("X_TABLE" )->asTable();
191
192 m_xSource = Parameters("X_SOURCE" )->asInt();
193 m_vMethod = Parameters("V_METHOD" )->asInt();
194
195 m_Linear_bSorted = Parameters("LINEAR_SORTED")->asBool();
196 m_Spline_bAll = Parameters("SPLINE_ALL" )->asBool() == false;
197 m_Trend_Order = Parameters("TREND_ORDER" )->asInt();
198
199 switch( Parameters("H_METHOD")->asInt() )
200 {
201 default: m_hMethod = GRID_RESAMPLING_NearestNeighbour; break;
202 case 1: m_hMethod = GRID_RESAMPLING_Bilinear ; break;
203 case 2: m_hMethod = GRID_RESAMPLING_BicubicSpline ; break;
204 case 3: m_hMethod = GRID_RESAMPLING_BSpline ; break;
205 }
206
207 //-----------------------------------------------------
208 if( m_pVariables->Get_Grid_Count() != (m_xSource == 0 ? m_pXTable->Get_Count() : m_pXGrids->Get_Grid_Count()) )
209 {
210 Error_Set(_TL("variable and height levels have to be of same number"));
211
212 return( false );
213 }
214
215 if( m_vMethod == 2 && m_pVariables->Get_Grid_Count() <= m_Trend_Order )
216 {
217 Error_Set(_TL("fitting a polynom of ith order needs at least i + 1 samples"));
218
219 return( false );
220 }
221
222 if( !Extent.Intersects(Get_System().Get_Extent(true)) )
223 {
224 Error_Set(_TL("target area is distinct from levels area "));
225
226 return( false );
227 }
228
229 //-----------------------------------------------------
230 CSG_Grid *pHeight_Min = m_xSource == 1 && Parameters("X_GRIDS_CHECK") ? Parameters("X_GRIDS_CHECK")->asGrid() : NULL;
231
232 if( pHeight_Min )
233 {
234 if( !Get_Parameters("INTERNAL")->Get_Parameter("X_GRIDS") )
235 {
236 Get_Parameters("INTERNAL")->Add_Grid_List("", "X_GRIDS", "", "", PARAMETER_INPUT_OPTIONAL);
237 }
238
239 CSG_Parameter_Grid_List *pXGrids = Get_Parameters("INTERNAL")->Get_Parameter("X_GRIDS")->asGridList();
240
241 for(int i=0; i<m_pXGrids->Get_Grid_Count(); i++)
242 {
243 CSG_Grid *pHeight = SG_Create_Grid(*m_pXGrids->Get_Grid(i));
244
245 #pragma omp parallel for
246 for(int y=0; y<Get_NY(); y++)
247 {
248 for(int x=0; x<Get_NX(); x++)
249 {
250 if( pHeight->asDouble(x, y) < pHeight_Min->asDouble(x, y) )
251 {
252 pHeight->Set_NoData(x, y);
253 }
254 }
255 }
256
257 pXGrids->Add_Item(pHeight);
258 }
259
260 m_pXGrids = pXGrids;
261 }
262
263 //-----------------------------------------------------
264 if( m_vMethod == 0 && Parameters("COEFFICIENTS")->asBool() && Parameters("SURFACE") ) // linear coefficients interpolation
265 {
266 CSG_Grid Surface(Get_System());
267
268 Surface.Assign(Parameters("SURFACE")->asGrid(), GRID_RESAMPLING_Mean_Cells);
269
270 m_Trend_Order = 1; // linear
271
272 m_Coeff = new CSG_Grid[1 + m_Trend_Order];
273
274 for(int i=0; i<=m_Trend_Order; i++)
275 {
276 if( !m_Coeff[i].Create(Get_System()) )
277 {
278 return( false );
279 }
280 }
281
282 #pragma omp parallel
283 for(int y=0; y<Get_NY(); y++)
284 {
285 double p_y = Get_YMin() + y * Get_Cellsize();
286
287 for(int x=0; x<Get_NX(); x++)
288 {
289 double p_x = Get_XMin() + x * Get_Cellsize();
290
291 double zz[2], vv[2];
292
293 CSG_Trend_Polynom Trend; Trend.Set_Order(m_Trend_Order);
294
295 if( Get_Linear_Coeff(p_x, p_y, Surface.asDouble(x, y), vv, zz) )
296 {
297 Trend.Add_Data(zz[0], vv[0]);
298 Trend.Add_Data(zz[1], vv[1]);
299 }
300
301 if( Trend.Get_Trend() )
302 {
303 m_Coeff[0].Set_Value(x, y, Trend.Get_Coefficient(0));
304 m_Coeff[1].Set_Value(x, y, Trend.Get_Coefficient(1));
305 }
306 else
307 {
308 m_Coeff[0].Set_NoData(x, y);
309 m_Coeff[1].Set_NoData(x, y);
310 }
311 }
312 }
313 }
314
315 //-----------------------------------------------------
316 if( m_vMethod == 2 && Parameters("COEFFICIENTS")->asBool() ) // polynomial trend coefficients interpolation
317 {
318 int i;
319
320 m_Coeff = new CSG_Grid[1 + m_Trend_Order];
321
322 for(i=0; i<=m_Trend_Order; i++)
323 {
324 if( !m_Coeff[i].Create(Get_System()) )
325 {
326 return( false );
327 }
328 }
329
330 #pragma omp parallel for private(i)
331 for(int y=0; y<Get_NY(); y++)
332 {
333 double p_y = Get_YMin() + y * Get_Cellsize();
334
335 for(int x=0; x<Get_NX(); x++)
336 {
337 double p_x = Get_XMin() + x * Get_Cellsize();
338
339 CSG_Trend_Polynom Trend; Trend.Set_Order(m_Trend_Order);
340
341 for(i=0; i<m_pVariables->Get_Grid_Count(); i++)
342 {
343 double Height, Variable;
344
345 if( Get_Height(p_x, p_y, i, Height) && Get_Variable(p_x, p_y, i, Variable) )
346 {
347 Trend.Add_Data(Height, Variable);
348 }
349 }
350
351 if( Trend.Get_Trend() )
352 {
353 for(i=0; i<=m_Trend_Order; i++)
354 {
355 m_Coeff[i].Set_Value(x, y, Trend.Get_Coefficient(i));
356 }
357 }
358 else
359 {
360 for(i=0; i<=m_Trend_Order; i++)
361 {
362 m_Coeff[i].Set_NoData(x, y);
363 }
364 }
365 }
366 }
367 }
368
369 //-----------------------------------------------------
370 return( true );
371 }
372
373 //---------------------------------------------------------
Finalize(void)374 bool CGrid_Levels_Interpolation::Finalize(void)
375 {
376 if( Get_Parameters("INTERNAL")->Get_Parameter("X_GRIDS")
377 && Get_Parameters("INTERNAL")->Get_Parameter("X_GRIDS")->asGridList() == m_pXGrids )
378 {
379 for(int i=0; i<m_pXGrids->Get_Grid_Count(); i++)
380 {
381 delete(m_pXGrids->Get_Grid(i));
382 }
383
384 m_pXGrids->Del_Items();
385 }
386
387 if( m_Coeff )
388 {
389 delete[](m_Coeff);
390
391 m_Coeff = NULL;
392 }
393
394 return( true );
395 }
396
397
398 ///////////////////////////////////////////////////////////
399 // //
400 ///////////////////////////////////////////////////////////
401
402 //---------------------------------------------------------
Get_Variable(double x,double y,int iLevel)403 inline double CGrid_Levels_Interpolation::Get_Variable(double x, double y, int iLevel)
404 {
405 return( m_pVariables->Get_Grid(iLevel)->Get_Value(x, y, m_hMethod) );
406 }
407
408 //---------------------------------------------------------
Get_Variable(double x,double y,int iLevel,double & Variable)409 inline bool CGrid_Levels_Interpolation::Get_Variable(double x, double y, int iLevel, double &Variable)
410 {
411 return( m_pVariables->Get_Grid(iLevel)->Get_Value(x, y, Variable, m_hMethod) );
412 }
413
414 //---------------------------------------------------------
Get_Height(double x,double y,int iLevel)415 inline double CGrid_Levels_Interpolation::Get_Height(double x, double y, int iLevel)
416 {
417 if( m_xSource == 0 )
418 {
419 return( m_pXTable->Get_Record(iLevel)->asDouble(0) );
420 }
421
422 return( m_pXGrids->Get_Grid(iLevel)->Get_Value(x, y, m_hMethod) );
423 }
424
425 //---------------------------------------------------------
Get_Height(double x,double y,int iLevel,double & Height)426 inline bool CGrid_Levels_Interpolation::Get_Height(double x, double y, int iLevel, double &Height)
427 {
428 if( m_xSource == 0 )
429 {
430 Height = m_pXTable->Get_Record(iLevel)->asDouble(0);
431
432 return( true );
433 }
434
435 return( m_pXGrids->Get_Grid(iLevel)->Get_Value(x, y, Height, m_hMethod) );
436 }
437
438
439 ///////////////////////////////////////////////////////////
440 // //
441 ///////////////////////////////////////////////////////////
442
443 //---------------------------------------------------------
Get_Values(double x,double y,double z,int & iLevel,CSG_Table & Values)444 bool CGrid_Levels_Interpolation::Get_Values(double x, double y, double z, int &iLevel, CSG_Table &Values)
445 {
446 Values.Destroy();
447
448 Values.Add_Field("h", SG_DATATYPE_Double);
449 Values.Add_Field("v", SG_DATATYPE_Double);
450
451 for(int i=0; i<m_pVariables->Get_Grid_Count(); i++)
452 {
453 double Height, Variable;
454
455 if( Get_Height(x, y, i, Height) && Get_Variable(x, y, i, Variable) )
456 {
457 CSG_Table_Record *pRecord = Values.Add_Record();
458
459 pRecord->Set_Value(0, Height);
460 pRecord->Set_Value(1, Variable);
461 }
462 }
463
464 if( Values.Get_Count() < 2 || !Values.Set_Index(0, TABLE_INDEX_Ascending) )
465 {
466 return( false );
467 }
468
469 for(iLevel=1; iLevel<Values.Get_Count()-1; iLevel++)
470 {
471 if( Values[iLevel][0].asDouble() > z )
472 {
473 return( true );
474 }
475 }
476
477 return( true );
478 }
479
480
481 ///////////////////////////////////////////////////////////
482 // //
483 ///////////////////////////////////////////////////////////
484
485 //---------------------------------------------------------
Get_Value(double x,double y,double z,double & Value)486 inline bool CGrid_Levels_Interpolation::Get_Value(double x, double y, double z, double &Value)
487 {
488 switch( m_vMethod )
489 {
490 default: // linear
491 if( m_Coeff )
492 {
493 return( Get_Trend_Coeff (x, y, z, Value) );
494 }
495 else
496 {
497 return( Get_Linear (x, y, z, Value) );
498 }
499
500 case 1: // spline
501 if( m_Spline_bAll )
502 {
503 return( Get_Spline_All (x, y, z, Value) );
504 }
505 else
506 {
507 return( Get_Spline (x, y, z, Value) );
508 }
509
510 case 2: // polynomial trend
511 if( m_Coeff )
512 {
513 return( Get_Trend_Coeff (x, y, z, Value) );
514 }
515 else
516 {
517 return( Get_Trend (x, y, z, Value) );
518 }
519 }
520 }
521
522 //---------------------------------------------------------
Get_Value(const TSG_Point & p,double z,double & Value)523 inline bool CGrid_Levels_Interpolation::Get_Value(const TSG_Point &p, double z, double &Value)
524 {
525 return( Get_Value(p.x, p.y, z, Value) );
526 }
527
528
529 ///////////////////////////////////////////////////////////
530 // //
531 ///////////////////////////////////////////////////////////
532
533 //---------------------------------------------------------
Get_Linear_Coeff(double x,double y,double z,double vv[2],double hh[2])534 bool CGrid_Levels_Interpolation::Get_Linear_Coeff(double x, double y, double z, double vv[2], double hh[2])
535 {
536 int iLevel;
537
538 if( m_Linear_bSorted )
539 {
540 for(iLevel=1; iLevel<m_pVariables->Get_Grid_Count()-1; iLevel++)
541 {
542 if( Get_Height(x, y, iLevel) > z )
543 {
544 break;
545 }
546 }
547
548 hh[0] = Get_Height (x, y, iLevel - 1);
549 vv[0] = Get_Variable(x, y, iLevel - 1);
550 hh[1] = Get_Height (x, y, iLevel );
551 vv[1] = Get_Variable(x, y, iLevel );
552 }
553 else
554 {
555 CSG_Table Values;
556
557 if( !Get_Values(x, y, z, iLevel, Values) )
558 {
559 return( false );
560 }
561
562 hh[0] = Values[iLevel - 1][0];
563 vv[0] = Values[iLevel - 1][1];
564 hh[1] = Values[iLevel ][0];
565 vv[1] = Values[iLevel ][1];
566 }
567
568 return( hh[0] < hh[1] );
569 }
570
571 //---------------------------------------------------------
Get_Linear(double x,double y,double z,double & Value)572 bool CGrid_Levels_Interpolation::Get_Linear(double x, double y, double z, double &Value)
573 {
574 double zz[2], vv[2];
575
576 if( Get_Linear_Coeff(x, y, z, vv, zz) )
577 {
578 Value = vv[0] + (z - zz[0]) * (vv[1] - vv[0]) / (zz[1] - zz[0]);
579
580 return( true );
581 }
582
583 return( false );
584 }
585
586
587 ///////////////////////////////////////////////////////////
588 // //
589 ///////////////////////////////////////////////////////////
590
591 //---------------------------------------------------------
Get_Spline_All(double x,double y,double z,double & Value)592 bool CGrid_Levels_Interpolation::Get_Spline_All(double x, double y, double z, double &Value)
593 {
594 CSG_Spline Spline;
595
596 for(int i=0; i<m_pVariables->Get_Grid_Count(); i++)
597 {
598 double Height, Variable;
599
600 if( Get_Height(x, y, i, Height) && Get_Variable(x, y, i, Variable) )
601 {
602 Spline.Add(Height, Variable);
603 }
604 }
605
606 if( Spline.Get_Value(z, Value) )
607 {
608 return( true );
609 }
610
611 return( false );
612 }
613
614 //---------------------------------------------------------
Get_Spline(double x,double y,double z,double & Value)615 bool CGrid_Levels_Interpolation::Get_Spline(double x, double y, double z, double &Value)
616 {
617 int iLevel;
618 CSG_Table Values;
619
620 if( !Get_Values(x, y, z, iLevel, Values) )
621 {
622 return( false );
623 }
624
625 if( Values.Get_Count() < 3 )
626 {
627 return( Get_Linear(x, y, z, Value) );
628 }
629
630 if( iLevel >= Values.Get_Count() - 1 )
631 {
632 iLevel--;
633 }
634
635 CSG_Spline Spline;
636
637 if( iLevel > 1 )
638 {
639 Spline.Add(Values[iLevel - 2][0], Values[iLevel - 2][1]);
640 }
641
642 Spline.Add(Values[iLevel - 1][0], Values[iLevel - 1][1]);
643 Spline.Add(Values[iLevel ][0], Values[iLevel ][1]);
644 Spline.Add(Values[iLevel + 1][0], Values[iLevel + 1][1]);
645
646 return( Spline.Get_Value(z, Value) );
647 }
648
649
650 ///////////////////////////////////////////////////////////
651 // //
652 ///////////////////////////////////////////////////////////
653
654 //---------------------------------------------------------
Get_Trend(double x,double y,double z,double & Value)655 bool CGrid_Levels_Interpolation::Get_Trend(double x, double y, double z, double &Value)
656 {
657 CSG_Trend_Polynom Trend;
658
659 Trend.Set_Order(m_Trend_Order);
660
661 for(int i=0; i<m_pVariables->Get_Grid_Count(); i++)
662 {
663 double Height, Variable;
664
665 if( Get_Height(x, y, i, Height) && Get_Variable(x, y, i, Variable) )
666 {
667 Trend.Add_Data(Height, Variable);
668 }
669 }
670
671 if( Trend.Get_Trend() )
672 {
673 Value = Trend.Get_Value(z);
674
675 return( true );
676 }
677
678 return( false );
679 }
680
681
682 ///////////////////////////////////////////////////////////
683 // //
684 ///////////////////////////////////////////////////////////
685
686 //---------------------------------------------------------
Get_Trend_Coeff(double x,double y,double z,double & Value)687 bool CGrid_Levels_Interpolation::Get_Trend_Coeff(double x, double y, double z, double &Value)
688 {
689 double Coeff, zPower = 1.0;
690
691 Value = 0.0;
692
693 for(int i=0; i<=m_Trend_Order; i++)
694 {
695 if( !m_Coeff[i].Get_Value(x, y, Coeff, m_hMethod) )
696 {
697 return( false );
698 }
699
700 Value += Coeff * zPower;
701 zPower *= z;
702 }
703
704 return( true );
705 }
706
707
708 ///////////////////////////////////////////////////////////
709 // //
710 // //
711 // //
712 ///////////////////////////////////////////////////////////
713
714 //---------------------------------------------------------
CGrid_Levels_to_Surface(void)715 CGrid_Levels_to_Surface::CGrid_Levels_to_Surface(void)
716 {
717 //-----------------------------------------------------
718 Set_Name (_TL("Multi Level to Surface Interpolation"));
719
720 Set_Description (_TW(
721 ""
722 ));
723
724 //-----------------------------------------------------
725 Parameters.Add_Grid_System("",
726 "SYSTEM" , _TL("Grid system"),
727 _TL("")
728 );
729
730 Parameters.Add_Grid("SYSTEM",
731 "SURFACE" , _TL("Surface"),
732 _TL(""),
733 PARAMETER_INPUT
734 );
735
736 Parameters.Add_Grid("SYSTEM",
737 "RESULT" , _TL("Interpolation"),
738 _TL(""),
739 PARAMETER_OUTPUT
740 );
741 }
742
743 //---------------------------------------------------------
On_Execute(void)744 bool CGrid_Levels_to_Surface::On_Execute(void)
745 {
746 //-----------------------------------------------------
747 CSG_Grid *pSurface = Parameters("SURFACE")->asGrid();
748 CSG_Grid *pResult = Parameters("RESULT" )->asGrid();
749
750 if( !(pSurface->Get_System() == pResult->Get_System()) )
751 {
752 Error_Set(_TL("surface and result grids have to share the same grid system"));
753
754 return( false );
755 }
756
757 if( !Initialize(pSurface->Get_Extent()) )
758 {
759 Finalize();
760
761 return( false );
762 }
763
764 //-----------------------------------------------------
765 for(int y=0; y<pSurface->Get_NY() && Set_Progress(y, pSurface->Get_NY()); y++)
766 {
767 double p_y = pSurface->Get_YMin() + y * pSurface->Get_Cellsize();
768
769 #pragma omp parallel for
770 for(int x=0; x<pSurface->Get_NX(); x++)
771 {
772 double Value, p_x = pSurface->Get_XMin() + x * pSurface->Get_Cellsize();
773
774 if( !pSurface->is_NoData(x, y) && Get_Value(p_x, p_y, pSurface->asDouble(x, y), Value) )
775 {
776 pResult->Set_Value(x, y, Value);
777 }
778 else
779 {
780 pResult->Set_NoData(x, y);
781 }
782 }
783 }
784
785 //-----------------------------------------------------
786 Finalize();
787
788 return( true );
789 }
790
791
792 ///////////////////////////////////////////////////////////
793 // //
794 // //
795 // //
796 ///////////////////////////////////////////////////////////
797
798 //---------------------------------------------------------
CGrid_Levels_to_Points(void)799 CGrid_Levels_to_Points::CGrid_Levels_to_Points(void)
800 {
801 //-----------------------------------------------------
802 Set_Name (_TL("Multi Level to Points Interpolation"));
803
804 Set_Description (_TW(
805 ""
806 ));
807
808 //-----------------------------------------------------
809 Parameters.Add_Shapes("",
810 "POINTS" , _TL("Points"),
811 _TL(""),
812 PARAMETER_INPUT, SHAPE_TYPE_Point
813 );
814
815 Parameters.Add_Table_Field("POINTS",
816 "ZFIELD" , _TL("Height"),
817 _TL("")
818 );
819
820 Parameters.Add_Shapes("POINTS",
821 "RESULT" , _TL("Result"),
822 _TL(""),
823 PARAMETER_OUTPUT_OPTIONAL, SHAPE_TYPE_Point
824 );
825
826 Parameters.Add_String("",
827 "NAME" , _TL("Field Name"),
828 _TL(""),
829 _TL("Variable")
830 );
831 }
832
833 //---------------------------------------------------------
On_Execute(void)834 bool CGrid_Levels_to_Points::On_Execute(void)
835 {
836 //-----------------------------------------------------
837 CSG_Shapes *pPoints = Parameters("POINTS")->asShapes();
838
839 if( !Initialize(pPoints->Get_Extent()) )
840 {
841 Finalize();
842
843 return( false );
844 }
845
846 //-----------------------------------------------------
847 if( Parameters("RESULT")->asShapes() && Parameters("RESULT")->asShapes() != pPoints )
848 {
849 Parameters("RESULT")->asShapes()->Create(*pPoints);
850
851 pPoints = Parameters("RESULT")->asShapes();
852 }
853
854 //-----------------------------------------------------
855 int zField = Parameters("ZFIELD")->asInt();
856 int vField = pPoints->Get_Field_Count();
857
858 CSG_String Name = Parameters("NAME")->asString(); if( Name.is_Empty() ) Name = _TL("Variable");
859
860 pPoints->Add_Field(Name, SG_DATATYPE_Double);
861
862 //-----------------------------------------------------
863 // #pragma omp parallel for
864 for(int iPoint=0; iPoint<pPoints->Get_Count() && Set_Progress(iPoint, pPoints->Get_Count()); iPoint++)
865 {
866 CSG_Shape *pPoint = pPoints->Get_Shape(iPoint);
867
868 double Value;
869
870 if( !pPoint->is_NoData(zField) && Get_Value(pPoint->Get_Point(0), pPoint->asDouble(zField), Value) )
871 {
872 pPoint->Set_Value(vField, Value);
873 }
874 else
875 {
876 pPoint->Set_NoData(vField);
877 }
878 }
879
880 //-----------------------------------------------------
881 if( Parameters("RESULT")->asShapes() == NULL )
882 {
883 DataObject_Update(pPoints);
884 }
885
886 Finalize();
887
888 return( true );
889 }
890
891
892 ///////////////////////////////////////////////////////////
893 // //
894 // //
895 // //
896 ///////////////////////////////////////////////////////////
897
898 //---------------------------------------------------------
899