📄 grid.cpp
字号:
// //
///////////////////////////////////////////////////////////
//---------------------------------------------------------
void CSG_Grid::Set_NoData_Value(double Value)
{
Set_NoData_Value_Range(Value, Value);
}
void CSG_Grid::Set_NoData_Value_Range(double loValue, double hiValue)
{
double d;
if( loValue > hiValue )
{
d = loValue;
loValue = hiValue;
hiValue = d;
}
if( !m_bUpdate )
{
m_bUpdate = loValue != m_NoData_Value || hiValue != m_NoData_hiValue;
}
m_NoData_Value = loValue;
m_NoData_hiValue = hiValue;
}
///////////////////////////////////////////////////////////
// //
// Checks //
// //
///////////////////////////////////////////////////////////
//---------------------------------------------------------
bool CSG_Grid::is_Valid(void) const
{
if( m_System.is_Valid() && m_Type > 0 && m_Type < GRID_TYPE_Count )
{
switch( m_Memory_Type )
{
default:
return( m_Values != NULL );
case GRID_MEMORY_Cache:
return( Cache_Stream.is_Open() );
}
}
return( false );
}
//---------------------------------------------------------
TSG_Intersection CSG_Grid::is_Intersecting(const CSG_Rect &Extent) const
{
return( Get_Extent().Intersects(Extent.m_rect) );
}
TSG_Intersection CSG_Grid::is_Intersecting(const TSG_Rect &Extent) const
{
return( Get_Extent().Intersects(Extent) );
}
TSG_Intersection CSG_Grid::is_Intersecting(double xMin, double yMin, double xMax, double yMax) const
{
return( is_Intersecting(CSG_Rect(xMin, yMin, xMax, yMax)) );
}
//---------------------------------------------------------
bool CSG_Grid::is_Compatible(CSG_Grid *pGrid) const
{
return( pGrid && is_Compatible(pGrid->Get_System()) );
}
bool CSG_Grid::is_Compatible(const CSG_Grid_System &System) const
{
return( m_System == System );
}
bool CSG_Grid::is_Compatible(int NX, int NY, double Cellsize, double xMin, double yMin) const
{
return( is_Compatible(CSG_Grid_System(Cellsize, xMin, yMin, NX, NY)) );
}
///////////////////////////////////////////////////////////
// //
// Value access by Position (-> Interpolation) //
// //
///////////////////////////////////////////////////////////
//---------------------------------------------------------
double CSG_Grid::Get_Value(TSG_Point Position, int Interpolation, bool bZFactor, bool bByteWise) const
{
double Value;
return( Get_Value(Position.x, Position.y, Value, Interpolation, bZFactor, bByteWise) ? Value : m_NoData_Value );
}
double CSG_Grid::Get_Value(double xPosition, double yPosition, int Interpolation, bool bZFactor, bool bByteWise) const
{
double Value;
return( Get_Value(xPosition, yPosition, Value, Interpolation, bZFactor, bByteWise) ? Value : m_NoData_Value );
}
bool CSG_Grid::Get_Value(TSG_Point Position, double &Value, int Interpolation, bool bZFactor, bool bByteWise) const
{
return( Get_Value(Position.x, Position.y, Value, Interpolation, bZFactor, bByteWise) );
}
//---------------------------------------------------------
inline bool CSG_Grid::Get_Value(double xPosition, double yPosition, double &Value, int Interpolation, bool bZFactor, bool bByteWise) const
{
int x, y;
double dx, dy;
if( Get_Extent().Contains(xPosition, yPosition) )
{
x = (int)(xPosition = (xPosition - Get_XMin()) / Get_Cellsize());
y = (int)(yPosition = (yPosition - Get_YMin()) / Get_Cellsize());
dx = xPosition - x;
dy = yPosition - y;
switch( Interpolation )
{
default:
return( false );
case GRID_INTERPOLATION_NearestNeighbour:
Value = _Get_ValAtPos_NearestNeighbour(x, y, dx, dy);
break;
case GRID_INTERPOLATION_Bilinear:
Value = _Get_ValAtPos_BiLinear (x, y, dx, dy, bByteWise);
break;
case GRID_INTERPOLATION_InverseDistance:
Value = _Get_ValAtPos_InverseDistance (x, y, dx, dy, bByteWise);
break;
case GRID_INTERPOLATION_BicubicSpline:
Value = _Get_ValAtPos_BiCubicSpline (x, y, dx, dy, bByteWise);
break;
case GRID_INTERPOLATION_BSpline:
Value = _Get_ValAtPos_BSpline (x, y, dx, dy, bByteWise);
break;
}
if( Value != m_NoData_Value )
{
if( bZFactor )
{
Value *= m_zFactor;
}
return( true );
}
}
return( false );
}
//---------------------------------------------------------
inline double CSG_Grid::_Get_ValAtPos_NearestNeighbour(int x, int y, double dx, double dy) const
{
x += (int)(0.5 + dx);
y += (int)(0.5 + dy);
if( is_InGrid(x, y) )
{
return( asDouble(x, y) );
}
return( m_NoData_Value );
}
//---------------------------------------------------------
#define BILINEAR_ADD(ix, iy, d) if( is_InGrid(ix, iy) ) { n += d;\
z += d * asDouble(ix, iy); }
#define BILINEAR_ADD_BYTE(ix, iy, d) if( is_InGrid(ix, iy) ) { n += d; v = asInt(ix, iy);\
z[0] += d * SG_GET_BYTE_0(v);\
z[1] += d * SG_GET_BYTE_1(v);\
z[2] += d * SG_GET_BYTE_2(v);\
z[3] += d * SG_GET_BYTE_3(v); }
inline double CSG_Grid::_Get_ValAtPos_BiLinear(int x, int y, double dx, double dy, bool bByteWise) const
{
if( !bByteWise )
{
double z = 0.0, n = 0.0;
BILINEAR_ADD(x , y , (1.0 - dx) * (1.0 - dy));
BILINEAR_ADD(x + 1, y , ( dx) * (1.0 - dy));
BILINEAR_ADD(x , y + 1, (1.0 - dx) * ( dy));
BILINEAR_ADD(x + 1, y + 1, ( dx) * ( dy));
if( n > 0.0 )
{
return( z / n );
}
}
else
{
long v;
double z[4], n = 0.0;
z[0] = z[1] = z[2] = z[3] = 0.0;
BILINEAR_ADD_BYTE(x , y , (1.0 - dx) * (1.0 - dy));
BILINEAR_ADD_BYTE(x + 1, y , ( dx) * (1.0 - dy));
BILINEAR_ADD_BYTE(x , y + 1, (1.0 - dx) * ( dy));
BILINEAR_ADD_BYTE(x + 1, y + 1, ( dx) * ( dy));
if( n > 0.0 )
{
z[0] /= n;
z[1] /= n;
z[2] /= n;
z[3] /= n;
return( SG_GET_LONG(z[0], z[1], z[2], z[3]) );
}
}
return( m_NoData_Value );
}
//---------------------------------------------------------
#define INVERSEDIST_ADD(ix, iy, dsx, dsy) if( is_InGrid(ix, iy) ) { d = 1.0 / sqrt((dsx)*(dsx) + (dsy)*(dsy)); n += d;\
z += d * asDouble(ix, iy); }
#define INVERSEDIST_ADD_BYTE(ix, iy, dsx, dsy) if( is_InGrid(ix, iy) ) { d = 1.0 / sqrt((dsx)*(dsx) + (dsy)*(dsy)); n += d; v = asInt(ix, iy);\
z[0] += d * SG_GET_BYTE_0(v);\
z[1] += d * SG_GET_BYTE_1(v);\
z[2] += d * SG_GET_BYTE_2(v);\
z[3] += d * SG_GET_BYTE_3(v); }
inline double CSG_Grid::_Get_ValAtPos_InverseDistance(int x, int y, double dx, double dy, bool bByteWise) const
{
if( dx > 0.0 || dy > 0.0 )
{
if( !bByteWise )
{
double z = 0.0, n = 0.0, d;
INVERSEDIST_ADD(x , y , dx, dy);
INVERSEDIST_ADD(x + 1, y , 1.0 - dx, dy);
INVERSEDIST_ADD(x , y + 1, dx, 1.0 - dy);
INVERSEDIST_ADD(x + 1, y + 1, 1.0 - dx, 1.0 - dy);
if( n > 0.0 )
{
return( z / n );
}
}
else
{
long v;
double z[4], n = 0.0, d;
z[0] = z[1] = z[2] = z[3] = 0.0;
INVERSEDIST_ADD_BYTE(x , y , dx, dy);
INVERSEDIST_ADD_BYTE(x + 1, y , 1.0 - dx, dy);
INVERSEDIST_ADD_BYTE(x , y + 1, dx, 1.0 - dy);
INVERSEDIST_ADD_BYTE(x + 1, y + 1, 1.0 - dx, 1.0 - dy);
if( n > 0.0 )
{
z[0] /= n;
z[1] /= n;
z[2] /= n;
z[3] /= n;
return( SG_GET_LONG(z[0], z[1], z[2], z[3]) );
}
}
}
else
{
return( asDouble(x, y) );
}
return( m_NoData_Value );
}
//---------------------------------------------------------
inline double CSG_Grid::_Get_ValAtPos_BiCubicSpline(double dx, double dy, double z_xy[4][4]) const
{
double a0, a2, a3, b1, b2, b3, c[4];
for(int i=0; i<4; i++)
{
a0 = z_xy[0][i] - z_xy[1][i];
a2 = z_xy[2][i] - z_xy[1][i];
a3 = z_xy[3][i] - z_xy[1][i];
b1 = -a0 / 3.0 + a2 - a3 / 6.0;
b2 = a0 / 2.0 + a2 / 2.0;
b3 = -a0 / 6.0 - a2 / 2.0 + a3 / 6.0;
c[i] = z_xy[1][i] + b1 * dx + b2 * dx*dx + b3 * dx*dx*dx;
}
a0 = c[0] - c[1];
a2 = c[2] - c[1];
a3 = c[3] - c[1];
b1 = -a0 / 3.0 + a2 - a3 / 6.0;
b2 = a0 / 2.0 + a2 / 2.0;
b3 = -a0 / 6.0 - a2 / 2.0 + a3 / 6.0;
return( c[1] + b1 * dy + b2 * dy*dy + b3 * dy*dy*dy );
}
inline double CSG_Grid::_Get_ValAtPos_BiCubicSpline(int x, int y, double dx, double dy, bool bByteWise) const
{
if( !bByteWise )
{
double z_xy[4][4];
if( _Get_ValAtPos_Fill4x4Submatrix(x, y, z_xy) )
{
return( _Get_ValAtPos_BiCubicSpline(dx, dy, z_xy) );
}
}
else
{
double z_xy[4][4][4], z[4];
if( _Get_ValAtPos_Fill4x4Submatrix(x, y, z_xy) )
{
z[0] = _Get_ValAtPos_BiCubicSpline(dx, dy, z_xy[0]);
z[1] = _Get_ValAtPos_BiCubicSpline(dx, dy, z_xy[1]);
z[2] = _Get_ValAtPos_BiCubicSpline(dx, dy, z_xy[2]);
z[3] = _Get_ValAtPos_BiCubicSpline(dx, dy, z_xy[3]);
return( SG_GET_LONG(z[0], z[1], z[2], z[3]) );
}
}
return( _Get_ValAtPos_BiLinear(x, y, dx, dy, bByteWise) );
}
//---------------------------------------------------------
inline double CSG_Grid::_Get_ValAtPos_BSpline(double dx, double dy, double z_xy[4][4]) const
{
int i, ix, iy;
double z, px, py, Rx[4], Ry[4];
for(i=0, px=-1.0-dx, py=-1.0-dy; i<4; i++, px++, py++)
{
Rx[i] = 0.0;
Ry[i] = 0.0;
if( (z = px + 2.0) > 0.0 )
Rx[i] += z*z*z;
if( (z = px + 1.0) > 0.0 )
Rx[i] += -4.0 * z*z*z;
if( (z = px + 0.0) > 0.0 )
Rx[i] += 6.0 * z*z*z;
if( (z = px - 1.0) > 0.0 )
Rx[i] += -4.0 * z*z*z;
if( (z = py + 2.0) > 0.0 )
Ry[i] += z*z*z;
if( (z = py + 1.0) > 0.0 )
Ry[i] += -4.0 * z*z*z;
if( (z = py + 0.0) > 0.0 )
Ry[i] += 6.0 * z*z*z;
if( (z = py - 1.0) > 0.0 )
Ry[i] += -4.0 * z*z*z;
Rx[i] /= 6.0;
Ry[i] /= 6.0;
}
for(iy=0, z=0.0; iy<4; iy++)
{
for(ix=0; ix<4; ix++)
{
z += z_xy[ix][iy] * Rx[ix] * Ry[iy];
}
}
return( z );
}
inline double CSG_Grid::_Get_ValAtPos_BSpline(int x, int y, double dx, double dy, bool bByteWise) const
{
if( !bByteWise )
{
double z_xy[4][4];
if( _Get_ValAtPos_Fill4x4Submatrix(x, y, z_xy) )
{
return( _Get_ValAtPos_BSpline(dx, dy, z_xy) );
}
}
else
{
double z_xy[4][4][4], z[4];
if( _Get_ValAtPos_Fill4x4Submatrix(x, y, z_xy) )
{
z[0] = _Get_ValAtPos_BSpline(dx, dy, z_xy[0]);
z[1] = _Get_ValAtPos_BSpline(dx, dy, z_xy[1]);
z[2] = _Get_ValAtPos_BSpline(dx, dy, z_xy[2]);
z[3] = _Get_ValAtPos_BSpline(dx, dy, z_xy[3]);
return( SG_GET_LONG(z[0], z[1], z[2], z[3]) );
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -