Skip to content

Commit

Permalink
Browse files Browse the repository at this point in the history
Added cubic resampler (not working yet)
  • Loading branch information
mhugent committed Dec 21, 2011
1 parent e9b8cc0 commit 3d0b85a
Show file tree
Hide file tree
Showing 6 changed files with 289 additions and 30 deletions.
1 change: 1 addition & 0 deletions src/core/CMakeLists.txt
Expand Up @@ -166,6 +166,7 @@ SET(QGIS_CORE_SRCS

raster/qgsrasterrenderer.cpp
raster/qgsbilinearrasterresampler.cpp
raster/qgscubicrasterresampler.cpp
raster/qgspalettedrasterrenderer.cpp
raster/qgsmultibandcolorrenderer.cpp

Expand Down
51 changes: 22 additions & 29 deletions src/core/raster/qgsbilinearrasterresampler.cpp
Expand Up @@ -43,27 +43,33 @@ void QgsBilinearRasterResampler::resample( const QImage& srcImage, QImage& dstIm
{
currentSrcRow += nSrcPerDstY;
}
else if ( currentSrcRow > srcImage.height() - 2 )
{
//todo: resample in one direction only
break;
}

currentSrcCol = nSrcPerDstY / 2.0 - 0.5;
for ( int j = 0; j < dstImage.width(); ++j )
{
double u = currentSrcCol - ( int )currentSrcCol;
double v = currentSrcRow - ( int )currentSrcRow;
if ( currentSrcCol > srcImage.width() - 2 )
{
//todo: resample in one direction only
//resample in one direction only
px1 = srcImage.pixel( currentSrcCol, currentSrcRow );
px2 = srcImage.pixel( currentSrcCol, currentSrcRow + 1 );
dstImage.setPixel( j, i, resampleColorValue( v, px1, px2 ) );
currentSrcCol += nSrcPerDstX;
continue;
}
else if ( currentSrcRow > srcImage.height() - 2 )
{
px1 = srcImage.pixel( currentSrcCol, currentSrcRow );
px2 = srcImage.pixel( currentSrcCol + 1, currentSrcRow );
dstImage.setPixel( j, i, resampleColorValue( u, px1, px2 ) );
currentSrcCol += nSrcPerDstX;
continue;
}
px1 = srcImage.pixel( currentSrcCol, currentSrcRow );
px2 = srcImage.pixel( currentSrcCol + 1, currentSrcRow );
px3 = srcImage.pixel( currentSrcCol + 1, currentSrcRow + 1 );
px4 = srcImage.pixel( currentSrcCol, currentSrcRow + 1 );
double u = currentSrcCol - ( int )currentSrcCol;
double v = currentSrcRow - ( int )currentSrcRow;
dstImage.setPixel( j, i, resampleColorValue( u, v, px1, px2, px3, px4 ) );
currentSrcCol += nSrcPerDstX;
}
Expand All @@ -73,26 +79,13 @@ void QgsBilinearRasterResampler::resample( const QImage& srcImage, QImage& dstIm

QRgb QgsBilinearRasterResampler::resampleColorValue( double u, double v, QRgb col1, QRgb col2, QRgb col3, QRgb col4 ) const
{
double r1 = qRed( col1 );
double g1 = qGreen( col1 );
double b1 = qBlue( col1 );
double r2 = qRed( col2 );
double g2 = qGreen( col2 );
double b2 = qBlue( col2 );
double r3 = qRed( col3 );
double g3 = qGreen( col3 );
double b3 = qBlue( col3 );
double r4 = qRed( col4 );
double g4 = qGreen( col4 );
double b4 = qBlue( col4 );

double rt1 = u * r2 + ( 1 - u ) * r1;
double gt1 = u * g2 + ( 1 - u ) * g1;
double bt1 = u * b2 + ( 1 - u ) * b1;

double rt2 = u * r3 + ( 1 - u ) * r4;
double gt2 = u * g3 + ( 1 - u ) * g4;
double bt2 = u * b3 + ( 1 - u ) * b4;
int red = bilinearInterpolation( u, v, qRed( col1 ), qRed( col2 ), qRed( col3 ), qRed( col4 ) );
int green = bilinearInterpolation( u, v, qGreen( col1 ), qGreen( col2 ), qGreen( col3 ), qGreen( col4 ) );
int blue = bilinearInterpolation( u, v, qBlue( col1 ), qBlue( col2 ), qBlue( col3 ), qBlue( col4 ) );
return qRgb( red, green, blue );
}

return qRgb( v * rt2 + ( 1 - v ) * rt1, v * gt2 + ( 1 - v ) * gt1, v * bt2 + ( 1 - v ) * bt1 );
QRgb QgsBilinearRasterResampler::resampleColorValue( double u, QRgb col1, QRgb col2 ) const
{
return qRgb( qRed( col1 ) * ( 1 - u ) + qRed( col2 ) * u, qGreen( col1 ) * ( 1 - u ) + qGreen( col2 ) * u, qBlue( col1 ) * ( 1 - u ) + qBlue( col2 ) * u );
}
7 changes: 7 additions & 0 deletions src/core/raster/qgsbilinearrasterresampler.h
Expand Up @@ -31,6 +31,13 @@ class QgsBilinearRasterResampler: public QgsRasterResampler

private:
QRgb resampleColorValue( double u, double v, QRgb col1, QRgb col2, QRgb col3, QRgb col4 ) const;
QRgb resampleColorValue( double u, QRgb col1, QRgb col2 ) const;
inline double bilinearInterpolation( double u, double v, int val1, int val2, int val3, int val4 ) const;
};

double QgsBilinearRasterResampler::bilinearInterpolation( double u, double v, int val1, int val2, int val3, int val4 ) const
{
return ( val1 * ( 1 - u ) * ( 1 - v ) + val2 * ( 1 - v ) * u + val4 * ( 1 - u ) * v + val3 * u * v );
}

#endif // QGSBILINEARRASTERRESAMPLER_H
192 changes: 192 additions & 0 deletions src/core/raster/qgscubicrasterresampler.cpp
@@ -0,0 +1,192 @@
/***************************************************************************
qgscubicrasterresampler.cpp
----------------------------
begin : December 2011
copyright : (C) 2011 by Marco Hugentobler
email : marco at sourcepole dot ch
***************************************************************************/

/***************************************************************************
* *
* This program 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; either version 2 of the License, or *
* (at your option) any later version. *
* *
***************************************************************************/

#include "qgscubicrasterresampler.h"
#include <QImage>

QgsCubicRasterResampler::QgsCubicRasterResampler()
{
}

QgsCubicRasterResampler::~QgsCubicRasterResampler()
{
}

void QgsCubicRasterResampler::resample( const QImage& srcImage, QImage& dstImage ) const
{
int nCols = srcImage.width();
int nRows = srcImage.height();

int pos = 0;
QRgb px;
int* redMatrix = new int[ nCols * nRows ];
int* greenMatrix = new int[ nCols * nRows ];
int* blueMatrix = new int[ nCols * nRows ];

for( int i = 0; i < nRows; ++i )
{
for( int j = 0; j < nCols; ++j )
{
px = srcImage.pixel( j, i );
redMatrix[pos] = qRed( px );
greenMatrix[pos] = qGreen( px );
blueMatrix[pos] = qBlue( px );
++pos;
}
}

//derivative x
double* xDerivativeMatrixRed = new double[ nCols * nRows ];
xDerivativeMatrix( nCols, nRows, xDerivativeMatrixRed, redMatrix );
double* xDerivativeMatrixGreen = new double[ nCols * nRows ];
xDerivativeMatrix( nCols, nRows, xDerivativeMatrixGreen, greenMatrix );
double* xDerivativeMatrixBlue = new double[ nCols * nRows ];
xDerivativeMatrix( nCols, nRows, xDerivativeMatrixBlue, blueMatrix );

//derivative y
double* yDerivativeMatrixRed = new double[ nCols * nRows ];
yDerivativeMatrix<int>( nCols, nRows, yDerivativeMatrixRed, redMatrix );
double* yDerivativeMatrixGreen = new double[ nCols * nRows ];
yDerivativeMatrix<int>( nCols, nRows, yDerivativeMatrixGreen, greenMatrix );
double* yDerivativeMatrixBlue = new double[ nCols * nRows ];
yDerivativeMatrix<int>( nCols, nRows, yDerivativeMatrixBlue, blueMatrix );

//derivative xy
double* xyDerivativeMatrixRed = new double[ nCols * nRows ];
yDerivativeMatrix<double>( nCols, nRows, xyDerivativeMatrixRed, xDerivativeMatrixRed );
double* xyDerivativeMatrixGreen = new double[ nCols * nRows ];
yDerivativeMatrix<double>( nCols, nRows, xyDerivativeMatrixGreen, xDerivativeMatrixGreen );
double* xyDerivativeMatrixBlue = new double[ nCols * nRows ];
yDerivativeMatrix<double>( nCols, nRows, xyDerivativeMatrixBlue, xDerivativeMatrixBlue );

//compute output
double nSrcPerDstX = ( double ) srcImage.width() / ( double ) dstImage.width();
double nSrcPerDstY = ( double ) srcImage.height() / ( double ) dstImage.height();

double currentSrcRow = nSrcPerDstX / 2.0 - 0.5;
double currentSrcCol;
int currentSrcColInt;
int currentSrcRowInt;

int r, g, b;
int index;

for ( int i = 0; i < dstImage.height(); ++i )
{
currentSrcRowInt = (int) currentSrcRow;

currentSrcCol = nSrcPerDstY / 2.0 - 0.5;
for( int j = 0; j < dstImage.width(); ++j )
{
//out of bounds check
currentSrcColInt = (int)currentSrcCol;
if( currentSrcColInt < 0 || currentSrcColInt > srcImage.width() - 2
|| currentSrcRowInt < 0 || currentSrcRowInt > srcImage.height() - 2 )
{
currentSrcCol += nSrcPerDstX;
++index;
dstImage.setPixel( j, i, qRgb( 0, 0, 255 ) );
continue;
}


//interpolation
index = currentSrcRowInt * nCols + currentSrcColInt;
r = cubicInterpolation( redMatrix[index], redMatrix[index + 1], redMatrix[index + nCols + 1 ], redMatrix[ index + nCols ],
xDerivativeMatrixRed[index], xDerivativeMatrixRed[index + 1], xDerivativeMatrixRed[index + nCols + 1 ], xDerivativeMatrixRed[ index + nCols ],
yDerivativeMatrixRed[index], yDerivativeMatrixRed[index + 1], yDerivativeMatrixRed[index + nCols + 1 ], yDerivativeMatrixRed[ index + nCols ],
xyDerivativeMatrixRed[index], xyDerivativeMatrixRed[index + 1], xyDerivativeMatrixRed[index + nCols + 1 ], xyDerivativeMatrixRed[ index + nCols ] );
g = cubicInterpolation( greenMatrix[index], greenMatrix[index + 1], greenMatrix[index + nCols + 1 ], greenMatrix[ index + nCols ],
xDerivativeMatrixGreen[index], xDerivativeMatrixGreen[index + 1], xDerivativeMatrixGreen[index + nCols + 1 ], xDerivativeMatrixGreen[ index + nCols ],
yDerivativeMatrixGreen[index], yDerivativeMatrixGreen[index + 1], yDerivativeMatrixGreen[index + nCols + 1 ], yDerivativeMatrixGreen[ index + nCols ],
xyDerivativeMatrixGreen[index], xyDerivativeMatrixGreen[index + 1], xyDerivativeMatrixGreen[index + nCols + 1 ], xyDerivativeMatrixGreen[ index + nCols ] );
b = cubicInterpolation( blueMatrix[index], blueMatrix[index + 1], blueMatrix[index + nCols + 1 ], blueMatrix[ index + nCols ],
xDerivativeMatrixBlue[index], xDerivativeMatrixBlue[index + 1], xDerivativeMatrixBlue[index + nCols + 1 ], xDerivativeMatrixBlue[ index + nCols ],
yDerivativeMatrixBlue[index], yDerivativeMatrixBlue[index + 1], yDerivativeMatrixBlue[index + nCols + 1 ], yDerivativeMatrixBlue[ index + nCols ],
xyDerivativeMatrixBlue[index], xyDerivativeMatrixBlue[index + 1], xyDerivativeMatrixBlue[index + nCols + 1 ], xyDerivativeMatrixBlue[ index + nCols ] );

dstImage.setPixel( j, i, qRgb( 255, 0, 0 ) );
currentSrcCol += nSrcPerDstX;
}
currentSrcRow += nSrcPerDstY;
}


//cleanup memory
delete[] redMatrix;
delete[] greenMatrix;
delete[] blueMatrix;
delete[] xDerivativeMatrixRed;
delete[] xDerivativeMatrixGreen;
delete[] xDerivativeMatrixBlue;
delete[] yDerivativeMatrixRed;
delete[] yDerivativeMatrixGreen;
delete[] yDerivativeMatrixBlue;
delete[] xyDerivativeMatrixRed;
delete[] xyDerivativeMatrixGreen;
delete[] xyDerivativeMatrixBlue;
}

void QgsCubicRasterResampler::xDerivativeMatrix( int nCols, int nRows, double* matrix, const int* colorMatrix )
{
double val;
int index = 0;

for( int i = 0; i < nRows; ++i )
{
for( int j = 0; j < nCols; ++j )
{
if( j < 1 )
{
val = colorMatrix[index + 1] - colorMatrix[index];
}
else if( j == ( nCols - 1 ) )
{
val = colorMatrix[index] - colorMatrix[ index - 1 ];
}
else
{
val = ( colorMatrix[index + 1] - colorMatrix[index - 1] ) / 2.0;
}
matrix[index] = val;
++index;
}
}
}

int QgsCubicRasterResampler::cubicInterpolation( double p1, double p2, double p3, double p4, double p1x, double p2x, double p3x, double p4x,
double p1y, double p2y, double p3y, double p4y, double p1xy, double p2xy, double p3xy, double p4xy )
{
#if 0
//calculate 16 coefficients
double a00 = p1 - 3 * p4 + 2* p3 - 3 * p1y + 9 * p4y - 6 * p3y + 2 * p1xy - 6 * p4xy+ 4 * p3xy;
double a10 = 3 * p4 - 2 * p3 - 9* p4y + 6 * p3y + 6 * p4xy - 4 * p3xy;
double a20 = 3 * p1y - 9 * p4y + 6 * p3y - 2 * p1xy + 6 * p4xy - 4 * p3xy;
double a30 = 9 * p4y - 6 * p3y - 6 * p4xy + 4 * p3xy;
double a01 = p2 - 2 * p4 + p3 - 3 * p2y + 6 * p4y -3 * p3y + 2 * p2xy - 4 * p4xy + 2 * p3xy;
double a11 = - p4 + p3 + 3 * p4y + -3 * p3y - 2 * p4xy;
double a21 = 3 * p2y - 6 * p4y + 3 * p3y -2 * p2xy + 4 * p4xy - 2 * p3xy;
double a31 = - 3 * p4y + 3 * p3y + 2 * p4xy -2 * p3xy;
double a02 = p1x - 3 * p4x + 2 * p3x - 2 * p1y + 6 * p4y - 4 * p3y + p1xy - 3 * p4xy + 2 * p3xy;
double a12 = 3 * p4x - 2 * p3x - 6 * p4y + 4 * p3y + 3 * p4xy - 2 * p3xy;
double a22 = - p1y + 3 * p4y - 2 * p3y + p1xy - 3 * p4xy + 2 * p3xy;
double a32 = -3 * p4y + 2 * p3y + 3 * p4xy - 2 * p3xy;
double a03 = p2x - 2 * p4x + p3x - 2 * p4y + 4 * p4y - 2 * p3y + p2xy - 2 * p4xy + p3xy;
double a13 = - p4x + p3x + 2 * p4y - 2 * p3y - p4xy + p3xy;
#endif //0
return 0;
}
64 changes: 64 additions & 0 deletions src/core/raster/qgscubicrasterresampler.h
@@ -0,0 +1,64 @@
/***************************************************************************
qgscubicrasterresampler.h
----------------------------
begin : December 2011
copyright : (C) 2011 by Marco Hugentobler
email : marco at sourcepole dot ch
***************************************************************************/

/***************************************************************************
* *
* This program 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; either version 2 of the License, or *
* (at your option) any later version. *
* *
***************************************************************************/

#ifndef QGSCUBICRASTERRESAMPLER_H
#define QGSCUBICRASTERRESAMPLER_H

#include "qgsrasterresampler.h"

class QgsCubicRasterResampler: public QgsRasterResampler
{
public:
QgsCubicRasterResampler();
~QgsCubicRasterResampler();
void resample( const QImage& srcImage, QImage& dstImage ) const;

private:
static void xDerivativeMatrix( int nCols, int nRows, double* matrix, const int* colorMatrix );
template <class T> void yDerivativeMatrix( int nCols, int nRows, double* dstMatrix, const T* srcMatrix ) const;
static int cubicInterpolation( double p1, double p2, double p3, double p4, double p1x, double p2x, double p3x, double p4x,
double p1y, double p2y, double p3y, double p4y, double p1xy, double p2xy, double p3xy, double p4xy );
};

template <class T> void QgsCubicRasterResampler::yDerivativeMatrix( int nCols, int nRows, double* dstMatrix, const T* srcMatrix ) const
{
double val;
int index = 0;

for( int i = 0; i < nRows; ++i )
{
for( int j = 0; j < nCols; ++j )
{
if( i == 0 )
{
val = srcMatrix[ index + nRows ] - srcMatrix[ index ];
}
else if( i == nRows - 1 )
{
val = srcMatrix[ index ] - srcMatrix[ index - nRows ];
}
else
{
val = ( srcMatrix[ index + nRows ] - srcMatrix[ index - nRows ] ) / 2.0;
}
dstMatrix[index] = val;
++index;
}
}
}

#endif // QGSCUBICRASTERRESAMPLER_H
4 changes: 3 additions & 1 deletion src/core/raster/qgsrasterlayer.cpp
Expand Up @@ -39,6 +39,7 @@ email : tim at linfiniti.com
//renderers
#include "qgspalettedrasterrenderer.h"
#include "qgsbilinearrasterresampler.h"
#include "qgscubicrasterresampler.h"
#include "qgsmultibandcolorrenderer.h"

#include <cstdio>
Expand Down Expand Up @@ -848,7 +849,8 @@ void QgsRasterLayer::draw( QPainter * theQPainter,
colorArray[( int )colorIt->value] = colorIt->color;
}

QgsBilinearRasterResampler resampler;
//QgsBilinearRasterResampler resampler;
QgsCubicRasterResampler resampler;
QgsPalettedRasterRenderer renderer( mDataProvider, bNumber, colorArray, itemList.size(), &resampler );
renderer.draw( theQPainter, theRasterViewPort, theQgsMapToPixel );
#if 0
Expand Down

0 comments on commit 3d0b85a

Please sign in to comment.