00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028 #include "GaborJet.h"
00029 #include "CelesteGlobals.h"
00030
00031 #ifndef M_PI
00032 #define M_PI 3.1415926535897932384626433832795
00033 #endif
00034
00035 using namespace std;
00036
00037 namespace celeste
00038 {
00039
00040 GaborJet::GaborJet()
00041 {
00042 mHeight = 512;
00043 mWidth = 512;
00044 mX = 128;
00045 mY = 128;
00046 mShowFilter = false;
00047 mFilters = NULL;
00048 mFiducials = NULL;
00049 }
00050
00051
00052 GaborJet::~GaborJet()
00053 {
00054 if ( mFilters != NULL )
00055 {
00056 for ( int i = 0; i < mAngles; i++ ) delete[] mFilters[i];
00057 delete[] mFilters;
00058 }
00059 if ( mFiducials != NULL ) delete[] mFiducials;
00060 }
00061
00062
00063
00064 void GaborJet::Initialize( int y, int x, int x0, int y0, int r,
00065 float s, int f, float maxF, float minF, int a, bool save )
00066 {
00067 int i, j;
00068 float angle, freq;
00069
00070
00071 mHeight = y;
00072 mWidth = x;
00073 mX = x0;
00074 mY = y0;
00075 mSigma = (float)(s * M_PI * M_PI);
00076 mAngles = a;
00077 mFreqs = f;
00078 mRadius = r;
00079 mMinFreq = minF;
00080 mMaxFreq = maxF;
00081 mShowFilter = save;
00082 mFiducials = new float[mAngles * mFreqs];
00083
00084
00085 mFilters = new GaborFilter * [mAngles];
00086 for ( i = 0; i < mAngles; i++ )
00087 {
00088
00089 angle = (float)((float)i * M_PI / (float)mAngles);
00090
00091
00092 mFilters[i] = new GaborFilter[mFreqs];
00093
00094
00095 for ( j = 0; j < mFreqs; j++ )
00096 {
00097
00098 freq = minF + ( j * ( maxF - minF ) ) / (float)mFreqs;
00099
00100
00101 mFilters[i][j].Initialize( mRadius, angle, freq, mSigma );
00102 if ( mShowFilter ) mFilters[i][j].Save( mFile, i, j );
00103 }
00104 }
00105 }
00106
00107
00108
00109 void GaborJet::Filter( float** image, int* len )
00110 {
00111 int x, y;
00112 int gx, gy;
00113 int a, f;
00114 int h, i, j;
00115 float sumI, sumR;
00116
00117 if ( kVerbosity ) cerr << "convoluting..." << endl;
00118
00119
00120
00121 h = 0;
00122 for ( a = 0; a < mAngles; a++ )
00123 {
00124 for ( f = 0; f < mFreqs; f++ )
00125 {
00126 sumR = 0.0;
00127 sumI = 0.0;
00128
00129
00130 y = mY - mRadius;
00131 for ( gy = y; gy < y + 2 * mRadius; gy++ )
00132 {
00133
00134 if ( gy < 0 || gy >= mHeight ) break;
00135
00136
00137 i = gy - y;
00138
00139 x = mX - mRadius;
00140 for ( gx = x; gx < x + 2 * mRadius; gx++ )
00141 {
00142
00143 if ( gx < 0 || gx >= mWidth ) break;
00144
00145
00146 j = gx - x;
00147
00148 sumR += image[gy][gx] * mFilters[a][f].GetReal(i,j);
00149 sumI += image[gy][gx] * mFilters[a][f].GetImaginary(i,j);
00150 }
00151 }
00152 mFiducials[h] = sqrt( sumR*sumR + sumI*sumI );
00153 h++;
00154 }
00155 }
00156
00157 *len = mAngles * mFreqs;
00158 }
00159
00160 };