Could someone post some pseudo code on implementing CLAHE in Processing?

edited December 2015 in Questions about Code

I'm interested in creating a filter that implements Contrast Limited Adaptive Histogram Equalization(CLAHE). I think I understand how it works. You divide the image into a number of areas and perform as histogram equalization on each of those areas. Where I get lost is when you interpolate the histograms for each area.

I found this code for such a plugin that works with imageJ- and the filter works great. This is exactly what I want, but I get lost in it. Could someone give me a play by play of whats going on and some pseudo code of how to implemt this in processing. I do understand how do a histogram equalization.

The code I'm looking at is here : http://imagej.nih.gov/ij/plugins/clahe/CLAHE_.java

Tagged:

Answers

  • So, how did this become a library question? Its a question about how to implement some code, i.e , a programming question.

  • edited December 2015

    I got this working in processing, and actually I got a better , much faster version working in processing - so fast it can be used as a real time filter ( from 2480 ms for a pic to 186 ms for a pic). I didn't write it, but adapted it to processing . If anyone is interested, Ill post the code. Its a cool effect. And its FAST.

  • Here's the imageJ/fiji code: https://github.com/axtimwalde/mpicbg/tree/master/mpicbg/src/main/java/mpicbg/ij/clahe

    Here's the code adapted for processing.Run it and press any key to see the effect you can play with the values of slope, bins, and radius. Its real nice effect that can be used with other filters for painterly effects.

    PImage img =loadImage("http://ichef-1.bbci.co.uk/news/976/media/images/83351000/jpg/_83351965_explorer273lincolnshirewoldssouthpicturebynicholassilkstone.jpg","jpg"); ;
    
    void setup(){
         size(800, 800);
    
    
    }
    void draw(){
         background(200);
         image(img,0,0);
    
    }
    void keyPressed(){
              fastCLAHE(img);
    }
    public int roundPos( float a ){
              return ( int )( a + 0.5f );
    }
    public int[] createHistogram(int blockRadius,int bins,int blockXCenter,int blockYCenter, int[][]src, int w, int h ){
         int[] hist = new int[ bins + 1 ];          
         int xMin = Math.max( 0, blockXCenter - blockRadius );
         int yMin = Math.max( 0, blockYCenter - blockRadius );          
         int xMax = Math.min( w, blockXCenter + blockRadius + 1 );
         int yMax = Math.min( h, blockYCenter + blockRadius + 1 );         
         for ( int y = yMin; y < yMax; ++y ){
              for ( int x = xMin; x < xMax; ++x ){
                   ++hist[roundPos( src[x][y] / 255.0f * bins ) ];
              }
         }          
         return hist;
    }
    float[] createTransfer(int[] hist,int limit ){
         int[] cdfs = new int[ hist.length ];
         clipHistogram( hist, cdfs, limit );          
         int hMin = hist.length - 1;
         for ( int i = 0; i < hMin; ++i )
              if ( cdfs[ i ] != 0 ) hMin = i;          
         int cdf = 0;
         for ( int i = hMin; i < hist.length; ++i ){
              cdf += cdfs[ i ];
              cdfs[ i ] = cdf;
         }          
         int cdfMin = cdfs[ hMin ];
         int cdfMax = cdfs[ hist.length - 1 ];          
         float[] transfer = new float[ hist.length ];
         for ( int i = 0; i < transfer.length; ++i )
              transfer[ i ] = ( cdfs[ i ] - cdfMin ) / ( float )( cdfMax - cdfMin );          
         return transfer;
    }
    public void clipHistogram(int[] hist,int[] clippedHist, int limit ){
         System.arraycopy( hist, 0, clippedHist, 0, hist.length );
         int clippedEntries = 0, clippedEntriesBefore;
         do{
              clippedEntriesBefore = clippedEntries;
              clippedEntries = 0;
              for ( int i = 0; i < hist.length; ++i ){
                   final int d = clippedHist[ i ] - limit;
                   if ( d > 0 ){
                        clippedEntries += d;
                        clippedHist[ i ] = limit;
                   }
              }               
              int d = clippedEntries / ( hist.length );
              int m = clippedEntries % ( hist.length );
              for ( int i = 0; i < hist.length; ++i)
                   clippedHist[ i ] += d;               
              if ( m != 0 ){
                   final int s = ( hist.length - 1 ) / m;
                   for ( int i = s / 2; i < hist.length; i += s )
                        ++clippedHist[ i ];
              }
         }while ( clippedEntries != clippedEntriesBefore );
    }
    public void fastCLAHE(PImage myImage){
         println("clahe");
         int in = millis();
         int blockRadius = 63;
         int bins = 255;
         float slope = 6.0;
         myImage.loadPixels();
         int [][]src = new int[myImage.width][myImage.height];
         int [][]dst = new int[myImage.width][myImage.height];
         for(int y = 0,pix = 0; y < myImage.height;y++){
              for(int x = 0;x < myImage.width; x++,pix++){
                   float r = (myImage.pixels[pix]>> 16 & 0xFF) * 0.21;
                   float g = (myImage.pixels[pix]>>8 & 0xFF) * 0.72;
                   float b = (myImage.pixels[pix] & 0xFF) * 0.07;
                   src[x][y] = (int)(r +g +b +0.5);
              }
         }
    
         int[] hist;
         float[] tl;
         float[] tr;
         float[] bl;
         float[] br;
    
         int blockSize = 2 * blockRadius + 1;
         int limit = ( int )( slope * blockSize * blockSize / bins + 0.5f );
    
         /* div */
         int nc = myImage.width / blockSize;
         int nr = myImage.height / blockSize;
    
         /* % */
         int cm = myImage.width - nc * blockSize; //column modulus
         int[] cs;  //column sections
         switch ( cm )
         {
         case 0:
              cs = new int[ nc ];
              for ( int i = 0; i < nc; ++i )
                   cs[ i ] = i * blockSize + blockRadius + 1;
              break;
         case 1:
              cs = new int[ nc + 1 ];
              for ( int i = 0; i < nc; ++i )
                   cs[ i ] = i * blockSize + blockRadius + 1;
              cs[ nc ] = myImage.width - blockRadius - 1;
              break;
         default:
              cs = new int[ nc + 2 ];
              cs[ 0 ] = blockRadius + 1;
              for ( int i = 0; i < nc; ++i )
                   cs[ i + 1 ] = i * blockSize + blockRadius + 1 + cm / 2;
              cs[ nc + 1 ] = myImage.width - blockRadius - 1;
         }
    
         int rm = myImage.height - nr * blockSize; // row modulus
         int[] rs;  // row sections
         switch ( rm )
         {
         case 0:
              rs = new int[ nr ];
              for ( int i = 0; i < nr; ++i )
                   rs[ i ] = i * blockSize + blockRadius + 1;
              break;
         case 1:
              rs = new int[ nr + 1 ];
              for ( int i = 0; i < nr; ++i )
                   rs[ i ] = i * blockSize + blockRadius + 1;
              rs[ nr ] = myImage.height - blockRadius - 1;
              break;
         default:
              rs = new int[ nr + 2 ];
              rs[ 0 ] = blockRadius + 1;
              for ( int i = 0; i < nr; ++i )
                   rs[ i + 1 ] = i * blockSize + blockRadius + 1 + rm / 2;
              rs[ nr + 1 ] = myImage.height - blockRadius - 1;
         }
    
         for ( int r = 0; r <= rs.length; ++r )
         {
              int r0 = Math.max( 0, r - 1 );
              int r1 = Math.min( rs.length - 1, r );
              int dr = rs[ r1 ] - rs[ r0 ];
    
              hist = createHistogram( blockRadius, bins, cs[ 0 ], rs[ r0 ], src, myImage.width, myImage.height );
              tr = createTransfer( hist, limit );
              if ( r0 == r1 )
                   br = tr;
              else
              {
                   hist = createHistogram( blockRadius, bins, cs[ 0 ], rs[ r1 ], src,myImage.width,myImage.height );
                   br = createTransfer( hist, limit );
              }
    
              int yMin = ( r == 0 ? 0 : rs[ r0 ] );
              int yMax = ( r < rs.length ? rs[ r1 ] : myImage.height - 1 );
    
              for ( int c = 0; c <= cs.length; ++c )
              {
                   int c0 = Math.max( 0, c - 1 );
                   int c1 = Math.min( cs.length - 1, c );
                   int dc = cs[ c1 ] - cs[ c0 ];
    
                   tl = tr;
                   bl = br;
    
                   if ( c0 != c1 )
                   {
                        hist = createHistogram( blockRadius, bins, cs[ c1 ], rs[ r0 ], src,myImage.width,myImage.height );
                        tr = createTransfer( hist, limit );
                        if ( r0 == r1 )
                             br = tr;
                        else
                        {
                             hist = createHistogram( blockRadius, bins, cs[ c1 ], rs[ r1 ], src,myImage.width,myImage.height );
                             br = createTransfer( hist, limit );
                        }
                   }
    
                   int xMin = ( c == 0 ? 0 : cs[ c0 ] );
                   int xMax = ( c < cs.length ? cs[ c1 ] : myImage.width - 1 );
    
                   for ( int y = yMin; y < yMax; ++y )
                   {
                        int o = y * myImage.width;
                        float wy = ( float )( rs[ r1 ] - y ) / dr;
    
                        for ( int x = xMin; x < xMax; ++x )
                        {
                             float wx = ( float )( cs[ c1 ] - x ) / dc;
                             int v = roundPos( src[x][y] / 255.0f * bins );
    
                             float t00 = tl[ v ];
                             float t01 = tr[ v ];
                             float t10 = bl[ v ];
                             float t11 = br[ v ];                              
                             float t0, t1;
                             if ( c0 == c1 )
                             {
                                  t0 = t00;
                                  t1 = t10;
                             }
                             else
                             {
                                  t0 = wx * t00 + ( 1.0f - wx ) * t01;
                                  t1 = wx * t10 + ( 1.0f - wx ) * t11;
                             }
    
                             float t;
                             if ( r0 == r1 )
                                  t = t0;
                             else
                                  t = wy * t0 + ( 1.0f - wy ) * t1;
    
                             dst[x][y] = Math.max( 0, Math.min( 255, roundPos( t * 255.0f ) ) ) ;
                        }
                   }
    
              }
    
         }
         for( int y = 0; y < myImage.height; y++){
              int loc = y * myImage.width;
              for ( int x = 0; x < width; ++x,loc++ ){
                   float a;
                   int argb = myImage.pixels[loc];
                   int test = src[x][y];
                   if(test == 0)
                        a = 0.0f;
                   else a = ( float )dst[x][y] / src[x][y];               
                   int r = Math.max( 0, Math.min( 255,(int)(a * ( ( argb >> 16 ) & 0xff )+0.5) ) );  
                   int g = Math.max( 0, Math.min( 255, (int)(a * ( ( argb >> 8 ) & 0xff )+0.5) ) );
                   int b = Math.max( 0, Math.min( 255,  (int)(a * ( argb & 0xff ) +0.5) ) );
                   myImage.pixels[loc] = 0xFF000000| ( r << 16 ) | ( g << 8 ) | b ;
              }
    
         }
         myImage.updatePixels();
         int out = millis();
         println(out-in);
    
    
    
    }
    
Sign In or Register to comment.