📄 mrvbf.cpp
字号:
///////////////////////////////////////////////////////////
// //
// SAGA //
// //
// System for Automated Geoscientific Analyses //
// //
// Module Library: //
// ta_morphometry //
// //
//-------------------------------------------------------//
// //
// MRVBF.cpp //
// //
// Copyright (C) 2006 by //
// Olaf Conrad //
// //
//-------------------------------------------------------//
// //
// This file is part of 'SAGA - System for Automated //
// Geoscientific Analyses'. SAGA is free software; you //
// can redistribute it and/or modify it under the terms //
// of the GNU General Public License as published by the //
// Free Software Foundation; version 2 of the License. //
// //
// SAGA is distributed in the hope that it will be //
// useful, but WITHOUT ANY WARRANTY; without even the //
// implied warranty of MERCHANTABILITY or FITNESS FOR A //
// PARTICULAR PURPOSE. See the GNU General Public //
// License for more details. //
// //
// You should have received a copy of the GNU General //
// Public License along with this program; if not, //
// write to the Free Software Foundation, Inc., //
// 59 Temple Place - Suite 330, Boston, MA 02111-1307, //
// USA. //
// //
//-------------------------------------------------------//
// //
// e-mail: oconrad@saga-gis.org //
// //
// contact: Olaf Conrad //
// Institute of Geography //
// University of Goettingen //
// Goldschmidtstr. 5 //
// 37077 Goettingen //
// Germany //
// //
///////////////////////////////////////////////////////////
///////////////////////////////////////////////////////////
// //
// //
// //
///////////////////////////////////////////////////////////
//---------------------------------------------------------
#include "mrvbf.h"
///////////////////////////////////////////////////////////
// //
// //
// //
///////////////////////////////////////////////////////////
//---------------------------------------------------------
CMRVBF::CMRVBF(void)
{
Set_Name (_TL("Multiresolution Index of Valley Bottom Flatness (MRVBF)"));
Set_Author (_TL("Copyrights (c) 2006 by Olaf Conrad"));
Set_Description (_TW(
"Calculation of the 'multiresolution index of valley bottom flatness' (MRVBF) and "
"the complementary 'multiresolution index of the ridge top flatness' (MRRTF). "
"\n\n"
"References:\n"
"- Gallant, J.C., Dowling, T.I. (2003): "
" 'A multiresolution index of valley bottom flatness for mapping depositional areas', "
" Water Resources Research, 39/12:1347-1359\n"
));
//-----------------------------------------------------
Parameters.Add_Grid(
NULL , "DEM" , _TL("Elevation"),
_TL(""),
PARAMETER_INPUT
);
Parameters.Add_Grid(
NULL , "MRVBF" , _TL("MRVBF"),
_TL(""),
PARAMETER_OUTPUT
);
Parameters.Add_Grid(
NULL , "MRRTF" , _TL("MRRTF"),
_TL(""),
PARAMETER_OUTPUT
);
Parameters.Add_Value(
NULL , "T_SLOPE" , _TL("Initial Threshold for Slope"),
_TL(""),
PARAMETER_TYPE_Double , 16.0, 0.0, true, 100.0, true
);
Parameters.Add_Value(
NULL , "T_PCTL_V" , _TL("Threshold for Elevation Percentile (Lowness)"),
_TL(""),
PARAMETER_TYPE_Double , 0.40, 0.0, true, 1.0, true
);
Parameters.Add_Value(
NULL , "T_PCTL_R" , _TL("Threshold for Elevation Percentile (Upness)"),
_TL(""),
PARAMETER_TYPE_Double , 0.35, 0.0, true, 1.0, true
);
Parameters.Add_Value(
NULL , "P_SLOPE" , _TL("Shape Parameter for Slope"),
_TL(""),
PARAMETER_TYPE_Double , 4.0
);
Parameters.Add_Value(
NULL , "P_PCTL" , _TL("Shape Parameter for Elevation Percentile"),
_TL(""),
PARAMETER_TYPE_Double , 3.0
);
Parameters.Add_Value(
NULL , "UPDATE" , _TL("Update Views"),
_TL(""),
PARAMETER_TYPE_Bool , true
);
Parameters.Add_Value(
NULL , "CLASSIFY" , _TL("Classify"),
_TL(""),
PARAMETER_TYPE_Bool , false
);
Parameters.Add_Value(
NULL , "MAX_RES" , _TL("Maximum Resolution (Percentage)"),
_TL("Maximum resolution as percentage of the diameter of the DEM."),
PARAMETER_TYPE_Double , 100.0, 0.0, true, 100.0, true
);
}
//---------------------------------------------------------
CMRVBF::~CMRVBF(void)
{}
///////////////////////////////////////////////////////////
// //
// //
// //
///////////////////////////////////////////////////////////
//---------------------------------------------------------
#define UPDATE_VIEWS(b) if( bUpdate ) { DataObject_Update(pMRRTF, b); DataObject_Update(pMRVBF, b); }
//---------------------------------------------------------
bool CMRVBF::On_Execute(void)
{
bool bUpdate;
int Level;
double T_Slope, Resolution, max_Resolution;
CSG_Grid *pDEM, *pMRVBF, *pMRRTF, CF, VF, RF, DEM, Slope, Pctl;
//-----------------------------------------------------
pDEM = Parameters("DEM") ->asGrid();
pMRVBF = Parameters("MRVBF") ->asGrid();
pMRRTF = Parameters("MRRTF") ->asGrid();
T_Slope = Parameters("T_SLOPE") ->asDouble();
m_T_Pctl_V = Parameters("T_PCTL_V") ->asDouble();
m_T_Pctl_R = Parameters("T_PCTL_R") ->asDouble();
m_P_Slope = Parameters("P_SLOPE") ->asDouble();
m_P_Pctl = Parameters("P_PCTL") ->asDouble();
bUpdate = Parameters("UPDATE") ->asBool();
max_Resolution = Parameters("MAX_RES") ->asDouble() / 100.0;
Resolution = SG_Get_Length(Get_System()->Get_XRange(), Get_System()->Get_YRange());
max_Resolution = max_Resolution * Resolution;
//-----------------------------------------------------
if( 1 )
{
// DataObject_Set_Colors(pMRVBF, 100, SG_COLORS_WHITE_BLUE , false);
DataObject_Set_Colors(pMRVBF, 100, SG_COLORS_RED_GREY_BLUE , false);
// DataObject_Set_Colors(pMRRTF, 100, SG_COLORS_RED_GREY_BLUE , true);
DataObject_Set_Colors(pMRRTF, 100, SG_COLORS_WHITE_RED , false);
CSG_Grid CF, VF, RF, DEM, Slopes, Percentiles;
VF.Create(*Get_System(), GRID_TYPE_Float);
RF.Create(*Get_System(), GRID_TYPE_Float);
CF.Create(*Get_System(), GRID_TYPE_Float);
CF.Assign(1.0);
DEM.Create(*pDEM);
//-------------------------------------------------
Level = 1;
Resolution = Get_Cellsize();
Process_Set_Text(CSG_String::Format(SG_T("%d. %s"), Level, _TL("step")));
Message_Add(CSG_String::Format(SG_T("%s: %d, %s: %.2f, %s %.2f"), _TL("step"), Level, _TL("resolution"), Resolution, _TL("threshold slope"), T_Slope));
Get_Slopes (&DEM, &Slopes);
Get_Percentiles (&DEM, &Percentiles, 3);
Get_Flatness (&Slopes, &Percentiles, &CF, pMRVBF, pMRRTF, T_Slope);
UPDATE_VIEWS (true);
//-------------------------------------------------
T_Slope /= 2.0;
Level++;
Process_Set_Text(CSG_String::Format(SG_T("%d. %s"), Level, _TL("step")));
Message_Add(CSG_String::Format(SG_T("%s: %d, %s: %.2f, %s %.2f"), _TL("step"), Level, _TL("resolution"), Resolution, _TL("threshold slope"), T_Slope));
Get_Percentiles (&DEM, &Percentiles, 6);
Get_Flatness (&Slopes, &Percentiles, &CF, &VF, &RF, T_Slope);
Get_MRVBF (Level, pMRVBF, &VF, pMRRTF, &RF);
UPDATE_VIEWS (false);
//-------------------------------------------------
while( Process_Get_Okay(false) && Resolution < max_Resolution )
{
Resolution *= 3.0;
T_Slope /= 2.0;
Level++;
Process_Set_Text(CSG_String::Format(SG_T("%d. %s"), Level, _TL("step")));
Message_Add(CSG_String::Format(SG_T("%s: %d, %s: %.2f, %s %.2f"), _TL("step"), Level, _TL("resolution"), Resolution, _TL("threshold slope"), T_Slope));
Get_Values (&DEM, &Slopes, &Percentiles, Resolution);
Get_Flatness (&Slopes, &Percentiles, &CF, &VF, &RF, T_Slope);
Get_MRVBF (Level, pMRVBF, &VF, pMRRTF, &RF);
UPDATE_VIEWS (false);
}
if( Parameters("CLASSIFY")->asBool() )
{
Get_Classified(pMRVBF);
Get_Classified(pMRRTF);
}
return( true );
}
return( false );
}
///////////////////////////////////////////////////////////
// //
// //
// //
///////////////////////////////////////////////////////////
//---------------------------------------------------------
inline double CMRVBF::Get_Transformation(double x, double t, double p)
{
return( 1.0 / (1.0 + pow(x / t, p)) );
}
///////////////////////////////////////////////////////////
// //
// //
// //
///////////////////////////////////////////////////////////
//---------------------------------------------------------
bool CMRVBF::Get_Smoothed(CSG_Grid *pDEM, CSG_Grid *pSmoothed, int Radius, double Smoothing)
{
if( pDEM && pSmoothed )
{
int x, y, ix, iy, jx, jy, kx, ky;
double d, s;
CSG_Grid Kernel(GRID_TYPE_Double, 1 + 2 * Radius, 1 + 2 * Radius);
for(iy=-Radius, ky=0; iy<=Radius; iy++, ky++)
{
for(ix=-Radius, kx=0; ix<=Radius; ix++, kx++)
{
Kernel.Set_Value(kx, ky, 4.3565 * exp(-SG_Get_Square(sqrt((double)ix*ix + iy*iy) / 3.0)) );
// Kernel.Set_Value(kx, ky, exp(-(ix*ix + iy*iy) / (2.0 * s)) / (2.0 * M_PI * s));
}
}
pSmoothed->Create(pDEM, GRID_TYPE_Float);
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -