Main Page | Namespace List | Class Hierarchy | Class List | Directories | File List | Class Members | File Members

edgeDetect.cpp

Go to the documentation of this file.
00001 //==============================================
00002 //  copyright            : (C) 2003-2005 by Will Stokes
00003 //==============================================
00004 //  This program is free software; you can redistribute it 
00005 //  and/or modify it under the terms of the GNU General 
00006 //  Public License as published by the Free Software 
00007 //  Foundation; either version 2 of the License, or  
00008 //  (at your option) any later version.         
00009 //==============================================
00010 
00011 //Systemwide includes
00012 #include <qpainter.h>
00013 #include <qimage.h>
00014 #include <math.h>
00015 
00016 //Projectwide includes
00017 #include "edgeDetect.h"
00018 #include "blur.h"
00019 #include "../enhancements/contrast.h"
00020 
00021 //==============================================
00022 EdgeDetect::EdgeDetect( QImage* image )
00023 {                     
00024   //load image
00025   this->image = image;
00026 
00027   //allocate and initialize objects used for edge detection
00028   allocateAndInitObjects();
00029 
00030   //fill lum map and lum histogram
00031   fillLumMapAndLumHistogram();
00032   
00033   //fill smoothed histogram
00034   smoothLumHistogram();
00035   
00036   //compute edge magnitude and GSLC maps
00037   computeEdgeMagAndGSLCmaps();
00038   
00039   //determine pixel clusters
00040   findPixelClusters();
00041   
00042   computeClusterStatistics();  
00043  
00044   computeClusterThresholds();
00045   
00046   constructEdgeImage();
00047 }
00048 //==============================================
00049 EdgeDetect::~EdgeDetect()
00050 {
00051   deallocateObjects();
00052 }
00053 //==============================================
00054 int EdgeDetect::getNumClusters()
00055 { return numClusters; }
00056 //==============================================
00057 PixelCluster* EdgeDetect::getClusters()
00058 { return clusters; }
00059 //==============================================
00060 int* EdgeDetect::getPeaks()
00061 { return clusterPeaks; }
00062 //==============================================
00063 int* EdgeDetect::getSmoothHist()
00064 { return smoothLumHist; }
00065 //==============================================
00066 QImage* EdgeDetect::getEdgeImage()
00067 {
00068   return image; 
00069 }
00070 //==============================================
00071 int* EdgeDetect::getClusterMap()
00072 {
00073   //construct map
00074   int* clusterMap = new int[image->width() * image->height()];
00075   
00076   //iterate over all pixels, determine cluster each pixel belongs to
00077   int i, cluster;
00078   for(i=0; i<image->width()*image->height(); i++)
00079   {
00080     for(cluster=0; cluster<numClusters; cluster++)
00081     {
00082       if( lumMap[i] >= clusters[cluster].minLuminance &&
00083           lumMap[i] <= clusters[cluster].maxLuminance )
00084       {
00085         clusterMap[i] = cluster;
00086         break;
00087       }
00088     } //cluster
00089   } //pixel
00090 
00091   return clusterMap;
00092 }
00093 //==============================================
00094 void EdgeDetect::allocateAndInitObjects()
00095 {
00096   //initialize: 
00097   //-luminosity histogram
00098   //-smoothed luminosity histogram
00099   //-identified peak regions
00100   int i;
00101   for(i=0; i<256; i++)
00102   { 
00103     lumHist[i] = 0; 
00104     smoothLumHist[i] = 0;
00105     clusterPeaks[i] = 0;
00106   }
00107   
00108   //allocate luminance map
00109   lumMap = new int[image->width() * image->height()];
00110   
00111   //allocate edge magnitude map
00112   edgeMagMap = new float[image->width() * image->height()];
00113   
00114   //allocate GSLC map
00115   GSLCmap = new int[image->width() * image->height()];
00116   
00117   //construct LUT
00118   constructGSLClut();
00119 }
00120 //==============================================
00121 void EdgeDetect::fillLumMapAndLumHistogram()
00122 {
00123   int x, y;
00124   QRgb* rgb;
00125   uchar* scanLine;
00126   int lumVal;
00127   for( y=0; y<image->height(); y++)
00128   {   
00129     scanLine = image->scanLine(y);
00130     for( x=0; x<image->width(); x++)
00131     {
00132       //get lum value for this pixel
00133       rgb = ((QRgb*)scanLine+x);
00134       lumVal = qGray(*rgb);
00135       
00136       //store in lum map
00137       lumMap[x + y*image->width()] = lumVal;
00138       
00139       //update lum histogram
00140       lumHist[ lumVal ]++;
00141     }
00142   }
00143 }
00144 //==============================================
00145 void EdgeDetect::smoothLumHistogram()
00146 {
00147   #define FILTER_SIZE 5
00148   int filter[FILTER_SIZE] = {2, 5, 8, 5, 2};
00149   
00150   int i,j;
00151   int filterIndex, sum, total;
00152   for(i = 0; i<256; i++)
00153   {
00154     sum = 0;
00155     total = 0;
00156     
00157     for( j= -FILTER_SIZE/2; j <= FILTER_SIZE/2; j++)
00158     {
00159       if( i+j > 0 && i+j < 256 )
00160       {
00161         filterIndex = j+ FILTER_SIZE/2;
00162         total+= filter[filterIndex] * lumHist[i+j];
00163         sum  += filter[filterIndex];
00164       }
00165     }
00166     
00167     smoothLumHist[i] = total / sum;
00168   }
00169 }
00170 //==============================================
00171 void EdgeDetect::computeEdgeMagAndGSLCmaps()
00172 {
00173   int x, y;
00174   int idealPattern[9];
00175   int pixelLums[9];
00176   
00177   //-------  
00178   //iterate over all pixels
00179   for( y=0; y<image->height(); y++)
00180   {   
00181     for( x=0; x<image->width(); x++)
00182     {
00183       //compute pixel luminances for entire grid
00184       pixelLums[0] = pixelLum(x-1,y-1);
00185       pixelLums[1] = pixelLum(x  ,y-1);
00186       pixelLums[2] = pixelLum(x+1,y-1);
00187       pixelLums[3] = pixelLum(x-1,y  );
00188       pixelLums[4] = pixelLum(x  ,y  );
00189       pixelLums[5] = pixelLum(x+1,y  );
00190       pixelLums[6] = pixelLum(x-1,y+1);
00191       pixelLums[7] = pixelLum(x  ,y+1);
00192       pixelLums[8] = pixelLum(x+1,y+1);
00193       
00194       //compute average
00195       float avg = 0;
00196       int i;
00197       for(i=0; i<=8; i++)
00198       {
00199         avg+= pixelLums[i];
00200       }
00201       avg = avg / 9;
00202       
00203       //determine ideal pattern and I0 and I1 averages
00204       int centerPixelLum = pixelLums[4];
00205       float centerDiff = centerPixelLum - avg;
00206       
00207       float I0avg = 0;
00208       int I0count = 0;
00209       
00210       float I1avg = 0;
00211       int I1count = 0;
00212       
00213       for(i=0; i<=8; i++)
00214       {
00215         if( centerDiff * (pixelLums[i]-avg) >=0 )
00216         { 
00217           I1avg+=pixelLums[i];
00218           I1count++;
00219           idealPattern[i] = 1; 
00220         }
00221         else 
00222         { 
00223           I0avg+=pixelLums[i];
00224           I0count++;
00225           idealPattern[i] = 0; 
00226         }
00227       }
00228 
00229       //compute and store edge magnitude
00230       if(I0count > 0) I0avg = I0avg/I0count;
00231       if(I1count > 0) I1avg = I1avg/I1count;     
00232       edgeMagMap[x + y*image->width()] = QABS( I1avg - I0avg );
00233       
00234       //compute and store GSLC
00235       int GSLC=0;
00236       int weight = 1;
00237       for(i=0; i<9; i++)
00238       {
00239         //skip center
00240         if(i == 4) continue;
00241         
00242         if(idealPattern[i] == 1)
00243         { GSLC+=weight; }
00244         
00245         weight = weight*2;
00246       }
00247       GSLCmap[x + y*image->width()] = GSLC;
00248     } //x
00249   } //y
00250 }
00251 //==============================================
00252 int EdgeDetect::pixelLum(int x, int y)
00253 {
00254   int clampedX = QMAX( QMIN( x, image->width()-1), 0);
00255   int clampedY = QMAX( QMIN( y, image->height()-1), 0);
00256   return lumMap[ clampedX + clampedY * image->width() ];
00257 }
00258 //==============================================
00259 void EdgeDetect::findPixelClusters()
00260 {
00261   //find max count
00262   int maxCount = 0;
00263   int i;
00264   for(i=0; i<256; i++)
00265   {
00266     if(smoothLumHist[i] > maxCount)
00267       maxCount = smoothLumHist[i];
00268   }
00269 
00270   //compute JND for histogram (2% of total spread)
00271   int histJND = maxCount/50;
00272 
00273   //construct temporary array for valley locations
00274   //1's will indicate a valley midpoint
00275   int tmpValleyArray[256];
00276   for(i=0; i<256; i++) { tmpValleyArray[i] = 0; }
00277   
00278   //move across histogram finding valley midpoints
00279   int curTrackedMin = smoothLumHist[0];
00280   
00281   //first and last indices tracked min was observed
00282   int firstMinIndex = 0;
00283   int lastMinIndex = 0;
00284   
00285   //only add valley midpoint if finished tracking a descent
00286   bool slopeNeg = false;
00287   
00288   for(i = 1; i<256; i++ )
00289   {
00290     if( smoothLumHist[i] < curTrackedMin - histJND )
00291     {
00292       //found a descent!
00293       slopeNeg = true;
00294       curTrackedMin = smoothLumHist[i];
00295       firstMinIndex = i;
00296     }
00297     //starting to go up again, add last min to list
00298     else if( smoothLumHist[i] > curTrackedMin + histJND )
00299     {
00300       //if finished tracing a negative slope find midpoint and set location to true
00301       if(slopeNeg)
00302       {
00303         tmpValleyArray[ (firstMinIndex + lastMinIndex)/2 ] = 1;
00304       }
00305       
00306       curTrackedMin = smoothLumHist[i];
00307       slopeNeg = false;
00308     }
00309     else
00310     {
00311       //still tracking a min, update the right 
00312       //hand index. center of valley is found
00313       //by averaging first and last min index
00314       lastMinIndex = i;
00315     }
00316   }
00317   
00318   //count valleys
00319   int numValleys = 0;
00320   for(i=0; i<256; i++)
00321   {
00322     if(tmpValleyArray[i] == 1 ) numValleys++;
00323   }
00324 
00325   //determine number of clusters
00326   numClusters = numValleys-1;
00327   if(tmpValleyArray[0] != 1)
00328     numClusters++;
00329   if(tmpValleyArray[255] != 1)
00330     numClusters++;
00331   
00332   //allocate clusters
00333   clusters = new PixelCluster[numClusters];
00334   
00335   //automatically start first cluster
00336   int cluster=0;
00337   clusters[cluster].minLuminance = 0;
00338   
00339   //initialize left and right boundaries of all clusters
00340   for(i=1; i<256; i++)
00341   {
00342     //reached next valley, end cluster
00343     if( tmpValleyArray[i] == 1)
00344     {
00345       clusters[cluster].maxLuminance = i-1;
00346       cluster++;
00347       clusters[cluster].minLuminance = i;
00348     }
00349     //end last cluster automatically at end
00350     else if(i == 255)
00351     {
00352       clusters[cluster].maxLuminance = i;
00353     }
00354   }
00355   
00356   //determine cluster peaks
00357   for(cluster=0; cluster<numClusters; cluster++)
00358   {
00359     //find max for current cluster
00360     int maxIndex = clusters[cluster].minLuminance;
00361     for(i=clusters[cluster].minLuminance; i<=clusters[cluster].maxLuminance; i++)
00362     {
00363       if(smoothLumHist[i] > smoothLumHist[maxIndex])
00364          maxIndex = i;
00365     }
00366     
00367     //mark peaks  
00368     int lumJND = 255/50;
00369     for(i=QMAX(0, maxIndex-lumJND); i<QMIN(256, maxIndex+lumJND); i++)
00370     { 
00371       clusterPeaks[i] = 1; 
00372     }
00373   }
00374 }
00375 //==============================================
00376 void EdgeDetect::computeClusterStatistics()
00377 {
00378   //initialize cluster stats
00379   int cluster;
00380   for(cluster=0; cluster<numClusters; cluster++)
00381   {
00382     int i;
00383     for(i=0; i<256; i++)
00384     {
00385       clusters[cluster].edgeMagHistogram[i] = 0;
00386     }
00387     clusters[cluster].totalEdgeMagnitude=0.0f;
00388     clusters[cluster].numPixels = 0;
00389   }
00390   
00391   //iterate over all pixels
00392   int i;
00393   for(i=0; i<image->width()*image->height(); i++)
00394   {
00395     //skip pixels that don't belong to peaks
00396     if( clusterPeaks[ lumMap[i] ] != 1)
00397       continue;
00398     
00399     //determine cluster pixel belongs to
00400     int cluster;
00401     for(cluster=0; cluster<numClusters; cluster++)
00402     {
00403       if( lumMap[i] >= clusters[cluster].minLuminance &&
00404           lumMap[i] <= clusters[cluster].maxLuminance )
00405       {      
00406         clusters[cluster].totalEdgeMagnitude+= edgeMagMap[i]; 
00407         clusters[cluster].numPixels++;
00408         clusters[cluster].edgeMagHistogram[ QMIN( QMAX( (int)edgeMagMap[i], 0), 255) ]++;
00409         break;
00410       }
00411     } //cluster
00412   } //pixel i
00413   
00414   //iterate over clusters to determine min and max peak cluster sizes
00415   minClusterSize = clusters[0].numPixels;
00416   maxClusterSize = clusters[0].numPixels;
00417   for(cluster=1; cluster<numClusters; cluster++)
00418   {
00419     if(clusters[cluster].numPixels < minClusterSize)
00420       minClusterSize = clusters[cluster].numPixels;
00421 
00422     if(clusters[cluster].numPixels > maxClusterSize)
00423       maxClusterSize = clusters[cluster].numPixels;
00424   }
00425 
00426   //iterate over clusters one final time to deduce normalized inputs to fuzzy logic process
00427   int JND = 255/50;
00428   for(cluster=0; cluster<numClusters; cluster++)
00429   {
00430     clusters[cluster].meanMode = QMIN( clusters[cluster].totalEdgeMagnitude / clusters[cluster].numPixels,
00431                                        3*JND );
00432     
00433     int i;
00434     int mode = 0;
00435     for(i=1; i<256; i++)
00436     {
00437       if( clusters[cluster].edgeMagHistogram[i] > clusters[cluster].edgeMagHistogram[ mode ] )
00438         mode = i;
00439     }
00440     clusters[cluster].mode = QMIN( mode, 2*JND );
00441         
00442     clusters[cluster].pixelCount = ((float)(clusters[cluster].numPixels - minClusterSize)) / 
00443                                    (maxClusterSize - minClusterSize);  
00444   }
00445 }
00446 //==============================================
00447 //compute edge thresholds for each cluster using 18-rule fuzzy logic approach
00448 void EdgeDetect::computeClusterThresholds()
00449 {
00450   //iterate over each cluster
00451   int cluster;
00452   float S1,M1,L1;
00453   float S2,M2,L2;
00454   float S3,L3;
00455   float outS, outM, outL;
00456   
00457   int JND = 255/50;
00458   
00459   for(cluster=0; cluster<numClusters; cluster++)
00460   {
00461     //----
00462     //compute S,M, and L values for each input
00463     //----
00464     S1 = QMAX( 1.0f - ((clusters[cluster].meanMode/JND) / 1.5f), 0 );
00465 
00466     if( (clusters[cluster].meanMode/JND) <= 1.5f )
00467       M1 = QMAX( (clusters[cluster].meanMode/JND) - 0.5f, 0 );
00468     else
00469       M1 = QMAX( 2.5f - (clusters[cluster].meanMode/JND), 0 );
00470     
00471     L1 = QMAX( ((clusters[cluster].meanMode/JND) - 1.5f) / 1.5f, 0 );
00472     //----
00473     S2 = QMAX( 1.0f - (clusters[cluster].mode/JND), 0 );
00474     
00475     if( (clusters[cluster].mode/JND) <= 1.0f )
00476       M2 = QMAX( -1.0f + 2*(clusters[cluster].mode/JND), 0 );
00477     else
00478       M2 = QMAX( 3.0f - 2*(clusters[cluster].mode/JND), 0 );
00479     
00480     L2 = QMAX( (clusters[cluster].mode/JND) - 1.0, 0 );
00481     //----
00482     S3 = QMAX( 1.0f - 2*clusters[cluster].pixelCount, 0 );
00483     L3 = QMAX( -1.0f + 2*clusters[cluster].pixelCount, 0 );
00484     //----
00485     
00486     //Compute M,L for outputs using set of 18 rules.
00487     //outS is inherantly S given the ruleset provided
00488     outM = 0.0f;
00489     outL = 0.0f;
00490     //Out 1
00491     if( numClusters > 2 )
00492     {
00493       outM += S1*S2*S3;   //rule 1
00494       
00495       //rule 2
00496       if( clusters[cluster].meanMode < clusters[cluster].mode )
00497         outS += S1*S2*L3;   
00498       else
00499         outM += S1*S2*L3;
00500 
00501       outM += S1*M2*S3;   //rule 3
00502       outM += S1*M2*L3;   //rule 4
00503       outM += S1*L2*S3;   //rule 5
00504       outM += S1*L2*L3;   //rule 6
00505       outM += M1*S2*S3;   //rule 7
00506       outM += M1*S2*L3;   //rule 8
00507       outM += M1*M2*S3;   //rule 9
00508       outL += M1*M2*L3;   //rule 10
00509       outM += M1*L2*S3;   //rule 11
00510       outL += M1*L2*L3;   //rule 12
00511       outM += L1*S2*S3;   //rule 13
00512       outL += L1*S2*L3;   //rule 14
00513       outM += L1*M2*S3;   //rule 15
00514       outL += L1*M2*L3;   //rule 16
00515       outL += L1*L2*S3;   //rule 17
00516       outL += L1*L2*L3;   //rule 18
00517     }
00518     //Out 2
00519     else
00520     {
00521       outL += S1*S2*S3;   //rule 1
00522       outL += S1*S2*L3;   //rule 2
00523       outM += S1*M2*S3;   //rule 3
00524       outL += S1*M2*L3;   //rule 4
00525       outM += S1*L2*S3;   //rule 5
00526       outM += S1*L2*L3;   //rule 6
00527       outL += M1*S2*S3;   //rule 7
00528       outL += M1*S2*L3;   //rule 8
00529       outL += M1*M2*S3;   //rule 9
00530       outL += M1*M2*L3;   //rule 10
00531       outL += M1*L2*S3;   //rule 11
00532       outL += M1*L2*L3;   //rule 12
00533       outL += L1*S2*S3;   //rule 13
00534       outL += L1*S2*L3;   //rule 14
00535       outL += L1*M2*S3;   //rule 15
00536       outL += L1*M2*L3;   //rule 16
00537       outL += L1*L2*S3;   //rule 17
00538       outL += L1*L2*L3;   //rule 18
00539     }
00540 
00541     //find centroid - Beta[k]
00542     float A = outM + 0.5f;
00543     float B = 2.5f - outM;
00544     float C = 1.5f * (outL + 1);
00545     float D = 1.5f * (outM + 1);
00546     float E = 2.5f - outL;
00547 
00548     //---------------------------------------------------------------
00549     //Case 1: Both outM and outL are above intersection point of diagonals
00550     if( outM > 0.5f && outL > 0.5f )
00551     {
00552       //find area of 7 subregions
00553       float area1 = ((A-0.5f)*outM)/2;
00554       float area2 = outM * (B-A);
00555       float area3 = ((2.1f-B) * (outM - 0.5)) / 2;
00556       float area4 = (2.1 - B) * 0.5f;
00557       float area5 = ((C - 2.1f) * (outL - 0.5)) / 2;
00558       float area6 = (C - 2.1f) * 0.5f;
00559       float area7 = (3.0f - C) * outL;
00560      
00561       //find half of total area
00562       float halfArea = (area1 + area2 + area3 + area4 + area5 + area6 + area7) / 2;
00563       
00564       //determine which region split will be within and resulting horizontal midpoint
00565       
00566       //Within area 1
00567       if( area1 > halfArea )
00568       {
00569         clusters[cluster].beta = 0.5f + (float)sqrt(2*halfArea);       
00570       }
00571       //Within area 2
00572       else if( area1 + area2 > halfArea )
00573       {
00574         clusters[cluster].beta = ((halfArea - area1) / outM) + A;        
00575       }
00576       //Within area 3-4
00577       else if( area1 + area2 + area3 + area4 > halfArea )
00578       {
00579         float a = -0.5f;
00580         float b = 2.8f;
00581         float c = area1 + area2 + area3 - halfArea - B/2 - 2.625f;
00582         clusters[cluster].beta = (-b + (float)sqrt( b*b - 4*a*c )) / (2*a);
00583       }
00584       //Within area 5-6
00585       else if( area1 + area2 + area3 + area4 + area5 + area6 > halfArea )
00586       {
00587         float a = 1.0f/3;
00588         float b = -0.7f;
00589         float c = area1 + area2 + area3 + area4 - halfArea;
00590         clusters[cluster].beta = (-b + (float)sqrt( b*b - 4*a*c )) / (2*a);
00591       }
00592       //Within area 7
00593       else
00594       {
00595         clusters[cluster].beta = ((halfArea - (area1 + area2 + area3 + area4 + area5 + area6) ) / outL) + C;       
00596       }
00597     } //end case 1
00598     //---------------------------------------------------------------
00599     //Case 2
00600     else if ( outM < 0.5f && outL > outM )
00601     {
00602       //find area of 5 subregions
00603       float area1 = (outM*(A-0.5f)) / 2;
00604       float area2 = (D-A) * outM;
00605       float area3 = ((C-D) * (outL - outM)) / 2;
00606       float area4 = (C-D) * outM;
00607       float area5 = (3.0f - C) * outL;
00608       
00609       //find half of total area
00610       float halfArea = (area1 + area2 + area3 + area4 + area5) / 2;
00611       
00612       //determine which region split will be within and resulting horizontal midpoint
00613 
00614       //Within area 1
00615       if( area1 > halfArea )
00616       {
00617         clusters[cluster].beta = 0.5f + (float)sqrt(2*halfArea);
00618       }
00619       //Within area 2
00620       else if( area1 + area2 > halfArea )
00621       {
00622         clusters[cluster].beta = ((halfArea - area1) / outM) + A;
00623       }
00624       //Within area 3-4
00625       else if( area1 + area2 + area3 + area4 > halfArea )
00626       {
00627         float a = 1.0f/3.0f;
00628         float b = outM - 0.5f - D/3;
00629         float c = area1 + area2 - D*outM + D/2 - halfArea;
00630         clusters[cluster].beta = (-b + (float)sqrt( b*b - 4*a*c )) / (2*a);
00631       }
00632       //Within area5
00633       else
00634       {
00635         clusters[cluster].beta = ((halfArea - (area1 + area2 + area3 + area4) ) / outL) + C;
00636       }
00637     } //end case 2
00638     //---------------------------------------------------------------
00639     //Case 3
00640     else
00641     {
00642       //find area of 5 subregions
00643       float area1 = (outM*(A-0.5f)) / 2;
00644       float area2 = (B-A) * outM;
00645       float area3 = ((E-B) * (outM - outL)) / 2;
00646       float area4 = (E-B) * outL;
00647       float area5 = (3.0f - E) * outL;
00648       
00649       //find half of total area
00650       float halfArea = (area1 + area2 + area3 + area4 + area5) / 2;
00651       
00652       //determine which region split will be within and resulting horizontal midpoint
00653       
00654       //Within area 1
00655       if( area1 > halfArea )
00656       {
00657         clusters[cluster].beta = 0.5f + (float)sqrt(2*halfArea);
00658       }
00659       //Within area 2
00660       else if( area1 + area2 > halfArea )
00661       {
00662         clusters[cluster].beta = ((halfArea - area1) / outM) + A;
00663       }
00664       //Within area 3-4
00665       else if( area1 + area2 + area3 + area4 > halfArea )
00666       {
00667         float a = -0.5f;
00668         float b = E/2 + 2.5f/2; 
00669         float c = area3 - 2.5f*E/2;
00670         clusters[cluster].beta = (-b + (float)sqrt( b*b - 4*a*c )) / (2*a);
00671       }
00672       //Within area5
00673       else
00674       {
00675         clusters[cluster].beta = ((halfArea - (area1 + area2 + area3 + area4) ) / outL) + E;
00676       }
00677     } //end case 3
00678     //---------------------------------------------------------------
00679     
00680     //Compute edge threshold
00681     int lumJND = 255/50;
00682     clusters[cluster].edgeThreshold = clusters[cluster].mode + clusters[cluster].beta*lumJND;
00683 
00684   } //end for cluster
00685   
00686 }
00687 //==============================================
00688 void EdgeDetect::constructEdgeImage()
00689 {
00690   int x, y;
00691   QRgb* rgb;
00692   
00693   uchar* scanLine;
00694   for( y=0; y<image->height(); y++)
00695   {   
00696     scanLine = image->scanLine(y);
00697     for( x=0; x<image->width(); x++)
00698     {
00699       //initialize pixel to black 
00700       rgb = ((QRgb*)scanLine+x);
00701       *rgb = qRgb( 0, 0, 0 );
00702       
00703       //lookup ESF for this pixel
00704       float ESF = LUT[ GSLCmap[x + y*image->width()] ].ESF;
00705 
00706       //If ESF value for this pixel is 0 skip
00707       if( ESF == 0.0f ) continue;
00708 
00709       //lookup edge magnitude threshold
00710       float lum = lumMap[x + y*image->width()];
00711       float edgeMagThresh = -1.0f;
00712       int cluster;
00713       for(cluster=0; cluster<numClusters; cluster++)
00714       {
00715         if(lum >= clusters[cluster].minLuminance &&
00716            lum <= clusters[cluster].maxLuminance)
00717         {
00718           edgeMagThresh = clusters[cluster].edgeThreshold;
00719           break;
00720         }
00721       }
00722 
00723       //if cluster not found bail
00724       if( cluster >= numClusters )
00725       {
00726 //        cout << "Error! Could not find cluster pixel belonged to!\n";
00727         continue;
00728       }
00729       
00730       //if edge mag below thresh then skip
00731       if( edgeMagMap[x + y*image->width()] < edgeMagThresh ) continue;
00732         
00733       //ok, last checks implement NMS (non-maximum supression)
00734       int direction = LUT[ GSLCmap[x + y*image->width()] ].direction;
00735       int neighborIndex1 = -1;
00736       int neighborIndex2 = -1;
00737 
00738       if( direction == 0)
00739       {
00740         if( x > 0) 
00741           neighborIndex1 = x-1 + y*image->width();
00742         if( x < image->width() - 1 ) 
00743           neighborIndex2 = x+1 + y*image->width();
00744       }
00745       else if(direction == 1)
00746       {
00747         if( x > 0 && y < image->height() - 1 )
00748           neighborIndex1 = x-1 + (y+1)*image->width();
00749         if( x < image->width() - 1 && y > 0 )
00750           neighborIndex2 = x+1 + (y-1)*image->width();
00751       }
00752       else if(direction == 2)
00753       {
00754         if( y < image->height() - 1 ) 
00755           neighborIndex1 = x + (y+1)*image->width();
00756         if( y > 0) 
00757           neighborIndex2 = x + (y-1)*image->width();
00758       }
00759       else if(direction == 3)
00760       {
00761         if( x < image->width() - 1 && y < image->height() - 1 )
00762           neighborIndex1 = x+1 + (y+1)*image->width();
00763         if( x > 0 && y > 0 )
00764           neighborIndex2 = x-1 + (y-1)*image->width();
00765       }
00766 
00767       //neighbor 1 has higher confidence, skip!
00768       if( neighborIndex1 != -1 &&
00769           LUT[ GSLCmap[neighborIndex1] ].ESF * edgeMagMap[neighborIndex1] >
00770           ESF * edgeMagMap[x + y*image->width()] )
00771         continue;
00772       
00773       //neighbor 2 has higher confidence, skip!
00774       if( neighborIndex2 != -1 &&
00775           LUT[ GSLCmap[neighborIndex2] ].ESF * edgeMagMap[neighborIndex2] >
00776           ESF * edgeMagMap[x + y*image->width()] )
00777         continue;
00778       
00779       //All tests passed! Mark edge!
00780       *rgb = qRgb( 255, 255, 255 );
00781     } //x
00782   } //y
00783   
00784   //blur image - all of it
00785   blurImage( *image, 2.0f );
00786 
00787   //normalize image
00788   enhanceImageContrast( image );
00789   
00790 }
00791 //==============================================
00792 void EdgeDetect::deallocateObjects()
00793 {
00794   delete[] lumMap;      lumMap = NULL;
00795   delete[] edgeMagMap;  edgeMagMap = NULL;
00796   delete[] GSLCmap;     GSLCmap = NULL;
00797   delete[] clusters;    clusters = NULL;
00798 }
00799 //==============================================
00800 void EdgeDetect::constructGSLClut()
00801 {
00802   //----------------------
00803   //First fill entire table with 0 ESF's and invalid directions
00804   int i;
00805   for(i=0; i<256; i++)
00806   {
00807     LUT[i].ESF = 0.0f;
00808     LUT[i].direction = -1;
00809   }
00810   //----------------------
00811   //Next code in all pattern that are highly 
00812   //likely to be on edges as described in the paper
00813   //----------------------
00814   //Pattern (a)
00815 
00816   // ###
00817   // ##.
00818   // ...
00819   LUT[15].ESF = 0.179f;
00820   LUT[15].direction = 3;
00821 
00822   // ...
00823   // .##
00824   // ###
00825   LUT[240].ESF = 0.179f;
00826   LUT[240].direction = 3;
00827 
00828   // ###
00829   // .##
00830   // ...
00831   LUT[23].ESF = 0.179f;
00832   LUT[23].direction = 1;
00833   
00834   // ...
00835   // ##.
00836   // ###
00837   LUT[232].ESF = 0.179f;
00838   LUT[232].direction = 1;
00839    
00840   // ##.
00841   // ##.
00842   // #..
00843   LUT[43].ESF = 0.179f;
00844   LUT[43].direction = 3;
00845   
00846   // ..#
00847   // .##
00848   // .##
00849   LUT[212].ESF = 0.179f;
00850   LUT[212].direction = 3;
00851 
00852   // #..
00853   // ##.
00854   // ##.
00855   LUT[105].ESF = 0.179f;
00856   LUT[105].direction = 1;
00857 
00858   // .##
00859   // .##
00860   // ..#
00861   LUT[150].ESF = 0.179f;
00862   LUT[150].direction = 1;
00863   //----------------------
00864   //Pattern (b)
00865 
00866   // ###
00867   // ###
00868   // ...
00869   LUT[31].ESF = 0.137f;
00870   LUT[31].direction = 2;
00871 
00872   // ...
00873   // ###
00874   // ###
00875   LUT[248].ESF = 0.137f;
00876   LUT[248].direction = 2;
00877   
00878   // ##.
00879   // ##.
00880   // ##.
00881   LUT[107].ESF = 0.137f;
00882   LUT[107].direction = 0;
00883   
00884   // .##
00885   // .##
00886   // .##
00887   LUT[214].ESF = 0.137f;
00888   LUT[214].direction = 0;
00889   //----------------------
00890   //Pattern (c)
00891   
00892   // ###
00893   // .#.
00894   // ...
00895   LUT[7].ESF = 0.126f;
00896   LUT[7].direction = 2;
00897   
00898   // ...
00899   // .#.
00900   // ###
00901   LUT[224].ESF = 0.126f;
00902   LUT[224].direction = 2;
00903 
00904   // #..
00905   // ##.
00906   // #..
00907   LUT[41].ESF = 0.126f;
00908   LUT[41].direction = 0; 
00909   
00910   // ..#
00911   // .##
00912   // ..#
00913   LUT[148].ESF = 0.126f;
00914   LUT[148].direction = 0;
00915   //----------------------
00916   //Pattern (d)
00917   
00918   // ###
00919   // ##.
00920   // #..
00921   LUT[47].ESF = 0.10f;
00922   LUT[47].direction = 3;
00923    
00924   // ..#
00925   // .##
00926   // ###
00927   LUT[244].ESF = 0.10f;
00928   LUT[244].direction = 3;
00929 
00930   // ###
00931   // .##
00932   // ..#
00933   LUT[151].ESF = 0.10f;
00934   LUT[151].direction = 1;
00935 
00936   // #..
00937   // ##.
00938   // ###
00939   LUT[233].ESF = 0.10f;
00940   LUT[233].direction = 1;
00941   //----------------------
00942   //Pattern (e)
00943   
00944   // ##.
00945   // ##.
00946   // ...
00947   LUT[11].ESF = 0.10f;
00948   LUT[11].direction = 3;
00949   
00950   // ...
00951   // .##
00952   // .##
00953   LUT[208].ESF = 0.10f;
00954   LUT[208].direction = 3;
00955 
00956   // .##
00957   // .##
00958   // ...
00959   LUT[22].ESF = 0.10f;
00960   LUT[22].direction = 1;
00961   
00962   // ...
00963   // ##.
00964   // ##.
00965   LUT[104].ESF = 0.10f;
00966   LUT[104].direction = 1; 
00967   //----------------------    
00968 }
00969 //==============================================

Generated on Sat Apr 2 05:44:04 2005 for AlbumShaper by  doxygen 1.3.9.1