ProgramGLSL.cpp
Go to the documentation of this file.
00001 
00002 //      File:           ProgramGLSL.cpp
00003 //      Author:         Changchang Wu
00004 //      Description : GLSL related classes
00005 //              class ProgramGLSL               A simple wrapper of GLSL programs
00006 //              class ShaderBagGLSL             GLSL shaders for SIFT
00007 //              class FilterGLSL                GLSL gaussian filters for SIFT
00008 //
00009 //      Copyright (c) 2007 University of North Carolina at Chapel Hill
00010 //      All Rights Reserved
00011 //
00012 //      Permission to use, copy, modify and distribute this software and its
00013 //      documentation for educational, research and non-profit purposes, without
00014 //      fee, and without a written agreement is hereby granted, provided that the
00015 //      above copyright notice and the following paragraph appear in all copies.
00016 //      
00017 //      The University of North Carolina at Chapel Hill make no representations
00018 //      about the suitability of this software for any purpose. It is provided
00019 //      'as is' without express or implied warranty. 
00020 //
00021 //      Please send BUG REPORTS to ccwu@cs.unc.edu
00022 //
00024 
00025   
00026 #include "GL/glew.h"
00027 #include <string.h>
00028 #include <stdio.h>
00029 #include <iomanip>
00030 #include <iostream>
00031 #include <strstream>
00032 #include <vector>
00033 #include <algorithm>
00034 #include <math.h>
00035 using namespace std;
00036 
00037 #include "GlobalUtil.h"
00038 #include "ProgramGLSL.h"
00039 #include "GLTexImage.h"
00040 #include "ShaderMan.h"
00041 #include "SiftGPU.h"
00042 
00043 ProgramGLSL::ShaderObject::ShaderObject(int shadertype, const char * source, int filesource)
00044 {
00045 
00046 
00047         _type = shadertype; 
00048         _compiled = 0;
00049 
00050 
00051         _shaderID = glCreateShader(shadertype);
00052         if(_shaderID == 0) return;
00053         
00054         if(source)
00055         {
00056 
00057                 GLint                           code_length;
00058                 if(filesource ==0)
00059                 {
00060                         const char* code  = source;
00061                         code_length = strlen(code);
00062                         glShaderSource(_shaderID, 1, (const char **) &code, &code_length);
00063                 }else
00064                 {
00065                         char * code;
00066                         if((code_length= ReadShaderFile(source, code)) ==0) return;
00067                         glShaderSource(_shaderID, 1, (const char **) &code, &code_length);
00068                         delete code;
00069                 }
00070 
00071                 glCompileShader(_shaderID);
00072 
00073                 CheckCompileLog();
00074 
00075                 if(!_compiled)          std::cout << source;
00076         }
00077 
00078 
00079 
00080 
00081 }
00082 
00083 int ProgramGLSL::ShaderObject::ReadShaderFile(const char *sourcefile,  char*& code )
00084 {
00085         code = NULL;
00086         FILE * file;
00087         int    len=0;
00088 
00089         if(sourcefile == NULL) return 0;
00090 
00091         file = fopen(sourcefile,"rt");
00092         if(file == NULL) return 0;
00093 
00094         
00095         fseek(file, 0, SEEK_END);
00096         len = ftell(file);
00097         rewind(file);
00098         if(len >1)
00099         {
00100                 code = new  char[len+1];
00101                 fread(code, sizeof( char), len, file);
00102                 code[len] = 0;
00103         }else
00104         {
00105                 len = 0;
00106         }
00107 
00108         fclose(file);
00109 
00110         return len;
00111         
00112 }
00113 
00114 void ProgramGLSL::ShaderObject::CheckCompileLog()
00115 {
00116 
00117         GLint status;
00118         glGetShaderiv(_shaderID, GL_COMPILE_STATUS, &status);
00119         _compiled = (status ==GL_TRUE);
00120 
00121         if(_compiled == 0)      PrintCompileLog(std::cout);
00122 
00123 
00124 }
00125 
00126 ProgramGLSL::ShaderObject::~ShaderObject()
00127 {
00128         if(_shaderID)   glDeleteShader(_shaderID);
00129 
00130 }
00131 
00132 int ProgramGLSL::ShaderObject::IsValidFragmentShader()
00133 {
00134         return _type == GL_FRAGMENT_SHADER && _shaderID && _compiled;
00135 }
00136 
00137 int  ProgramGLSL::ShaderObject::IsValidVertexShader()
00138 {
00139         return _type == GL_VERTEX_SHADER && _shaderID && _compiled;
00140 }
00141 
00142 
00143 void ProgramGLSL::ShaderObject::PrintCompileLog(ostream&os)
00144 {
00145         GLint len = 0;  
00146 
00147         glGetShaderiv(_shaderID, GL_INFO_LOG_LENGTH , &len);
00148         if(len <=1) return;
00149         
00150         char * compileLog = new char[len+1];
00151         if(compileLog == NULL) return;
00152 
00153         glGetShaderInfoLog(_shaderID, len, &len, compileLog);
00154         
00155 
00156         os<<"Compile Log\n"<<compileLog<<"\n";
00157 
00158         delete compileLog;
00159 }
00160 
00161 
00162 ProgramGLSL::ProgramGLSL()
00163 {
00164         _linked = 0;
00165         _programID = glCreateProgram();
00166 }
00167 ProgramGLSL::~ProgramGLSL()
00168 {
00169         if(_programID)glDeleteProgram(_programID);
00170 }
00171 void ProgramGLSL::AttachShaderObject(ShaderObject &shader)
00172 {
00173         if(_programID  && shader.IsValidShaderObject()) 
00174                 glAttachShader(_programID, shader.GetShaderID());
00175 }
00176 void ProgramGLSL::DetachShaderObject(ShaderObject &shader)
00177 {
00178         if(_programID  && shader.IsValidShaderObject()) 
00179                 glDetachShader(_programID, shader.GetShaderID());
00180 }
00181 int ProgramGLSL::LinkProgram()
00182 {
00183         _linked = 0;
00184 
00185         if(_programID==0) return 0;
00186 
00187         glLinkProgram(_programID);
00188 
00189         CheckLinkLog();
00190 
00191 //      GlobalUtil::StartTimer("100 link test");
00192 //      for(int i = 0; i<100; i++) glLinkProgram(_programID);
00193 //      GlobalUtil::StopTimer();
00194 
00195         return _linked;
00196 }
00197 
00198 void ProgramGLSL::CheckLinkLog()
00199 {
00200         GLint status;
00201         glGetProgramiv(_programID, GL_LINK_STATUS, &status);
00202 
00203         _linked = (status == GL_TRUE);
00204 
00205 }
00206 
00207 
00208 int ProgramGLSL::ValidateProgram()
00209 {
00210         if(_programID && _linked)
00211         {
00213 //              glValidateProgram(_programID);
00214 //              glGetProgramiv(_programID, GL_VALIDATE_STATUS, &status);
00215 //              return status == GL_TRUE;
00216                 return 1;
00217         }
00218         else
00219                 return 0;
00220 }
00221 
00222 void ProgramGLSL::PrintLinkLog(std::ostream &os)
00223 {
00224         GLint len = 0;  
00225 
00226         glGetProgramiv(_programID, GL_INFO_LOG_LENGTH , &len);
00227         if(len <=1) return;
00228         
00229         char* linkLog = new char[len+1];
00230         if(linkLog == NULL) return;
00231 
00232         glGetProgramInfoLog(_programID, len, &len, linkLog);
00233         
00234         linkLog[len] = 0;
00235 
00236         if(strstr(linkLog, "failed"))
00237         {
00238                 os<<linkLog + (linkLog[0] == ' '? 1:0)<<"\n";
00239                 _linked = 0;
00240         }
00241 
00242         delete linkLog;
00243 }
00244 
00245 int ProgramGLSL::UseProgram()
00246 {
00247         if(ValidateProgram())
00248         {
00249                 glUseProgram(_programID);
00250                 return true;
00251         }
00252         else
00253         {
00254                 return false;
00255         }
00256 }
00257 
00258 
00259 ProgramGLSL::ProgramGLSL(const char *frag_source)
00260 {
00261         _linked = 0;
00262         _programID = glCreateProgram();
00263         ShaderObject shader(GL_FRAGMENT_SHADER, frag_source);
00264 
00265         if(shader.IsValidFragmentShader())
00266         {
00267                 AttachShaderObject(shader);
00268                 LinkProgram();
00269 
00270                 if(!_linked)
00271                 {
00272                         //shader.PrintCompileLog(std::cout);
00273                         PrintLinkLog(std::cout);
00274                 }
00275         }else
00276         {
00277                 _linked = 0;
00278         }
00279         
00280 }
00281 
00282 /*
00283 ProgramGLSL::ProgramGLSL(char*frag_source, char * vert_source)
00284 {
00285         _used = 0;
00286         _linked = 0;
00287         _programID = glCreateProgram();
00288         ShaderObject shader(GL_FRAGMENT_SHADER, frag_source);
00289         ShaderObject vertex_shader(GL_VERTEX_SHADER, vert_source);
00290         AttachShaderObject(shader);
00291         AttachShaderObject(vertex_shader);
00292         LinkProgram();
00293         if(!_linked)
00294         {
00295                 shader.PrintCompileLog(std::cout);
00296                 vertex_shader.PrintCompileLog(std::cout);
00297                 PrintLinkLog(std::cout);
00298                 std::cout<<vert_source;
00299                 std::cout<<frag_source;
00300         }
00301 
00302 }
00303 */
00304 
00305 
00306 
00307 void ProgramGLSL::ReLink()
00308 {
00309         glLinkProgram(_programID);
00310 }
00311 
00312 int ProgramGLSL::IsNative()
00313 {
00314         return _linked;
00315 }
00316 
00317 FilterGLSL::FilterGLSL(float sigma) 
00318 {
00319         //pixel inside 3*sigma box
00320         int sz = int( ceil( GlobalUtil::_FilterWidthFactor * sigma -0.5) ) ;//
00321         int width = 2*sz + 1;
00322 
00323         //filter size truncation
00324         if(GlobalUtil::_MaxFilterWidth >0 && width > GlobalUtil::_MaxFilterWidth)
00325         {
00326                 std::cout<<"Filter size truncated from "<<width<<" to "<<GlobalUtil::_MaxFilterWidth<<endl;
00327                 sz = GlobalUtil::_MaxFilterWidth>>1;
00328                 width = 2 * sz + 1;
00329         }
00330 
00331         int i;
00332         float * kernel = new float[width];
00333         float   rv = 1.0f/(sigma*sigma);
00334         float   v, ksum =0; 
00335 
00336         // pre-compute filter
00337         for( i = -sz ; i <= sz ; ++i) 
00338         {
00339                 kernel[i+sz] =  v = exp(-0.5f * i * i *rv) ;
00340                 ksum += v;
00341         }
00342 
00343         //normalize the kernel
00344         rv = 1.0f / ksum;
00345         for(i = 0; i< width ;i++) kernel[i]*=rv;
00346         //
00347         
00348     MakeFilterProgram(kernel, width);
00349 
00350         _size = sz;
00351 
00352         delete[] kernel; 
00353     if(GlobalUtil::_verbose && GlobalUtil::_timingL) std::cout<<"Filter: sigma = "<<sigma<<", size = "<<width<<"x"<<width<<endl;
00354 }
00355 
00356 
00357 void FilterGLSL::MakeFilterProgram(float kernel[], int width)
00358 {
00359         if(GlobalUtil::_usePackedTex)
00360         {
00361                 s_shader_h = CreateFilterHPK(kernel, width);
00362                 s_shader_v = CreateFilterVPK(kernel, width);
00363         }else
00364         {
00365                 s_shader_h = CreateFilterH(kernel, width);
00366                 s_shader_v = CreateFilterV(kernel, width);
00367         }
00368 }
00369 
00370 ProgramGPU* FilterGLSL::CreateFilterH(float kernel[], int width)
00371 {
00372 
00373         char buffer[10240];
00374         ostrstream out(buffer, 10240);
00375         out<<setprecision(8);
00376 
00377         out<<  "uniform sampler2DRect tex;";
00378         out<< "\nvoid main(void){ float intensity = 0.0 ;  vec2 pos;\n";
00379 
00380     int half_width = width / 2; 
00381         for(int i = 0; i< width; i++)
00382         {
00383                 if(i == half_width)
00384                 {
00385 
00386                         out<<"float or = texture2DRect(tex, gl_TexCoord[0].st).r;\n";
00387                         out<<"intensity+= or * "<<kernel[i]<<";\n";
00388                 }else
00389                 {
00390                         out<<"pos = gl_TexCoord[0].st + vec2(float("<< (i - half_width) <<") , 0);\n";
00391                         out<<"intensity+= "<<kernel[i]<<"*texture2DRect(tex, pos).r;\n";
00392                 }
00393         }
00394 
00395         //copy original data to red channel
00396         out<<"gl_FragColor.r = or;\n"; 
00397         out<<"gl_FragColor.b  = intensity;}\n"<<'\0';
00398 
00399         return new ProgramGLSL( buffer);
00400 }
00401 
00402 
00403 ProgramGPU* FilterGLSL::CreateFilterV(float kernel[], int height)
00404 {
00405 
00406         char buffer[10240];
00407         ostrstream out(buffer, 10240);
00408         out<<setprecision(8);
00409 
00410         out<<  "uniform sampler2DRect tex;";
00411         out<< "\nvoid main(void){ float intensity = 0.0;vec2 pos; \n";
00412     int half_height = height / 2; 
00413         for(int i = 0; i< height; i++)
00414         {
00415 
00416                 if(i == half_height)
00417                 {
00418                         out<<"vec2 orb = texture2DRect(tex, gl_TexCoord[0].st).rb;\n";
00419                         out<<"intensity+= orb.y * "<<kernel[i]<<";\n";
00420 
00421                 }else
00422                 {
00423                         out<<"pos = gl_TexCoord[0].st + vec2(0, float("<<(i - half_height) <<") );\n";
00424                         out<<"intensity+= texture2DRect(tex, pos).b * "<<kernel[i]<<";\n";
00425                 }
00426                 
00427         }
00428 
00429         out<<"gl_FragColor.b = orb.y;\n";
00430         out<<"gl_FragColor.g = intensity - orb.x;\n"; // difference of gaussian..
00431         out<<"gl_FragColor.r = intensity;}\n"<<'\0';
00432         
00433 //      std::cout<<buffer<<endl;
00434         return new ProgramGLSL( buffer);
00435 }
00436 
00437 
00438 
00439 ProgramGPU* FilterGLSL::CreateFilterHPK(float kernel[], int width)
00440 {
00441         //both h and v are packed...
00442         int i, j , xw, xwn;
00443 
00444         int halfwidth  = width >>1;
00445         float * pf = kernel + halfwidth;
00446         int nhpixel = (halfwidth+1)>>1; //how many neighbour pixels need to be looked up
00447         int npixel  = (nhpixel<<1)+1;//
00448         char buffer[10240];
00449         float weight[3];
00450         ostrstream out(buffer, 10240);
00451         out<<setprecision(8);
00452 
00453         out<<  "uniform sampler2DRect tex;";
00454         out<< "\nvoid main(void){ vec4 result = vec4(0, 0, 0, 0);\n";
00456         out<<"vec4 pc; vec2 coord; \n";
00457         for( i = 0 ; i < npixel ; i++)
00458         {
00459                 out<<"coord = gl_TexCoord[0].xy + vec2(float("<<i-nhpixel<<"),0);\n";
00460                 out<<"pc=texture2DRect(tex, coord);\n";
00461                 if(GlobalUtil::_PreciseBorder)          out<<"if(coord.x < 0) pc = pc.rrbb;\n";
00462                 //for each sub-pixel j  in center, the weight of sub-pixel k 
00463                 xw = (i - nhpixel)*2;
00464                 for( j = 0; j < 3; j++)
00465                 {
00466                         xwn = xw  + j  -1;
00467                         weight[j] = xwn < -halfwidth || xwn > halfwidth? 0 : pf[xwn];
00468                 }
00469                 if(weight[1] == 0.0)
00470                 {
00471                         out<<"result += vec4("<<weight[2]<<","<<weight[0]<<","<<weight[2]<<","<<weight[0]<<")*pc.grab;\n";
00472                 }
00473                 else
00474                 {
00475                         out<<"result += vec4("<<weight[1]<<", "<<weight[0]<<", "<<weight[1]<<", "<<weight[0]<<")*pc.rrbb;\n";
00476                         out<<"result += vec4("<<weight[2]<<", "<<weight[1]<<", "<<weight[2]<<", "<<weight[1]<<")*pc.ggaa;\n";
00477                 }       
00478         
00479         }
00480         out<<"gl_FragColor = result;}\n"<<'\0';
00481         
00482         return new ProgramGLSL( buffer);
00483 
00484 
00485 }
00486 
00487 
00488 ProgramGPU* FilterGLSL::CreateFilterVPK(float kernel[], int height)
00489 {
00490 
00491         //both h and v are packed...
00492         int i, j, yw, ywn;
00493 
00494         int halfh  = height >>1;
00495         float * pf = kernel + halfh;
00496         int nhpixel = (halfh+1)>>1;     //how many neighbour pixels need to be looked up
00497         int npixel  = (nhpixel<<1)+1;//
00498         char buffer[10240];
00499         float weight[3];
00500         ostrstream out(buffer, 10240);
00501         out<<setprecision(8);
00502 
00503         out<<  "uniform sampler2DRect tex;";
00504         out<< "\nvoid main(void){ vec4 result = vec4(0, 0, 0, 0);\n";
00506         out<<"vec4 pc; vec2 coord;\n";
00507         for( i = 0 ; i < npixel ; i++)
00508         {
00509                 out<<"coord = gl_TexCoord[0].xy + vec2(0, float("<<i-nhpixel<<"));\n";
00510                 out<<"pc=texture2DRect(tex, coord);\n";
00511                 if(GlobalUtil::_PreciseBorder)  out<<"if(coord.y < 0) pc = pc.rgrg;\n";
00512 
00513                 //for each sub-pixel j  in center, the weight of sub-pixel k 
00514                 yw = (i - nhpixel)*2;
00515                 for( j = 0; j < 3; j++)
00516                 {
00517                         ywn = yw + j  -1;
00518                         weight[j] = ywn < -halfh || ywn > halfh? 0 : pf[ywn];
00519                 }
00520                 if(weight[1] == 0.0)
00521                 {
00522                         out<<"result += vec4("<<weight[2]<<","<<weight[2]<<","<<weight[0]<<","<<weight[0]<<")*pc.barg;\n";
00523                 }else
00524                 {
00525                         out<<"result += vec4("<<weight[1]<<","<<weight[1]<<","<<weight[0]<<","<<weight[0]<<")*pc.rgrg;\n";
00526                         out<<"result += vec4("<<weight[2]<<","<<weight[2]<<","<<weight[1]<<","<<weight[1]<<")*pc.baba;\n";
00527                 }
00528         }
00529         out<<"gl_FragColor = result;}\n"<<'\0';
00530 
00531         return new ProgramGLSL( buffer);
00532 }
00533 
00534 
00535 
00536 ShaderBag::ShaderBag()
00537 {
00538         s_debug = 0;
00539         s_orientation = 0;
00540         s_display_gaussian = 0;
00541         s_display_dog = 0;
00542         s_display_grad = 0;
00543         s_display_keys = 0;
00544         s_sampling = 0;
00545         s_grad_pass = 0;
00546         s_dog_pass = 0;
00547         s_keypoint = 0;
00548         s_genlist_init_tight = 0;
00549         s_genlist_init_ex = 0;
00550         s_genlist_histo = 0;
00551         s_genlist_start = 0;
00552         s_genlist_step = 0;
00553         s_genlist_end = 0;
00554         s_vertex_list = 0;
00555         s_descriptor_fp = 0;
00556         s_margin_copy = 0;
00558     f_gaussian_skip0 = NULL;
00559     f_gaussian_skip1 = NULL;
00560     f_gaussian_step = NULL;
00561     _gaussian_step_num = 0; 
00562 
00563 }
00564 
00565 ShaderBag::~ShaderBag()
00566 {
00567         if(s_debug)delete s_debug;
00568         if(s_orientation)delete s_orientation;
00569         if(s_display_gaussian)delete s_display_gaussian;
00570         if(s_display_dog)delete s_display_dog;
00571         if(s_display_grad)delete s_display_grad;
00572         if(s_display_keys)delete s_display_keys;
00573         if(s_sampling)delete s_sampling;
00574         if(s_grad_pass)delete s_grad_pass;
00575         if(s_dog_pass) delete s_dog_pass;
00576         if(s_keypoint)delete s_keypoint;
00577         if(s_genlist_init_tight)delete s_genlist_init_tight;
00578         if(s_genlist_init_ex)delete s_genlist_init_ex;
00579         if(s_genlist_histo)delete s_genlist_histo;
00580         if(s_genlist_start)delete s_genlist_start;
00581         if(s_genlist_step)delete s_genlist_step;
00582         if(s_genlist_end)delete s_genlist_end;
00583         if(s_vertex_list)delete s_vertex_list;
00584         if(s_descriptor_fp)delete s_descriptor_fp;
00585         if(s_margin_copy) delete s_margin_copy;
00586 
00588     if(f_gaussian_skip1) delete f_gaussian_skip1;
00589 
00590     for(unsigned int i = 0; i < f_gaussian_skip0_v.size(); i++)
00591     {
00592             if(f_gaussian_skip0_v[i]) delete f_gaussian_skip0_v[i];
00593     }
00594     if(f_gaussian_step && _gaussian_step_num > 0) 
00595     {
00596             for(int i = 0; i< _gaussian_step_num; i++)
00597             {
00598                     delete f_gaussian_step[i];
00599             }
00600             delete[] f_gaussian_step;
00601     }
00602 }
00603 
00604 
00605 void ShaderBag::SelectInitialSmoothingFilter(int octave_min, SiftParam&param)
00606 {
00607     float sigma = param.GetInitialSmoothSigma(octave_min);
00608     if(sigma == 0) 
00609     {
00610        f_gaussian_skip0 = NULL; 
00611     }else
00612     {
00613             for(unsigned int i = 0; i < f_gaussian_skip0_v.size(); i++)
00614             {
00615                     if(f_gaussian_skip0_v[i]->_id == octave_min)
00616                     {
00617                             f_gaussian_skip0 = f_gaussian_skip0_v[i];
00618                             return ;
00619                     }
00620             }
00621             FilterGLSL * filter = new FilterGLSL(sigma); 
00622             filter->_id = octave_min;
00623             f_gaussian_skip0_v.push_back(filter);
00624             f_gaussian_skip0 = filter; 
00625     }
00626 }
00627 
00628 void ShaderBag::CreateGaussianFilters(SiftParam&param)
00629 {
00630         if(param._sigma_skip0>0.0f) 
00631         {
00632         FilterGLSL * filter;
00633                 f_gaussian_skip0 = filter = new FilterGLSL(param._sigma_skip0);
00634                 filter->_id = GlobalUtil::_octave_min_default; 
00635                 f_gaussian_skip0_v.push_back(filter);
00636         }
00637         if(param._sigma_skip1>0.0f) 
00638         {
00639                 f_gaussian_skip1 = new FilterGLSL(param._sigma_skip1);
00640         }
00641 
00642         f_gaussian_step = new FilterProgram*[param._sigma_num];
00643         for(int i = 0; i< param._sigma_num; i++)
00644         {
00645                 f_gaussian_step[i] =  new FilterGLSL(param._sigma[i]);
00646         }
00647     _gaussian_step_num = param._sigma_num;
00648 }
00649 
00650 
00651 void ShaderBag::LoadDynamicShaders(SiftParam& param)
00652 {
00653     LoadKeypointShader(param._dog_threshold, param._edge_threshold);
00654     LoadGenListShader(param._dog_level_num, 0);
00655     CreateGaussianFilters(param);
00656 }
00657 
00658 
00659 void ShaderBagGLSL::LoadFixedShaders()
00660 {
00661 
00662 
00663         s_gray = new ProgramGLSL( 
00664                 "uniform sampler2DRect rgbTex; void main(void){\n"
00665                 "float intensity = dot(vec3(0.299, 0.587, 0.114), texture2DRect(rgbTex,gl_TexCoord[0].st ).rgb);\n"
00666                 "gl_FragColor = vec4(intensity, intensity, intensity, 1.0);}");
00667 
00668 
00669         s_debug = new ProgramGLSL( "void main(void){gl_FragColor.rg =  gl_TexCoord[0].st;}");
00670 
00671 
00672         s_sampling = new ProgramGLSL(
00673                 "uniform sampler2DRect tex; void main(void){gl_FragColor.rg= texture2DRect(tex, gl_TexCoord[0].st).rg;}");
00674 
00675         //
00676         s_grad_pass = new ProgramGLSL(
00677         "uniform sampler2DRect tex; void main ()\n"
00678         "{\n"
00679         "       vec4 v1, v2, gg;\n"
00680         "       vec4 cc  = texture2DRect(tex, gl_TexCoord[0].xy);\n"
00681         "       gg.x = texture2DRect(tex, gl_TexCoord[1].xy).r;\n"
00682         "       gg.y = texture2DRect(tex, gl_TexCoord[2].xy).r;\n"
00683         "       gg.z = texture2DRect(tex, gl_TexCoord[3].xy).r;\n"
00684         "       gg.w = texture2DRect(tex, gl_TexCoord[4].xy).r;\n"
00685         "       vec2 dxdy = (gg.yw - gg.xz); \n"
00686         "       float grad = 0.5*length(dxdy);\n"
00687         "       float theta = grad==0? 0.0: atan(dxdy.y, dxdy.x);\n"
00688         "       gl_FragData[0] = vec4(cc.rg, grad, theta);\n"
00689         "}\n\0");
00690 
00691         ProgramGLSL * program;
00692         s_margin_copy = program = new ProgramGLSL(
00693         "uniform sampler2DRect tex; uniform vec2 truncate;\n"
00694         "void main(){ gl_FragColor = texture2DRect(tex, min(gl_TexCoord[0].xy, truncate)); }");
00695 
00696         _param_margin_copy_truncate = glGetUniformLocation(*program, "truncate");
00697 
00698 
00699         GlobalUtil::_OrientationPack2 = 0;
00700         LoadOrientationShader();
00701 
00702         if(s_orientation == NULL)
00703         {
00704                 //Load a simplified version if the right version is not supported
00705                 s_orientation = program =  new ProgramGLSL(
00706                 "uniform sampler2DRect fTex; uniform sampler2DRect oTex;\n"
00707         "       uniform float size; void main(){\n"
00708         "       vec4 cc = texture2DRect(fTex, gl_TexCoord[0].st);\n"
00709         "       vec4 oo = texture2DRect(oTex, cc.rg);\n"
00710         "       gl_FragColor.rg = cc.rg;\n"
00711         "       gl_FragColor.b = oo.a;\n"
00712         "       gl_FragColor.a = size;}");  
00713 
00714                 _param_orientation_gtex = glGetUniformLocation(*program, "oTex");
00715                 _param_orientation_size = glGetUniformLocation(*program, "size");
00716 
00717                 GlobalUtil::_MaxOrientation = 0;
00718                 GlobalUtil::_FullSupported = 0;
00719                 std::cerr<<"Orientation simplified on this hardware"<<endl;
00720         }
00721 
00722 
00723         if(GlobalUtil::_DescriptorPPT) LoadDescriptorShader();
00724         if(s_descriptor_fp == NULL) 
00725         {
00726                 GlobalUtil::_DescriptorPPT = GlobalUtil::_FullSupported = 0; 
00727                 std::cerr<<"Descriptor ignored on this hardware"<<endl;
00728         }
00729 
00730         s_zero_pass = new ProgramGLSL("void main(){gl_FragColor = vec4(0.0);}");
00731 }
00732 
00733 
00734 void ShaderBagGLSL::LoadDisplayShaders()
00735 {
00736         s_copy_key = new ProgramGLSL(
00737                 "uniform sampler2DRect tex; void main(){\n"
00738         "gl_FragColor.rg= texture2DRect(tex, gl_TexCoord[0].st).rg; gl_FragColor.ba = vec2(0.0,1.0);    }");
00739         
00740 
00741         ProgramGLSL * program;
00742         s_vertex_list = program = new ProgramGLSL(
00743         "uniform vec4 sizes; uniform sampler2DRect tex;\n"
00744         "void main(void){\n"
00745         "float fwidth = sizes.y; float twidth = sizes.z; float rwidth = sizes.w; \n"
00746         "float index = 0.1*(fwidth*floor(gl_TexCoord[0].y) + gl_TexCoord[0].x);\n"
00747         "float px = mod(index, twidth);\n"
00748         "vec2 tpos= floor(vec2(px, index*rwidth))+0.5;\n"
00749         "vec4 cc = texture2DRect(tex, tpos );\n"
00750         "float size = 3.0f * cc.a; //sizes.x;// \n"
00751         "gl_FragColor.zw = vec2(0.0, 1.0);\n"
00752         "if(cc.x<=0 || cc.y <=0) {gl_FragColor.xy = cc.xy; }\n"
00753         "else {float type = fract(px);\n"
00754         "vec2 dxy; \n"
00755         "dxy.x = type < 0.1 ? 0.0 : ((type <0.5 || type > 0.9)? size : -size);\n"
00756         "dxy.y = type < 0.2 ? 0.0 : ((type < 0.3 || type > 0.7 )? -size :size); \n"
00757         "float s = sin(cc.b); float c = cos(cc.b); \n"
00758         "gl_FragColor.x = cc.x + c*dxy.x-s*dxy.y;\n"
00759         "gl_FragColor.y = cc.y + c*dxy.y+s*dxy.x;}\n}\n");
00760 
00761 
00762         _param_genvbo_size = glGetUniformLocation(*program, "sizes");
00763         s_display_gaussian =  new ProgramGLSL(
00764         "uniform sampler2DRect tex; void main(void){float r = texture2DRect(tex, gl_TexCoord[0].st).r;\n"
00765         "gl_FragColor = vec4(r, r, r, 1);}" );
00766 
00767         s_display_dog =  new ProgramGLSL(
00768         "uniform sampler2DRect tex; void main(void){float g = 0.5+(20.0*texture2DRect(tex, gl_TexCoord[0].st).g);\n"
00769         "gl_FragColor = vec4(g, g, g, 0.0);}" );
00770 
00771         s_display_grad = new ProgramGLSL(
00772                 "uniform sampler2DRect tex; void main(void){\n"
00773     "   vec4 cc = texture2DRect(tex, gl_TexCoord[0].st);gl_FragColor = vec4(5.0* cc.bbb, 1.0);}");
00774 
00775         s_display_keys= new ProgramGLSL(
00776                 "uniform sampler2DRect tex; void main(void){\n"
00777         "       vec4 cc = texture2DRect(tex, gl_TexCoord[0].st);\n"
00778         "       if(cc.r ==0.0) discard; gl_FragColor =  (cc.r==1.0? vec4(1.0, 0.0, 0,1.0):vec4(0.0,1.0,0.0,1.0));}");
00779 }
00780 
00781 void ShaderBagGLSL::LoadKeypointShader(float threshold, float edge_threshold)
00782 {
00783 
00784         char buffer[10240];
00785         float threshold0 = threshold* (GlobalUtil::_SubpixelLocalization?0.8f:1.0f);    
00786         float threshold1 = threshold;
00787         float threshold2 = (edge_threshold+1)*(edge_threshold+1)/edge_threshold;
00788         ostrstream out(buffer, 10240);
00789         streampos pos;
00790 
00791         //tex(X)(Y)
00792         //X: (CLR) (CENTER 0, LEFT -1, RIGHT +1)  
00793         //Y: (CDU) (CENTER 0, DOWN -1, UP    +1) 
00794 
00795         out <<  "#define THRESHOLD0 " << threshold0 << "\n"
00796                         "#define THRESHOLD1 " << threshold1 << "\n"
00797                         "#define THRESHOLD2 " << threshold2 << "\n";
00798 
00799         out<<
00800         "uniform sampler2DRect tex, texU, texD; void main ()\n"
00801         "{\n"
00802         "       vec4 v1, v2, gg, temp;\n"
00803         "       vec2 TexRU = vec2(gl_TexCoord[2].x, gl_TexCoord[4].y); \n"
00804         "       vec4 cc  = texture2DRect(tex, gl_TexCoord[0].xy);\n"
00805         "       temp =  texture2DRect(tex, gl_TexCoord[1].xy);\n"
00806         "       v1.x =  temp.g;                 gg.x = temp.r;\n"
00807         "       temp = texture2DRect(tex, gl_TexCoord[2].xy) ;\n"
00808         "       v1.y = temp.g;                  gg.y = temp.r;\n"
00809         "       temp = texture2DRect(tex, gl_TexCoord[3].xy) ;\n"
00810         "       v1.z = temp.g;                  gg.z = temp.r;\n"
00811         "       temp = texture2DRect(tex, gl_TexCoord[4].xy) ;\n"
00812         "       v1.w = temp.g;                  gg.w = temp.r;\n"
00813         "       v2.x = texture2DRect(tex, gl_TexCoord[5].xy).g;\n"
00814         "       v2.y = texture2DRect(tex, gl_TexCoord[6].xy).g;\n"
00815         "       v2.z = texture2DRect(tex, gl_TexCoord[7].xy).g;\n"
00816         "       v2.w = texture2DRect(tex, TexRU.xy).g;\n"
00817         "       vec2 dxdy = (gg.yw - gg.xz); \n"
00818         "       float grad = 0.5*length(dxdy);\n"
00819         "       float theta = grad==0? 0.0: atan(dxdy.y, dxdy.x);\n"
00820         "       gl_FragData[0] = vec4(cc.rg, grad, theta);\n"
00821 
00822         //test against 8 neighbours
00823         //use variable to identify type of extremum
00824         //1.0 for local maximum and 0.5 for minimum
00825         <<
00826         "       float dog = 0.0; \n"
00827         "       gl_FragData[1] = vec4(0, 0, 0, 0); \n"
00828         "       dog = cc.g > THRESHOLD0 && all(greaterThan(cc.gggg, max(v1, v2)))?1.0: 0.0;\n"
00829         "       dog = cc.g < -THRESHOLD0 && all(lessThan(cc.gggg, min(v1, v2)))?0.5: dog;\n"
00830         "       if(dog == 0.0) return;\n";
00831 
00832         pos = out.tellp();
00833         //do edge supression first.. 
00834         //vector v1 is < (-1, 0), (1, 0), (0,-1), (0, 1)>
00835         //vector v2 is < (-1,-1), (-1,1), (1,-1), (1, 1)>
00836 
00837         out<<
00838         "       float fxx, fyy, fxy; \n"
00839         "       vec4 D2 = v1.xyzw - cc.gggg;\n"
00840         "       vec2 D4 = v2.xw - v2.yz;\n"
00841         "       fxx = D2.x + D2.y;\n"
00842         "       fyy = D2.z + D2.w;\n"
00843         "       fxy = 0.25*(D4.x + D4.y);\n"
00844         "       float fxx_plus_fyy = fxx + fyy;\n"
00845         "       float score_up = fxx_plus_fyy*fxx_plus_fyy; \n"
00846         "       float score_down = (fxx*fyy - fxy*fxy);\n"
00847         "       if( score_down <= 0 || score_up > THRESHOLD2 * score_down)return;\n";
00848 
00849         //...
00850         out<<" \n"
00851         "       vec2 D5 = 0.5*(v1.yw-v1.xz); \n"
00852         "       float fx = D5.x, fy = D5.y ; \n"
00853         "       float fs, fss , fxs, fys ; \n"
00854         "       vec2 v3; vec4 v4, v5, v6;\n"
00855         //read 9 pixels of upper level
00856         <<
00857         "       v3.x = texture2DRect(texU, gl_TexCoord[0].xy).g;\n"
00858         "       v4.x = texture2DRect(texU, gl_TexCoord[1].xy).g;\n"
00859         "       v4.y = texture2DRect(texU, gl_TexCoord[2].xy).g;\n"
00860         "       v4.z = texture2DRect(texU, gl_TexCoord[3].xy).g;\n"
00861         "       v4.w = texture2DRect(texU, gl_TexCoord[4].xy).g;\n"
00862         "       v6.x = texture2DRect(texU, gl_TexCoord[5].xy).g;\n"
00863         "       v6.y = texture2DRect(texU, gl_TexCoord[6].xy).g;\n"
00864         "       v6.z = texture2DRect(texU, gl_TexCoord[7].xy).g;\n"
00865         "       v6.w = texture2DRect(texU, TexRU.xy).g;\n"
00866         //compare with 9 pixels of upper level
00867         //read and compare with 9 pixels of lower level
00868         //the maximum case
00869         <<
00870         "       if(dog == 1.0)\n"
00871         "       {\n"
00872         "               if(cc.g < v3.x || any(lessThan(cc.gggg, v4)) ||any(lessThan(cc.gggg, v6)))return; \n"
00873         "               v3.y = texture2DRect(texD, gl_TexCoord[0].xy).g;\n"
00874         "               v5.x = texture2DRect(texD, gl_TexCoord[1].xy).g;\n"
00875         "               v5.y = texture2DRect(texD, gl_TexCoord[2].xy).g;\n"
00876         "               v5.z = texture2DRect(texD, gl_TexCoord[3].xy).g;\n"
00877         "               v5.w = texture2DRect(texD, gl_TexCoord[4].xy).g;\n"
00878         "               v6.x = texture2DRect(texD, gl_TexCoord[5].xy).g;\n"
00879         "               v6.y = texture2DRect(texD, gl_TexCoord[6].xy).g;\n"
00880         "               v6.z = texture2DRect(texD, gl_TexCoord[7].xy).g;\n"
00881         "               v6.w = texture2DRect(texD, TexRU.xy).g;\n"
00882         "               if(cc.g < v3.y || any(lessThan(cc.gggg, v5)) ||any(lessThan(cc.gggg, v6)))return; \n"
00883         "       }\n"
00884         //the minimum case
00885         <<
00886         "       else{\n"
00887         "       if(cc.g > v3.x || any(greaterThan(cc.gggg, v4)) ||any(greaterThan(cc.gggg, v6)))return; \n"
00888         "               v3.y = texture2DRect(texD, gl_TexCoord[0].xy).g;\n"
00889         "               v5.x = texture2DRect(texD, gl_TexCoord[1].xy).g;\n"
00890         "               v5.y = texture2DRect(texD, gl_TexCoord[2].xy).g;\n"
00891         "               v5.z = texture2DRect(texD, gl_TexCoord[3].xy).g;\n"
00892         "               v5.w = texture2DRect(texD, gl_TexCoord[4].xy).g;\n"
00893         "               v6.x = texture2DRect(texD, gl_TexCoord[5].xy).g;\n"
00894         "               v6.y = texture2DRect(texD, gl_TexCoord[6].xy).g;\n"
00895         "               v6.z = texture2DRect(texD, gl_TexCoord[7].xy).g;\n"
00896         "               v6.w = texture2DRect(texD, TexRU.xy).g;\n"
00897         "               if(cc.g > v3.y || any(greaterThan(cc.gggg, v5)) ||any(greaterThan(cc.gggg, v6)))return; \n"
00898         "       }\n";
00899 
00900         if(GlobalUtil::_SubpixelLocalization)
00901 
00902         // sub-pixel localization FragData1 = vec4(dog, 0, 0, 0); return;
00903         out <<
00904         "       fs = 0.5*( v3.x - v3.y );  \n"
00905         "       fss = v3.x + v3.y - cc.g - cc.g;\n"
00906         "       fxs = 0.25 * ( v4.y + v5.x - v4.x - v5.y);\n"
00907         "       fys = 0.25 * ( v4.w + v5.z - v4.z - v5.w);\n"
00908         
00909         // 
00910         // let dog difference be quatratic function  of dx, dy, ds; 
00911         // df(dx, dy, ds) = fx * dx + fy*dy + fs * ds + 
00912         //                                + 0.5 * ( fxx * dx * dx + fyy * dy * dy + fss * ds * ds)
00913         //                                + (fxy * dx * dy + fxs * dx * ds + fys * dy * ds)
00914         // (fx, fy, fs, fxx, fyy, fss, fxy, fxs, fys are the derivatives)
00915         
00916         //the local extremum satisfies
00917         // df/dx = 0, df/dy = 0, df/dz = 0
00918         
00919         //that is 
00920         // |-fx|     | fxx fxy fxs |   |dx|
00921         // |-fy|  =  | fxy fyy fys | * |dy|
00922         // |-fs|     | fxs fys fss |   |ds|
00923         // need to solve dx, dy, ds
00924 
00925         // Use Gauss elimination to solve the linear system
00926     <<
00927         "       vec3 dxys = vec3(0.0);          \n"
00928         "       vec4 A0, A1, A2 ;                       \n"
00929         "       A0 = vec4(fxx, fxy, fxs, -fx);  \n"
00930         "       A1 = vec4(fxy, fyy, fys, -fy);  \n"
00931         "       A2 = vec4(fxs, fys, fss, -fs);  \n"
00932         "       vec3 x3 = abs(vec3(fxx, fxy, fxs));             \n"
00933         "       float maxa = max(max(x3.x, x3.y), x3.z);        \n"
00934         "       if(maxa >= 1e-10 ) {                                            \n"
00935         "               if(x3.y ==maxa )                                                        \n"
00936         "               {                                                                                       \n"
00937         "                       vec4 TEMP = A1; A1 = A0; A0 = TEMP;     \n"
00938         "               }else if( x3.z == maxa )                                        \n"
00939         "               {                                                                                       \n"
00940         "                       vec4 TEMP = A2; A2 = A0; A0 = TEMP;     \n"
00941         "               }                                                                                       \n"
00942         "               A0 /= A0.x;                                                                     \n"
00943         "               A1 -= A1.x * A0;                                                        \n"
00944         "               A2 -= A2.x * A0;                                                        \n"
00945         "               vec2 x2 = abs(vec2(A1.y, A2.y));                \n"
00946         "               if( x2.y > x2.x )                                                       \n"
00947         "               {                                                                                       \n"
00948         "                       vec3 TEMP = A2.yzw;                                     \n"
00949         "                       A2.yzw = A1.yzw;                                                \n"
00950         "                       A1.yzw = TEMP;                                                  \n"
00951         "                       x2.x = x2.y;                                                    \n"
00952         "               }                                                                                       \n"
00953         "               if(x2.x >= 1e-10) {                                             \n"
00954         "                       A1.yzw /= A1.y;                                                         \n"
00955         "                       A2.yzw -= A2.y * A1.yzw;                                        \n"
00956         "                       if(abs(A2.z) >= 1e-10) {                \n"
00957         // compute dx, dy, ds: 
00958         <<
00959         "                               \n"
00960         "                               dxys.z = A2.w /A2.z;                                \n"
00961         "                               dxys.y = A1.w - dxys.z*A1.z;                        \n"
00962         "                               dxys.x = A0.w - dxys.z*A0.z - dxys.y*A0.y;      \n"
00963 
00964         //one more threshold which I forgot in versions prior to 286
00965         <<
00966         "                               bool dog_test = (abs(cc.g + 0.5*dot(vec3(fx, fy, fs), dxys ))<= THRESHOLD1) ;\n"
00967         "                               if(dog_test || any(greaterThan(abs(dxys), vec3(1.0)))) dog = 0;\n"
00968         "                       }\n"
00969         "               }\n"
00970         "       }\n"
00971     //keep the point when the offset is less than 1
00972         <<
00973         "       gl_FragData[1] = vec4( dog, dxys); \n";
00974         else
00975 
00976         out<<
00977         "       gl_FragData[1] =  vec4( dog, 0, 0, 0) ; \n";
00978 
00979         out<<
00980         "}\n" <<'\0';
00981 
00982 
00983 
00984         ProgramGLSL * program = new ProgramGLSL(buffer); 
00985         if(program->IsNative())
00986         {
00987                 s_keypoint = program ;
00988                 //parameter
00989         }else
00990         {
00991                 delete program;
00992                 out.seekp(pos);
00993                 out << 
00994         "       gl_FragData[1] =  vec4(dog, 0, 0, 0) ;  \n"
00995         "}\n" <<'\0';
00996                 s_keypoint = program = new ProgramGLSL(buffer);
00997                 GlobalUtil::_SubpixelLocalization = 0;
00998                 std::cerr<<"Detection simplified on this hardware"<<endl;
00999         }
01000 
01001         _param_dog_texu = glGetUniformLocation(*program, "texU");
01002         _param_dog_texd = glGetUniformLocation(*program, "texD");
01003 }
01004 
01005 
01006 void ShaderBagGLSL::SetDogTexParam(int texU, int texD)
01007 {
01008         glUniform1i(_param_dog_texu, 1);
01009         glUniform1i(_param_dog_texd, 2);
01010 }
01011 
01012 void ShaderBagGLSL::SetGenListStepParam(int tex, int tex0)
01013 {
01014         glUniform1i(_param_genlist_step_tex0, 1);       
01015 }
01016 void ShaderBagGLSL::SetGenVBOParam( float width, float fwidth,  float size)
01017 {
01018         float sizes[4] = {size*3.0f, fwidth, width, 1.0f/width};
01019         glUniform4fv(_param_genvbo_size, 1, sizes);
01020 
01021 }
01022 
01023 
01024 
01025 void ShaderBagGLSL::UnloadProgram()
01026 {
01027         glUseProgram(0);
01028 } 
01029 
01030 
01031 
01032 void ShaderBagGLSL::LoadGenListShader(int ndoglev, int nlev)
01033 {
01034         ProgramGLSL * program;
01035 
01036         s_genlist_init_tight = new ProgramGLSL(
01037         "uniform sampler2DRect tex; void main (void){\n"
01038         "vec4 helper = vec4( texture2DRect(tex, gl_TexCoord[0].xy).r,  texture2DRect(tex, gl_TexCoord[1].xy).r,\n"
01039         "texture2DRect(tex, gl_TexCoord[2].xy).r, texture2DRect(tex, gl_TexCoord[3].xy).r);\n"
01040         "gl_FragColor = vec4(greaterThan(helper, vec4(0.0,0.0,0.0,0.0)));\n"
01041         "}");
01042 
01043         
01044         s_genlist_init_ex = program = new ProgramGLSL(
01045         "uniform sampler2DRect tex;uniform vec2 bbox;\n"
01046         "void main (void ){\n"
01047         "vec4 helper = vec4( texture2DRect(tex, gl_TexCoord[0].xy).r,  texture2DRect(tex, gl_TexCoord[1].xy).r,\n"
01048         "texture2DRect(tex, gl_TexCoord[2].xy).r, texture2DRect(tex, gl_TexCoord[3].xy).r);\n"
01049         "bvec4 helper2 = bvec4( \n"
01050         "all(lessThan(gl_TexCoord[0].xy , bbox)) && helper.x >0,\n"
01051         "all(lessThan(gl_TexCoord[1].xy , bbox)) && helper.y >0,\n"
01052         "all(lessThan(gl_TexCoord[2].xy , bbox)) && helper.z >0,\n"
01053         "all(lessThan(gl_TexCoord[3].xy , bbox)) && helper.w >0);\n"
01054         "gl_FragColor = vec4(helper2);\n"
01055         "}");
01056         _param_genlist_init_bbox = glGetUniformLocation( *program, "bbox");
01057 
01058 
01059         //reduction ...
01060         s_genlist_histo = new ProgramGLSL(
01061         "uniform sampler2DRect tex; void main (void){\n"
01062         "vec4 helper; vec4 helper2; \n"
01063         "helper = texture2DRect(tex, gl_TexCoord[0].xy); helper2.xy = helper.xy + helper.zw; \n"
01064         "helper = texture2DRect(tex, gl_TexCoord[1].xy); helper2.zw = helper.xy + helper.zw; \n"
01065         "gl_FragColor.rg = helper2.xz + helper2.yw;\n"
01066         "helper = texture2DRect(tex, gl_TexCoord[2].xy); helper2.xy = helper.xy + helper.zw; \n"
01067         "helper = texture2DRect(tex, gl_TexCoord[3].xy); helper2.zw = helper.xy + helper.zw; \n"
01068         "gl_FragColor.ba= helper2.xz+helper2.yw;\n"
01069         "}");
01070 
01071 
01072         //read of the first part, which generates tex coordinates 
01073         s_genlist_start= program =  LoadGenListStepShader(1, 1);
01074         _param_ftex_width= glGetUniformLocation(*program, "width");
01075         _param_genlist_start_tex0 = glGetUniformLocation(*program, "tex0");
01076         //stepping
01077         s_genlist_step = program = LoadGenListStepShader(0, 1);
01078         _param_genlist_step_tex= glGetUniformLocation(*program, "tex");
01079         _param_genlist_step_tex0= glGetUniformLocation(*program, "tex0");
01080 
01081 }
01082 
01083 void ShaderBagGLSL::SetMarginCopyParam(int xmax, int ymax)
01084 {
01085         float truncate[2] = {xmax - 0.5f , ymax - 0.5f};
01086         glUniform2fv(_param_margin_copy_truncate, 1, truncate);
01087 }
01088 
01089 void ShaderBagGLSL::SetGenListInitParam(int w, int h)
01090 {
01091         float bbox[2] = {w - 1.0f, h - 1.0f};
01092         glUniform2fv(_param_genlist_init_bbox, 1, bbox);
01093 }
01094 void ShaderBagGLSL::SetGenListStartParam(float width, int tex0)
01095 {
01096         glUniform1f(_param_ftex_width, width);
01097 }
01098 
01099 
01100 ProgramGLSL* ShaderBagGLSL::LoadGenListStepShader(int start, int step)
01101 {
01102         int i;
01103         char buffer[10240];
01104         // char chanels[5] = "rgba";
01105         ostrstream out(buffer, 10240);
01106 
01107         for(i = 0; i < step; i++) out<<"uniform sampler2DRect tex"<<i<<";\n";
01108         if(start)
01109         {
01110                 out<<"uniform float width;\n";
01111                 out<<"void main(void){\n";
01112                 out<<"float  index = floor(gl_TexCoord[0].y) * width + floor(gl_TexCoord[0].x);\n";
01113                 out<<"vec2 pos = vec2(0.5, 0.5);\n";
01114         }else
01115         {
01116                 out<<"uniform sampler2DRect tex;\n";
01117                 out<<"void main(void){\n";
01118                 out<<"vec4 tc = texture2DRect( tex, gl_TexCoord[0].xy);\n";
01119                 out<<"vec2 pos = tc.rg; float index = tc.b;\n";
01120         }
01121         out<<"vec2 sum;         vec4 cc;\n";
01122 
01123 
01124         if(step>0)
01125         {
01126                 out<<"vec2 cpos = vec2(-0.5, 0.5);\t vec2 opos;\n";
01127                 for(i = 0; i < step; i++)
01128                 {
01129 
01130                         out<<"cc = texture2DRect(tex"<<i<<", pos);\n";
01131                         out<<"sum.x = cc.r + cc.g; sum.y = sum.x + cc.b;  \n";
01132                         out<<"if (index <cc.r){ opos = cpos.xx;}\n";
01133                         out<<"else if(index < sum.x ) {opos = cpos.yx; index -= cc.r;}\n";
01134                         out<<"else if(index < sum.y ) {opos = cpos.xy; index -= sum.x;}\n";
01135                         out<<"else {opos = cpos.yy; index -= sum.y;}\n";
01136                         out<<"pos = (pos + pos + opos);\n";
01137                 }
01138         }
01139         out<<"gl_FragColor = vec4(pos, index, 1.0);\n";
01140         out<<"}\n"<<'\0';
01141         return new ProgramGLSL(buffer);
01142 }
01143 
01144 
01145 void ShaderBagGLSL::LoadOrientationShader()
01146 {
01147         char buffer[10240];
01148         ostrstream out(buffer,10240);
01149 
01150         if(GlobalUtil::_IsNvidia)
01151         {
01152         out <<  "#pragma optionNV(ifcvt none)\n"
01153                         "#pragma optionNV(unroll all)\n";
01154         }
01155 
01156         out<<"\n"
01157         "#define GAUSSIAN_WF "<<GlobalUtil::_OrientationGaussianFactor<<" \n"
01158         "#define SAMPLE_WF ("<<GlobalUtil::_OrientationWindowFactor<< " )\n"
01159         "#define ORIENTATION_THRESHOLD "<< GlobalUtil::_MulitiOrientationThreshold << "\n"
01160         "uniform sampler2DRect tex;                                     \n"
01161         "uniform sampler2DRect gradTex;                         \n"
01162         "uniform vec4 size;                                             \n"
01163         << (GlobalUtil::_SubpixelLocalization || GlobalUtil::_KeepExtremumSign? 
01164         "       uniform sampler2DRect texS;     \n" : " ")
01165         <<
01166         "void main()            \n"
01167         "{                                                                                                      \n"
01168         "       vec4 bins[10];                                                          \n"
01169         "       bins[0] = vec4(0.0);bins[1] = vec4(0.0);bins[2] = vec4(0.0);    \n"
01170         "       bins[3] = vec4(0.0);bins[4] = vec4(0.0);bins[5] = vec4(0.0);    \n"
01171         "       bins[6] = vec4(0.0);bins[7] = vec4(0.0);bins[8] = vec4(0.0);    \n"
01172         "       vec4 loc = texture2DRect(tex, gl_TexCoord[0].xy);       \n"
01173         "       vec2 pos = loc.xy;              \n"
01174         "       bool orientation_mode = (size.z != 0);                  \n"
01175         "       float sigma = orientation_mode? abs(size.z) : loc.w; \n";
01176         if(GlobalUtil::_SubpixelLocalization || GlobalUtil::_KeepExtremumSign)
01177         {
01178                 out<<
01179         "       if(orientation_mode){\n"
01180         "               vec4 offset = texture2DRect(texS, pos);\n"
01181         "               pos.xy = pos.xy + offset.yz; \n"
01182         "               sigma = sigma * pow(size.w, offset.w);\n"
01183         "               #if "<< GlobalUtil::_KeepExtremumSign << "\n"
01184         "                       if(offset.x < 0.6) sigma = -sigma; \n"
01185         "               #endif\n"
01186         "       }\n";
01187         }
01188         out<<
01189         "       //bool fixed_orientation = (size.z < 0);                \n"
01190         "       if(size.z < 0) {gl_FragData[0] = vec4(pos, 0, sigma); return;}"
01191         "       float gsigma = sigma * GAUSSIAN_WF;                             \n"
01192         "       vec2 win = abs(vec2(sigma)) * (SAMPLE_WF * GAUSSIAN_WF);        \n"
01193         "       vec2 dim = size.xy;                                                     \n"
01194         "       float dist_threshold = win.x*win.x+0.5;                 \n"
01195         "       float factor = -0.5/(gsigma*gsigma);                    \n"
01196         "       vec4 sz;        vec2 spos;                                              \n"
01197         "       //if(any(pos.xy <= 1)) discard;                                 \n"
01198         "       sz.xy = max( pos - win, vec2(1,1));                     \n"
01199         "       sz.zw = min( pos + win, dim-2);                         \n"
01200         "       sz = floor(sz)+0.5;";
01201         //loop to get the histogram
01202 
01203         out<<"\n"
01204         "       for(spos.y = sz.y; spos.y <= sz.w;      spos.y+=1.0)                            \n"
01205         "       {                                                                                                                               \n"
01206         "               for(spos.x = sz.x; spos.x <= sz.z;      spos.x+=1.0)                    \n"
01207         "               {                                                                                                                       \n"
01208         "                       vec2 offset = spos - pos;                                                               \n"
01209         "                       float sq_dist = dot(offset,offset);                                             \n"
01210         "                       if( sq_dist < dist_threshold){                                                  \n"
01211         "                               vec4 cc = texture2DRect(gradTex, spos);                         \n"
01212         "                               float grad = cc.b;      float theta = cc.a;                             \n"
01213         "                               float idx = floor(degrees(theta)*0.1);                          \n"
01214         "                               if(idx < 0 ) idx += 36;                                                                 \n"
01215         "                               float weight = grad*exp(sq_dist * factor);                              \n"
01216         "                               float vidx = fract(idx * 0.25) * 4.0;//mod(idx, 4.0) ;                                                  \n"
01217         "                               vec4 inc = weight*vec4(equal(vec4(vidx), vec4(0.0,1.0,2.0,3.0)));";
01218 
01219         if(GlobalUtil::_UseDynamicIndexing && GlobalUtil::_IsNvidia)
01220         {
01221                 //dynamic indexing may not be faster
01222                 out<<"\n"
01223         "                               int iidx = int((idx*0.25));     \n"
01224         "                               bins[iidx]+=inc;                                        \n"
01225         "                       }                                                                               \n"
01226         "               }                                                                                       \n"
01227         "       }";
01228 
01229         }else
01230         {
01231                 //nvfp40 still does not support dynamic array indexing
01232                 //unrolled binary search...
01233                 out<<"\n"
01234         "                               if(idx < 16)                                                    \n"
01235         "                               {                                                                               \n"
01236         "                                       if(idx < 8)                                                     \n"
01237         "                                       {                                                                       \n"
01238         "                                               if(idx < 4)     {       bins[0]+=inc;}  \n"
01239         "                                               else            {       bins[1]+=inc;}  \n"
01240         "                                       }else                                                           \n"
01241         "                                       {                                                                       \n"
01242         "                                               if(idx < 12){   bins[2]+=inc;}  \n"
01243         "                                               else            {       bins[3]+=inc;}  \n"
01244         "                                       }                                                                       \n"
01245         "                               }else if(idx < 32)                                              \n"
01246         "                               {                                                                               \n"
01247         "                                       if(idx < 24)                                            \n"
01248         "                                       {                                                                       \n"
01249         "                                               if(idx <20)     {       bins[4]+=inc;}  \n"
01250         "                                               else            {       bins[5]+=inc;}  \n"
01251         "                                       }else                                                           \n"
01252         "                                       {                                                                       \n"
01253         "                                               if(idx < 28){   bins[6]+=inc;}  \n"
01254         "                                               else            {       bins[7]+=inc;}  \n"
01255         "                                       }                                                                       \n"
01256         "                               }else                                           \n"
01257         "                               {                                                                               \n"
01258         "                                       bins[8]+=inc;                                           \n"
01259         "                               }                                                                               \n"
01260         "                       }                                                                               \n"
01261         "               }                                                                                       \n"
01262         "       }";
01263 
01264         }
01265 
01266         WriteOrientationCodeToStream(out);
01267 
01268         ProgramGLSL * program = new ProgramGLSL(buffer);
01269         if(program->IsNative())
01270         {
01271                 s_orientation = program ;
01272                 _param_orientation_gtex = glGetUniformLocation(*program, "gradTex");
01273                 _param_orientation_size = glGetUniformLocation(*program, "size");
01274                 _param_orientation_stex = glGetUniformLocation(*program, "texS");
01275         }else
01276         {
01277                 delete program;
01278         }
01279 }
01280 
01281 
01282 void ShaderBagGLSL::WriteOrientationCodeToStream(std::ostream& out)
01283 {
01284         //smooth histogram and find the largest
01285 /*
01286         smoothing kernel:        (1 3 6 7 6 3 1 )/27
01287         the same as 3 pass of (1 1 1)/3 averaging
01288         maybe better to use 4 pass on the vectors...
01289 */
01290 
01291 
01292         //the inner loop on different array numbers is always unrolled in fp40
01293 
01294         //bug fixed here:)
01295         out<<"\n"
01296         "       //mat3 m1 = mat3(1, 0, 0, 3, 1, 0, 6, 3, 1)/27.0;  \n"
01297         "       mat3 m1 = mat3(1, 3, 6, 0, 1, 3,0, 0, 1)/27.0;  \n"
01298         "       mat4 m2 = mat4(7, 6, 3, 1, 6, 7, 6, 3, 3, 6, 7, 6, 1, 3, 6, 7)/27.0;\n"
01299         "       #define FILTER_CODE(i) {                                                \\\n"
01300         "                       vec4 newb       =       (bins[i]* m2);                  \\\n"
01301         "                       newb.xyz        +=      ( prev.yzw * m1);               \\\n"
01302         "                       prev = bins[i];                                                 \\\n"
01303         "                       newb.wzy        +=      ( bins[i+1].zyx *m1);   \\\n"
01304         "                       bins[i] = newb;}\n"
01305         "       for (int j=0; j<2; j++)                                                         \n"
01306         "       {                                                                                               \n"
01307         "               vec4 prev  = bins[8];                                           \n"
01308         "               bins[9]          = bins[0];                                             \n";
01309 
01310         if(GlobalUtil::_IsNvidia)
01311         {
01312                 out<<
01313         "               for (int i=0; i<9; i++)                                                 \n"
01314         "               {                                                                                               \n"
01315         "                       FILTER_CODE(i);                                                         \n"
01316         "               }                                                                                               \n"
01317         "       }";
01318 
01319         }else
01320         {
01321                 //manually unroll the loop for ATI.
01322                 out << 
01323         "          FILTER_CODE(0);\n"
01324         "          FILTER_CODE(1);\n"
01325         "          FILTER_CODE(2);\n"
01326         "          FILTER_CODE(3);\n"
01327         "          FILTER_CODE(4);\n"
01328         "          FILTER_CODE(5);\n"
01329         "          FILTER_CODE(6);\n"
01330         "          FILTER_CODE(7);\n"
01331         "          FILTER_CODE(8);\n"
01332         "       }\n";
01333         }
01334         //find the maximum voting
01335         out<<"\n"
01336         "       vec4 maxh; vec2 maxh2;  \n"
01337         "       vec4 maxh4 = max(max(max(max(max(max(max(max(bins[0], bins[1]), bins[2]), \n"
01338         "                       bins[3]), bins[4]), bins[5]), bins[6]), bins[7]), bins[8]);\n"
01339         "       maxh2 = max(maxh4.xy, maxh4.zw); maxh = vec4(max(maxh2.x, maxh2.y));";
01340 
01341         char *testpeak_code;
01342         char *savepeak_code;
01343 
01344         //save two/three/four orientations with the largest votings?
01345 
01346         if(GlobalUtil::_MaxOrientation>1)
01347         {
01348                 out<<"\n"
01349                 "       vec4 Orientations = vec4(0, 0, 0, 0);                           \n"
01350                 "       vec4 weights = vec4(0,0,0,0);           ";      
01351                 
01352                 testpeak_code = "\\\n"
01353                 "       {test = greaterThan(bins[i], hh);";
01354 
01355                 //save the orientations in weight-decreasing order
01356                 if(GlobalUtil::_MaxOrientation ==2)
01357                 {
01358                 savepeak_code = "\\\n"
01359                 "                       if(weight <=weights.g){}\\\n"
01360                 "                       else if(weight >weights.r)\\\n"
01361                 "                       {weights.rg = vec2(weight, weights.r); Orientations.rg = vec2(th, Orientations.r);}\\\n"
01362                 "                       else {weights.g = weight; Orientations.g = th;}";
01363                 }else if(GlobalUtil::_MaxOrientation ==3)
01364                 {
01365                 savepeak_code = "\\\n"
01366                 "                       if(weight <=weights.b){}\\\n"
01367                 "                       else if(weight >weights.r)\\\n"
01368                 "                       {weights.rgb = vec3(weight, weights.rg); Orientations.rgb = vec3(th, Orientations.rg);}\\\n"
01369                 "                       else if(weight >weights.g)\\\n"
01370                 "                       {weights.gb = vec2(weight, weights.g); Orientations.gb = vec2(th, Orientations.g);}\\\n"
01371                 "                       else {weights.b = weight; Orientations.b = th;}";
01372                 }else
01373                 {
01374                 savepeak_code = "\\\n"
01375                 "                       if(weight <=weights.a){}\\\n"
01376                 "                       else if(weight >weights.r)\\\n"
01377                 "                       {weights = vec4(weight, weights.rgb); Orientations = vec4(th, Orientations.rgb);}\\\n"
01378                 "                       else if(weight >weights.g)\\\n"
01379                 "                       {weights.gba = vec3(weight, weights.gb); Orientations.gba = vec3(th, Orientations.gb);}\\\n"
01380                 "                       else if(weight >weights.b)\\\n"
01381                 "                       {weights.ba = vec2(weight, weights.b); Orientations.ba = vec2(th, Orientations.b);}\\\n"
01382                 "                       else {weights.a = weight; Orientations.a = th;}";
01383                 }
01384 
01385         }else
01386         {
01387                 out<<"\n"
01388                 "       float Orientation;                              ";
01389                 testpeak_code ="\\\n"
01390                 "       if(npeaks<=0.0){\\\n"
01391                 "       test = equal(bins[i], maxh)     ;";
01392                 savepeak_code="\\\n"
01393                 "                       npeaks++;       \\\n"
01394                 "                       Orientation = th;";
01395 
01396         }
01397         //find the peaks
01398         out <<"\n"
01399         "       #define FINDPEAK(i, k)" <<testpeak_code<<"\\\n"
01400         "       if( any ( test) )                                                       \\\n"
01401         "       {                                                                                       \\\n"
01402         "               if(test.r && bins[i].x > prevb && bins[i].x > bins[i].y )       \\\n"
01403         "               {                                                                                       \\\n"
01404         "                   float       di = -0.5 * (bins[i].y-prevb) / (bins[i].y+prevb-bins[i].x - bins[i].x) ; \\\n"
01405         "                   float       th = (k+di+0.5);        float weight = bins[i].x;"
01406                                 <<savepeak_code<<"\\\n"
01407         "               }\\\n"
01408         "               else if(test.g && all( greaterThan(bins[i].yy , bins[i].xz)) )  \\\n"
01409         "               {                                                                                       \\\n"
01410         "                   float       di = -0.5 * (bins[i].z-bins[i].x) / (bins[i].z+bins[i].x-bins[i].y- bins[i].y) ; \\\n"
01411         "                   float       th = (k+di+1.5);        float weight = bins[i].y;                               "
01412                                 <<savepeak_code<<"      \\\n"
01413         "               }\\\n"
01414         "               if(test.b && all( greaterThan( bins[i].zz , bins[i].yw)) )      \\\n"
01415         "               {                                                                                       \\\n"
01416         "                   float       di = -0.5 * (bins[i].w-bins[i].y) / (bins[i].w+bins[i].y-bins[i].z- bins[i].z) ; \\\n"
01417         "                   float       th = (k+di+2.5);        float weight = bins[i].z;                               "
01418                                 <<savepeak_code<<"      \\\n"
01419         "               }\\\n"
01420         "               else if(test.a && bins[i].w > bins[i].z && bins[i].w > bins[i+1].x )    \\\n"
01421         "               {                                                                                       \\\n"
01422         "                   float       di = -0.5 * (bins[i+1].x-bins[i].z) / (bins[i+1].x+bins[i].z-bins[i].w - bins[i].w) ; \\\n"
01423         "                   float       th = (k+di+3.5);        float weight = bins[i].w;                               "
01424                                 <<savepeak_code<<"      \\\n"
01425         "               }\\\n"
01426         "       }}\\\n"
01427         "       prevb = bins[i].w;";
01428         //the following loop will be unrolled anyway in fp40,
01429         //taking more than 1000 instrucsions..
01430         //....
01431         if(GlobalUtil::_IsNvidia)
01432         {
01433         out<<"\n"
01434         "       vec4 hh = maxh * ORIENTATION_THRESHOLD; bvec4 test;     \n"
01435         "       bins[9] = bins[0];                                                              \n"
01436         "       float npeaks = 0.0, k = 0;                                              \n"
01437         "       float prevb     = bins[8].w;                                            \n"
01438         "       for (int i = 0; i <9 ; i++)                                             \n"
01439         "       {\n"
01440         "               FINDPEAK(i, k);\n"
01441         "               k = k + 4.0;    \n"
01442         "       }";
01443         }else
01444         {
01445                 //loop unroll for ATI.
01446         out <<"\n"
01447         "       vec4 hh = maxh * ORIENTATION_THRESHOLD; bvec4 test;\n"
01448         "       bins[9] = bins[0];                                                              \n"
01449         "       float npeaks = 0.0;                                                             \n"
01450         "       float prevb     = bins[8].w;                                            \n"
01451         "       FINDPEAK(0, 0.0);\n"
01452         "       FINDPEAK(1, 4.0);\n"
01453         "       FINDPEAK(2, 8.0);\n"
01454         "       FINDPEAK(3, 12.0);\n"
01455         "       FINDPEAK(4, 16.0);\n"
01456         "       FINDPEAK(5, 20.0);\n"
01457         "       FINDPEAK(6, 24.0);\n"
01458         "       FINDPEAK(7, 28.0);\n"
01459         "       FINDPEAK(8, 32.0);\n";
01460         }
01461         //WRITE output
01462         if(GlobalUtil::_MaxOrientation>1)
01463         {
01464         out<<"\n"
01465         "       if(orientation_mode){\n"
01466         "               npeaks = dot(vec4(1,1,"
01467                         <<(GlobalUtil::_MaxOrientation>2 ? 1 : 0)<<","
01468                         <<(GlobalUtil::_MaxOrientation >3? 1 : 0)<<"), vec4(greaterThan(weights, hh)));\n"
01469         "               gl_FragData[1] = radians((Orientations )*10.0);\n"
01470         "               gl_FragData[0] = vec4(pos, npeaks, sigma);\n"
01471         "       }else{\n"
01472         "               gl_FragData[0] = vec4(pos, radians((Orientations.x)*10.0), sigma);\n"
01473         "       }\n";
01474         }else
01475         {
01476         out<<"\n"
01477         "        gl_FragData[0] = vec4(pos, radians((Orientation.x)*10.0), sigma);\n";
01478         }
01479         //end
01480         out<<"\n"
01481         "}\n"<<'\0';
01482 
01483 
01484 }
01485 
01486 void ShaderBagGLSL::SetSimpleOrientationInput(int oTex, float sigma, float sigma_step)
01487 {
01488         glUniform1i(_param_orientation_gtex, 1);
01489         glUniform1f(_param_orientation_size, sigma);
01490 }
01491 
01492 
01493 
01494 
01495 void ShaderBagGLSL::SetFeatureOrientationParam(int gtex, int width, int height, float sigma, int stex, float step)
01496 {
01498         glUniform1i(_param_orientation_gtex, 1);        
01499 
01500         if((GlobalUtil::_SubpixelLocalization || GlobalUtil::_KeepExtremumSign)&& stex)
01501         {
01502                 //specify texutre for subpixel subscale localization
01503                 glUniform1i(_param_orientation_stex, 2);
01504         }
01505 
01506         float size[4];
01507         size[0] = (float)width;
01508         size[1] = (float)height;
01509         size[2] = sigma;
01510         size[3] = step;
01511         glUniform4fv(_param_orientation_size, 1, size);
01512 
01513 }
01514 
01515 
01516 void ShaderBagGLSL::LoadDescriptorShaderF2()
01517 {
01518         //one shader outpout 128/8 = 16 , each fragout encodes 4
01519         //const double twopi = 2.0*3.14159265358979323846;
01520         //const double rpi  = 8.0/twopi;
01521         char buffer[10240];
01522         ostrstream out(buffer, 10240);
01523 
01524         out<<setprecision(8);
01525 
01526         out<<"\n"
01527         "#define M_PI 3.14159265358979323846\n"
01528         "#define TWO_PI (2.0*M_PI)\n"
01529         "#define RPI 1.2732395447351626861510701069801\n"
01530         "#define WF  size.z\n"
01531         "uniform sampler2DRect tex;                             \n"
01532         "uniform sampler2DRect gradTex;                 \n"
01533         "uniform vec4 dsize;                                            \n"
01534         "uniform vec3 size;                                             \n"
01535         "void main()            \n"
01536         "{\n"
01537         "       vec2 dim        = size.xy;      //image size                    \n"
01538         "       float index = dsize.x*floor(gl_TexCoord[0].y * 0.5) + gl_TexCoord[0].x;\n"
01539         "       float idx = 8.0 * fract(index * 0.125) + 8.0 * floor(2.0 * fract(gl_TexCoord[0].y * 0.5));              \n"
01540         "       index = floor(index*0.125) + 0.49;  \n"
01541         "       vec2 coord = floor( vec2( mod(index, dsize.z), index*dsize.w)) + 0.5 ;\n"
01542         "       vec2 pos = texture2DRect(tex, coord).xy;                \n"
01543         "       if(any(lessThanEqual(pos.xy,  vec2(1.0))) || any(greaterThanEqual(pos.xy, dim-1.0)))// discard; \n"
01544         "       { gl_FragData[0] = gl_FragData[1] = vec4(0.0); return; }\n"
01545         "       float  anglef = texture2DRect(tex, coord).z;\n"
01546         "       if(anglef > M_PI) anglef -= TWO_PI;\n"
01547         "       float sigma = texture2DRect(tex, coord).w; \n"
01548         "       float spt  = abs(sigma * WF);   //default to be 3*sigma \n";
01549 
01550         //rotation
01551         out<<
01552         "       vec4 cscs, rots;                                                                \n"
01553         "       cscs.y = sin(anglef);   cscs.x = cos(anglef);   \n"
01554         "       cscs.zw = - cscs.xy;                                                    \n"
01555         "       rots = cscs /spt;                                                               \n"
01556         "       cscs *= spt; \n";
01557 
01558         //here cscs is actually (cos, sin, -cos, -sin) * (factor: 3)*sigma
01559         //and rots is  (cos, sin, -cos, -sin ) /(factor*sigma)
01560         //devide the 4x4 sift grid into 16 1x1 block, and each corresponds to a shader thread
01561         //To use linear interoplation, 1x1 is increased to 2x2, by adding 0.5 to each side
01562 
01563         out<<
01564         "vec4 temp; vec2 pt, offsetpt;                          \n"
01565         "       /*the fraction part of idx is .5*/                      \n"
01566         "       offsetpt.x = 4.0* fract(idx*0.25) - 2.0;                                \n"
01567         "       offsetpt.y = floor(idx*0.25) - 1.5;                     \n"
01568         "       temp = cscs.xwyx*offsetpt.xyxy;                         \n"
01569         "       pt = pos + temp.xz + temp.yw;                           \n";
01570         
01571         //get a horizontal bounding box of the rotated rectangle
01572         out<<
01573         "       vec2 bwin = abs(cscs.xy);                                       \n"
01574         "       float bsz = bwin.x + bwin.y;                                    \n"
01575         "       vec4 sz;                                        \n"
01576         "       sz.xy = max(pt - vec2(bsz), vec2(1,1));\n"
01577         "       sz.zw = min(pt + vec2(bsz), dim - 2);           \n"
01578         "       sz = floor(sz)+0.5;"; //move sample point to pixel center
01579         //get voting for two box
01580 
01581         out<<"\n"
01582         "       vec4 DA, DB; vec2 spos;                 \n"
01583         "       DA = DB  = vec4(0, 0, 0, 0);            \n"
01584         "       for(spos.y = sz.y; spos.y <= sz.w;      spos.y+=1.0)                            \n"
01585         "       {                                                                                                                               \n"
01586         "               for(spos.x = sz.x; spos.x <= sz.z;      spos.x+=1.0)                    \n"
01587         "               {                                                                                                                       \n"
01588         "                       vec2 diff = spos - pt;                                                          \n"
01589         "                       temp = rots.xywx * diff.xyxy;\n"
01590         "                       vec2 nxy = (temp.xz + temp.yw); \n"
01591         "                       vec2 nxyn = abs(nxy);                   \n"
01592         "                       if(all( lessThan(nxyn, vec2(1.0)) ))\n"
01593         "                       {\n"
01594         "                               vec4 cc = texture2DRect(gradTex, spos);                                         \n"
01595         "                               float mod = cc.b;       float angle = cc.a;                                     \n"
01596         "                               float theta0 = RPI * (anglef - angle);                          \n"
01597         "                               float theta = theta0 < 0.0? theta0 + 8.0 : theta0;;\n"
01598         "                               diff = nxy + offsetpt.xy;                                                               \n"
01599         "                               float ww = exp(-0.125*dot(diff, diff));\n"
01600         "                               vec2 weights = 1 - nxyn;\n"
01601         "                               float weight = weights.x * weights.y *mod*ww; \n"
01602         "                               float theta1 = floor(theta); \n"
01603         "                               float weight2 = (theta - theta1) * weight;\n"
01604         "                               float weight1 = weight - weight2;\n"
01605         "                               DA += vec4(equal(vec4(theta1),  vec4(0, 1, 2, 3)))*weight1;\n"
01606         "                               DA += vec4(equal(vec4(theta1),  vec4(7, 0, 1, 2)))*weight2; \n"
01607         "                               DB += vec4(equal(vec4(theta1),  vec4(4, 5, 6, 7)))*weight1;\n"
01608         "                               DB += vec4(equal(vec4(theta1),  vec4(3, 4, 5, 6)))*weight2; \n"
01609         "                       }\n"
01610         "               }\n"
01611         "       }\n";
01612 
01613         out<<
01614         "        gl_FragData[0] = DA; gl_FragData[1] = DB;\n"
01615         "}\n"<<'\0';
01616 
01617         ProgramGLSL * program =  new ProgramGLSL(buffer); 
01618 
01619         if(program->IsNative())
01620         {
01621                 s_descriptor_fp = program ;
01622                 _param_descriptor_gtex = glGetUniformLocation(*program, "gradTex");
01623                 _param_descriptor_size = glGetUniformLocation(*program, "size");
01624                 _param_descriptor_dsize = glGetUniformLocation(*program, "dsize");
01625         }else
01626         {
01627                 delete program;
01628         }
01629 
01630 
01631 }
01632 
01633 void ShaderBagGLSL::LoadDescriptorShader()
01634 {
01635         GlobalUtil::_DescriptorPPT = 16;
01636         LoadDescriptorShaderF2();
01637 }
01638 
01639 
01640 void ShaderBagGLSL::SetFeatureDescirptorParam(int gtex, int otex, float dwidth, float fwidth,  float width, float height, float sigma)
01641 {
01643         glUniform1i(_param_descriptor_gtex, 1); 
01644 
01645         float dsize[4] ={dwidth, 1.0f/dwidth, fwidth, 1.0f/fwidth};
01646         glUniform4fv(_param_descriptor_dsize, 1, dsize);
01647         float size[3];
01648         size[0] = width;
01649         size[1] = height;
01650         size[2] = GlobalUtil::_DescriptorWindowFactor;
01651         glUniform3fv(_param_descriptor_size, 1, size);
01652 
01653 }
01654 
01656 
01657 void ShaderBagPKSL::LoadFixedShaders()
01658 {
01659         ProgramGLSL * program;
01660 
01661 
01662         s_gray = new ProgramGLSL( 
01663         "uniform sampler2DRect tex; void main(){\n"
01664         "float intensity = dot(vec3(0.299, 0.587, 0.114), texture2DRect(tex,gl_TexCoord[0].xy ).rgb);\n"
01665         "gl_FragColor= vec4(intensity, intensity, intensity, 1.0);}"    );
01666 
01667 
01668         s_sampling = new ProgramGLSL(
01669         "uniform sampler2DRect tex; void main(){\n"
01670         "gl_FragColor= vec4(    texture2DRect(tex,gl_TexCoord[0].st ).r,texture2DRect(tex,gl_TexCoord[1].st ).r,\n"
01671         "                                               texture2DRect(tex,gl_TexCoord[2].st ).r,texture2DRect(tex,gl_TexCoord[3].st ).r);}"     );
01672 
01673 
01674         s_margin_copy = program = new ProgramGLSL(
01675         "uniform sampler2DRect tex;  uniform vec4 truncate; void main(){\n"
01676         "vec4 cc = texture2DRect(tex, min(gl_TexCoord[0].xy, truncate.xy)); \n"
01677         "bvec2 ob = lessThan(gl_TexCoord[0].xy, truncate.xy);\n"
01678         "if(ob.y) { gl_FragColor = (truncate.z ==0 ? cc.rrbb : cc.ggaa); } \n"
01679         "else if(ob.x) {gl_FragColor = (truncate.w <1.5 ? cc.rgrg : cc.baba);} \n"
01680         "else { vec4 weights = vec4(vec4(0, 1, 2, 3) == truncate.w);\n"
01681         "float v = dot(weights, cc); gl_FragColor = vec4(v);}}");
01682 
01683         _param_margin_copy_truncate = glGetUniformLocation(*program, "truncate");
01684 
01685 
01686 
01687         s_zero_pass = new ProgramGLSL("void main(){gl_FragColor = vec4(0.0);}");
01688 
01689 
01690 
01691         s_grad_pass = program = new ProgramGLSL(
01692         "uniform sampler2DRect tex; uniform sampler2DRect texp; void main ()\n"
01693         "{\n"
01694         "       vec4 v1, v2, gg;\n"
01695         "       vec4 cc = texture2DRect(tex, gl_TexCoord[0].xy);\n"
01696         "       vec4 cp = texture2DRect(texp, gl_TexCoord[0].xy);\n"
01697         "       gl_FragData[0] = cc - cp; \n"
01698         "       vec4 cl = texture2DRect(tex, gl_TexCoord[1].xy); vec4 cr = texture2DRect(tex, gl_TexCoord[2].xy);\n"
01699         "       vec4 cd = texture2DRect(tex, gl_TexCoord[3].xy); vec4 cu = texture2DRect(tex, gl_TexCoord[4].xy);\n"
01700         "       vec4 dx = (vec4(cr.rb, cc.ga) - vec4(cc.rb, cl.ga)).zxwy;\n"
01701         "       vec4 dy = (vec4(cu.rg, cc.ba) - vec4(cc.rg, cd.ba)).zwxy;\n"
01702         "       vec4 grad = 0.5 * sqrt(dx*dx + dy * dy);\n"
01703         "       gl_FragData[1] = grad;\n"
01704         "       vec4 invalid = vec4(equal(grad, vec4(0.0)));    \n"
01705         "       vec4 ov = atan(dy, dx + invalid);               \n"
01706         "       gl_FragData[2] = ov; \n"
01707         "}\n\0"); //when 
01708 
01709         _param_grad_pass_texp = glGetUniformLocation(*program, "texp");
01710 
01711 
01712         GlobalUtil::_OrientationPack2 = 0;
01713         LoadOrientationShader();
01714 
01715         if(s_orientation == NULL)
01716         {
01717                 //Load a simplified version if the right version is not supported
01718                 s_orientation = program =  new ProgramGLSL(
01719                 "uniform sampler2DRect fTex; uniform sampler2DRect oTex; uniform vec2 size; void main(){\n"
01720                 "       vec4 cc = texture2DRect(fTex, gl_TexCoord[0].xy);\n"
01721                 "       vec2 co = cc.xy * 0.5; \n"
01722                 "       vec4 oo = texture2DRect(oTex, co);\n"
01723                 "       bvec2 bo = lessThan(fract(co), vec2(0.5)); \n"
01724                 "       float o = bo.y? (bo.x? oo.r : oo.g) : (bo.x? oo.b : oo.a); \n"
01725                 "       gl_FragColor = vec4(cc.rg, o, size.x * pow(size.y, cc.a));}");  
01726                 _param_orientation_gtex= glGetUniformLocation(*program, "oTex");
01727                 _param_orientation_size= glGetUniformLocation(*program, "size");
01728 
01729                 GlobalUtil::_MaxOrientation = 0;
01730                 GlobalUtil::_FullSupported = 0;
01731                 std::cerr<<"Orientation simplified on this hardware"<<endl;
01732         }
01733 
01734         if(GlobalUtil::_DescriptorPPT)
01735         {
01736                 LoadDescriptorShader();
01737                 if(s_descriptor_fp == NULL) 
01738                 {
01739                         GlobalUtil::_DescriptorPPT = GlobalUtil::_FullSupported = 0; 
01740                         std::cerr<<"Descriptor ignored on this hardware"<<endl;
01741                 }
01742         }
01743 }
01744 
01745 
01746 void ShaderBagPKSL::LoadDisplayShaders()
01747 {
01748         ProgramGLSL * program;
01749 
01750         s_copy_key = new ProgramGLSL(
01751         "uniform sampler2DRect tex;void main(){\n"
01752         "gl_FragColor= vec4(texture2DRect(tex, gl_TexCoord[0].xy).rg, 0,1);}");
01753 
01754         //shader used to write a vertex buffer object
01755         //which is used to draw the quads of each feature
01756         s_vertex_list = program = new ProgramGLSL(
01757         "uniform sampler2DRect tex; uniform vec4 sizes; void main(){\n"
01758         "float fwidth = sizes.y; \n"
01759         "float twidth = sizes.z; \n"
01760         "float rwidth = sizes.w; \n"
01761         "float index = 0.1*(fwidth*floor(gl_TexCoord[0].y) + gl_TexCoord[0].x);\n"
01762         "float px = mod(index, twidth);\n"
01763         "vec2 tpos= floor(vec2(px, index*rwidth))+0.5;\n"
01764         "vec4 cc = texture2DRect(tex, tpos );\n"
01765         "float size = 3.0f * cc.a; \n"
01766         "gl_FragColor.zw = vec2(0.0, 1.0);\n"
01767         "if(any(lessThan(cc.xy,vec2(0.0)))) {gl_FragColor.xy = cc.xy;}else \n"
01768         "{\n"
01769         "       float type = fract(px);\n"
01770         "       vec2 dxy; float s, c;\n"
01771         "       dxy.x = type < 0.1 ? 0.0 : ((type <0.5 || type > 0.9)? size : -size);\n"
01772         "       dxy.y = type < 0.2 ? 0.0 : ((type < 0.3 || type > 0.7 )? -size :size); \n"
01773         "       s = sin(cc.b); c = cos(cc.b); \n"
01774         "       gl_FragColor.x = cc.x + c*dxy.x-s*dxy.y;\n"
01775         "       gl_FragColor.y = cc.y + c*dxy.y+s*dxy.x;}\n"
01776         "}\n\0");
01777         /*gl_FragColor = vec4(tpos, 0.0, 1.0);}\n\0");*/
01778 
01779         _param_genvbo_size = glGetUniformLocation(*program, "sizes");
01780 
01781         s_display_gaussian = new ProgramGLSL(
01782         "uniform sampler2DRect tex; void main(){\n"
01783     "vec4 pc = texture2DRect(tex, gl_TexCoord[0].xy);   bvec2 ff = lessThan(fract(gl_TexCoord[0].xy), vec2(0.5));\n"
01784     "float v = ff.y?(ff.x? pc.r : pc.g):(ff.x?pc.b:pc.a); gl_FragColor = vec4(vec3(v), 1.0);}");
01785 
01786         s_display_dog =  new ProgramGLSL(
01787         "uniform sampler2DRect tex; void main(){\n"
01788         "vec4 pc = texture2DRect(tex, gl_TexCoord[0].xy); bvec2 ff = lessThan(fract(gl_TexCoord[0].xy), vec2(0.5));\n"
01789         "float v = ff.y ?(ff.x ? pc.r : pc.g):(ff.x ? pc.b : pc.a);float g = (0.5+20.0*v);\n"
01790         "gl_FragColor = vec4(g, g, g, 1.0);}" );
01791 
01792 
01793         s_display_grad = new ProgramGLSL(
01794         "uniform sampler2DRect tex; void main(){\n"
01795         "vec4 pc = texture2DRect(tex, gl_TexCoord[0].xy); bvec2 ff = lessThan(fract(gl_TexCoord[0].xy), vec2(0.5));\n"
01796         "float v = ff.y ?(ff.x ? pc.r : pc.g):(ff.x ? pc.b : pc.a); gl_FragColor = vec4(5.0 *vec3(v), 1.0); }");
01797 
01798         s_display_keys= new ProgramGLSL(
01799         "uniform sampler2DRect tex; void main(){\n"
01800         "vec4 oc = texture2DRect(tex, gl_TexCoord[0].xy); \n"
01801         "vec4 cc = vec4(equal(abs(oc.rrrr), vec4(1.0, 2.0, 3.0, 4.0))); \n"
01802         "bvec2 ff = lessThan(fract(gl_TexCoord[0].xy) , vec2(0.5));\n"
01803         "float v = ff.y ?(ff.x ? cc.r : cc.g):(ff.x ? cc.b : cc.a);\n"
01804         "if(v == 0) discard;    \n"
01805         "else if(oc.r > 0) gl_FragColor = vec4(1.0, 0, 0,1.0); \n"
01806         "else gl_FragColor = vec4(0.0,1.0,0.0,1.0);     }" );
01807 }
01808 
01809 void ShaderBagPKSL::LoadOrientationShader(void)
01810 {
01811         char buffer[10240];
01812         ostrstream out(buffer,10240);
01813         if(GlobalUtil::_IsNvidia)
01814         {
01815                 out <<  "#pragma optionNV(ifcvt none)\n"
01816                                 "#pragma optionNV(unroll all)\n";
01817         }
01818         out<<"\n"
01819         "#define GAUSSIAN_WF "<<GlobalUtil::_OrientationGaussianFactor<<" \n"
01820         "#define SAMPLE_WF ("<<GlobalUtil::_OrientationWindowFactor<< " )\n"
01821         "#define ORIENTATION_THRESHOLD "<< GlobalUtil::_MulitiOrientationThreshold << "\n"
01822         "uniform sampler2DRect tex;     uniform sampler2DRect gtex;\n"
01823         "uniform sampler2DRect otex; uniform vec4 size;\n"
01824         "void main()            \n"
01825         "{                                                                                                      \n"
01826         "       vec4 bins[10];                                                          \n"
01827         "       bins[0] = vec4(0.0);bins[1] = vec4(0.0);bins[2] = vec4(0.0);    \n"
01828         "       bins[3] = vec4(0.0);bins[4] = vec4(0.0);bins[5] = vec4(0.0);    \n"
01829         "       bins[6] = vec4(0.0);bins[7] = vec4(0.0);bins[8] = vec4(0.0);    \n"
01830         "       vec4 sift = texture2DRect(tex, gl_TexCoord[0].xy);      \n"
01831         "       vec2 pos = sift.xy; \n"
01832         "       bool orientation_mode = (size.z != 0);          \n"
01833         "       float sigma = orientation_mode? (abs(size.z) * pow(size.w, sift.w) * sift.z) : (sift.w); \n"
01834         "       //bool fixed_orientation = (size.z < 0);                \n"
01835         "       if(size.z < 0) {gl_FragData[0] = vec4(pos, 0, sigma); return;}"
01836         "       float gsigma = sigma * GAUSSIAN_WF;                             \n"
01837         "       vec2 win = abs(vec2(sigma)) * (SAMPLE_WF * GAUSSIAN_WF);        \n"
01838         "       vec2 dim = size.xy;                                                     \n"
01839         "       vec4 dist_threshold = vec4(win.x*win.x+0.5);                    \n"
01840         "       float factor = -0.5/(gsigma*gsigma);                    \n"
01841         "       vec4 sz;        vec2 spos;                                              \n"
01842         "       //if(any(pos.xy <= 1)) discard;                                 \n"
01843         "       sz.xy = max( pos - win, vec2(2,2));                     \n"
01844         "       sz.zw = min( pos + win, dim-3);                         \n"
01845         "       sz = floor(sz*0.5) + 0.5; ";
01846                 //loop to get the histogram
01847 
01848         out<<"\n"
01849         "       for(spos.y = sz.y; spos.y <= sz.w;      spos.y+=1.0)                            \n"
01850         "       {                                                                                                                               \n"
01851         "               for(spos.x = sz.x; spos.x <= sz.z;      spos.x+=1.0)                    \n"
01852         "               {                                                                                                                       \n"
01853         "                       vec2 offset = 2* spos - pos - 0.5;                                      \n"
01854         "                       vec4 off = vec4(offset, offset + 1);                            \n"
01855         "                       vec4 distsq = off.xzxz * off.xzxz + off.yyww * off.yyww;        \n"
01856         "                       bvec4 inside = lessThan(distsq, dist_threshold);                        \n"
01857         "                       if(any(inside))                                                                         \n"
01858         "                       {                                                                                                               \n"
01859         "                               vec4 gg = texture2DRect(gtex, spos);                            \n"
01860         "                               vec4 oo = texture2DRect(otex, spos);                            \n"
01861         "                               vec4 weight = gg * exp(distsq * factor);                        \n"
01862         "                               vec4 idxv  = floor(degrees(oo)*0.1);                            \n"
01863         "                               idxv+= (vec4(lessThan(idxv, vec4(0.0)))*36.0);                  \n"
01864         "                               vec4 vidx = fract(idxv * 0.25) * 4.0;//mod(idxv, 4.0);  \n";
01865         //
01866         if(GlobalUtil::_UseDynamicIndexing && GlobalUtil::_IsNvidia)
01867         {
01868                 // it might be slow on some GPUs
01869                 out<<"\n"
01870         "                               for(int i = 0 ; i < 4; i++)\n"
01871         "                               {\n"
01872         "                                       if(inside[i])\n"
01873         "                                       {\n"
01874         "                                               float idx = idxv[i];                                                            \n"
01875         "                                               vec4 inc = weight[i] * vec4(equal(vec4(vidx[i]), vec4(0,1,2,3)));       \n"
01876         "                                               int iidx = int(floor(idx*0.25));        \n"
01877         "                                               bins[iidx]+=inc;                                        \n"
01878         "                                       }                                                                               \n"
01879         "                               }                                                                                       \n"
01880         "                       }                                                                                               \n"
01881         "               }                                                                                                       \n"
01882         "       }";
01883 
01884         }else
01885         {
01886                 //nvfp40 still does not support dynamic array indexing
01887                 //unrolled binary search
01888                 //it seems to be faster than the dyanmic indexing version on some GPUs
01889                 out<<"\n"
01890         "                               for(int i = 0 ; i < 4; i++)\n"
01891         "                               {\n"
01892         "                                       if(inside[i])\n"
01893         "                                       {\n"
01894         "                                               float idx = idxv[i];                                                                            \n"
01895         "                                               vec4 inc = weight[i] * vec4(equal(vec4(vidx[i]), vec4(0,1,2,3)));       \n"
01896         "                                               if(idx < 16)                                                    \n"
01897         "                                               {                                                                               \n"
01898         "                                                       if(idx < 8)                                                     \n"
01899         "                                                       {                                                                       \n"
01900         "                                                               if(idx < 4)     {       bins[0]+=inc;}  \n"
01901         "                                                               else            {       bins[1]+=inc;}  \n"
01902         "                                                       }else                                                           \n"
01903         "                                                       {                                                                       \n"
01904         "                                                               if(idx < 12){   bins[2]+=inc;}  \n"
01905         "                                                               else            {       bins[3]+=inc;}  \n"
01906         "                                                       }                                                                       \n"
01907         "                                               }else if(idx < 32)                                              \n"
01908         "                                               {                                                                               \n"
01909         "                                                       if(idx < 24)                                            \n"
01910         "                                                       {                                                                       \n"
01911         "                                                               if(idx <20)     {       bins[4]+=inc;}  \n"
01912         "                                                               else            {       bins[5]+=inc;}  \n"
01913         "                                                       }else                                                           \n"
01914         "                                                       {                                                                       \n"
01915         "                                                               if(idx < 28){   bins[6]+=inc;}  \n"
01916         "                                                               else            {       bins[7]+=inc;}  \n"
01917         "                                                       }                                                                       \n"
01918         "                                               }else                                           \n"
01919         "                                               {                                                                               \n"
01920         "                                                       bins[8]+=inc;                                           \n"
01921         "                                               }                                                                               \n"
01922         "                                       }                                                                                       \n"
01923         "                               }                                                                                               \n"
01924         "                       }                                                                               \n"
01925         "               }                                                                                       \n"
01926         "       }";
01927 
01928         }
01929 
01930         //reuse the code from the unpacked version..
01931         ShaderBagGLSL::WriteOrientationCodeToStream(out);
01932 
01933 
01934 
01935         ProgramGLSL * program = new ProgramGLSL(buffer);
01936         if(program->IsNative())
01937         {
01938                 s_orientation = program ;
01939                 _param_orientation_gtex = glGetUniformLocation(*program, "gtex");
01940                 _param_orientation_otex = glGetUniformLocation(*program, "otex");
01941                 _param_orientation_size = glGetUniformLocation(*program, "size");
01942         }else
01943         {
01944                 delete program;
01945         }
01946 }
01947 
01948 void ShaderBagPKSL::SetGenListStartParam(float width, int tex0)
01949 {
01950         glUniform1f(_param_ftex_width, width);
01951 }
01952 
01953 void ShaderBagPKSL::LoadGenListShader(int ndoglev,int nlev)
01954 {
01955         ProgramGLSL * program;
01956 
01957         s_genlist_init_tight = new ProgramGLSL(
01958         "uniform sampler2DRect tex; void main ()\n"
01959         "{\n"
01960         "       vec4 key = vec4(texture2DRect(tex, gl_TexCoord[0].xy).r, \n"
01961         "                                       texture2DRect(tex, gl_TexCoord[1].xy).r, \n"
01962         "                                       texture2DRect(tex, gl_TexCoord[2].xy).r, \n"
01963         "                                       texture2DRect(tex, gl_TexCoord[3].xy).r); \n"
01964         "                                       gl_FragColor = vec4(notEqual(key, vec4(0.0))); \n"
01965         "}");
01966 
01967         s_genlist_init_ex = program = new ProgramGLSL(
01968         "uniform sampler2DRect tex; uniform vec4 bbox; void main ()\n"
01969         "{\n"
01970         "       vec4 helper1 = vec4(equal(vec4(abs(texture2DRect(tex, gl_TexCoord[0].xy).r)), vec4(1.0, 2.0, 3.0, 4.0)));\n"
01971         "       vec4 helper2 = vec4(equal(vec4(abs(texture2DRect(tex, gl_TexCoord[1].xy).r)), vec4(1.0, 2.0, 3.0, 4.0)));\n"
01972         "       vec4 helper3 = vec4(equal(vec4(abs(texture2DRect(tex, gl_TexCoord[2].xy).r)), vec4(1.0, 2.0, 3.0, 4.0)));\n"
01973         "       vec4 helper4 = vec4(equal(vec4(abs(texture2DRect(tex, gl_TexCoord[3].xy).r)), vec4(1.0, 2.0, 3.0, 4.0)));\n"
01974         "       vec4 bx1 = vec4(lessThan(gl_TexCoord[0].xxyy, bbox)); \n"
01975         "       vec4 bx4 = vec4(lessThan(gl_TexCoord[3].xxyy, bbox)); \n"
01976         "       vec4 bx2 = vec4(bx4.xy, bx1.zw); \n"
01977         "       vec4 bx3 = vec4(bx1.xy, bx4.zw);\n"
01978         "       helper1 = min(min(bx1.xyxy, bx1.zzww), helper1);\n"
01979         "       helper2 = min(min(bx2.xyxy, bx2.zzww), helper2);\n"
01980         "       helper3 = min(min(bx3.xyxy, bx3.zzww), helper3);\n"
01981         "       helper4 = min(min(bx4.xyxy, bx4.zzww), helper4);\n"
01982         "       gl_FragColor.r = float(any(greaterThan(max(helper1.xy, helper1.zw), vec2(0.0))));       \n"
01983         "       gl_FragColor.g = float(any(greaterThan(max(helper2.xy, helper2.zw), vec2(0.0))));       \n"
01984         "       gl_FragColor.b = float(any(greaterThan(max(helper3.xy, helper3.zw), vec2(0.0))));       \n"
01985         "       gl_FragColor.a = float(any(greaterThan(max(helper4.xy, helper4.zw), vec2(0.0))));       \n"
01986         "}");
01987         _param_genlist_init_bbox = glGetUniformLocation( *program, "bbox");
01988 
01989         s_genlist_end = program = new ProgramGLSL(
01990                 GlobalUtil::_KeepExtremumSign == 0 ? 
01991         
01992         "uniform sampler2DRect tex; uniform sampler2DRect ktex; void main()\n"
01993         "{\n"
01994         "       vec4 tc = texture2DRect( tex, gl_TexCoord[0].xy);\n"
01995         "       vec2 pos = tc.rg; float index = tc.b;\n"
01996         "       vec4 tk = texture2DRect( ktex, pos); \n"
01997         "       vec4 keys = vec4(equal(abs(tk.rrrr), vec4(1.0, 2.0, 3.0, 4.0))); \n"
01998         "       vec2 opos; \n"
01999         "       opos.x = dot(keys, vec4(-0.5, 0.5, -0.5, 0.5));\n"
02000         "       opos.y = dot(keys, vec4(-0.5, -0.5, 0.5, 0.5));\n"
02001         "       gl_FragColor = vec4(opos + pos * 2.0 + tk.yz, 1.0, tk.w);\n"
02002         "}" : 
02003         
02004         "uniform sampler2DRect tex; uniform sampler2DRect ktex; void main()\n"
02005         "{\n"
02006         "       vec4 tc = texture2DRect( tex, gl_TexCoord[0].xy);\n"
02007         "       vec2 pos = tc.rg; float index = tc.b;\n"
02008         "       vec4 tk = texture2DRect( ktex, pos); \n"
02009         "       vec4 keys = vec4(equal(abs(tk.rrrr), vec4(1.0, 2.0, 3.0, 4.0))) \n"
02010         "       vec2 opos; \n"
02011         "       opos.x = dot(keys, vec4(-0.5, 0.5, -0.5, 0.5));\n"
02012         "       opos.y = dot(keys, vec4(-0.5, -0.5, 0.5, 0.5));\n"
02013         "       gl_FragColor = vec4(opos + pos * 2.0 + tk.yz, sign(tk.r), tk.w);\n"
02014         "}"     
02015         );
02016         _param_genlist_end_ktex = glGetUniformLocation(*program, "ktex");
02017 
02018         //reduction ...
02019         s_genlist_histo = new ProgramGLSL(
02020         "uniform sampler2DRect tex; void main ()\n"
02021         "{\n"
02022         "       vec4 helper; vec4 helper2; \n"
02023         "       helper = texture2DRect(tex, gl_TexCoord[0].xy); helper2.xy = helper.xy + helper.zw; \n"
02024         "       helper = texture2DRect(tex, gl_TexCoord[1].xy); helper2.zw = helper.xy + helper.zw; \n"
02025         "       gl_FragColor.rg = helper2.xz + helper2.yw;\n"
02026         "       helper = texture2DRect(tex, gl_TexCoord[2].xy); helper2.xy = helper.xy + helper.zw; \n"
02027         "       helper = texture2DRect(tex, gl_TexCoord[3].xy); helper2.zw = helper.xy + helper.zw; \n"
02028         "       gl_FragColor.ba= helper2.xz+helper2.yw;\n"
02029         "}");
02030 
02031 
02032         //read of the first part, which generates tex coordinates 
02033 
02034         s_genlist_start= program =  ShaderBagGLSL::LoadGenListStepShader(1, 1);
02035         _param_ftex_width= glGetUniformLocation(*program, "width");
02036         _param_genlist_start_tex0 = glGetUniformLocation(*program, "tex0");
02037         //stepping
02038         s_genlist_step = program = ShaderBagGLSL::LoadGenListStepShader(0, 1);
02039         _param_genlist_step_tex= glGetUniformLocation(*program, "tex");
02040         _param_genlist_step_tex0= glGetUniformLocation(*program, "tex0");
02041 
02042 }
02043 void ShaderBagPKSL::UnloadProgram(void)
02044 {
02045         glUseProgram(0);
02046 }
02047 void ShaderBagPKSL::LoadKeypointShader(float dog_threshold, float edge_threshold)
02048 {
02049         //
02050         char buffer[10240];
02051         float threshold0 = dog_threshold* (GlobalUtil::_SubpixelLocalization?0.8f:1.0f);
02052         float threshold1 = dog_threshold;
02053         float threshold2 = (edge_threshold+1)*(edge_threshold+1)/edge_threshold;
02054         ostrstream out(buffer, 10240);
02055         out<<setprecision(8);
02056 
02057         if(GlobalUtil::_IsNvidia)
02058         {
02059                 out << "#pragma optionNV(ifcvt none)\n"
02060                                 "#pragma optionNV(unroll all)\n"
02061                                 "#define REPEAT4(FUNCTION)\\\n"
02062                                 "for(int i = 0; i < 4; ++i)\\\n"
02063                                 "{\\\n"
02064                                 "       FUNCTION(i);\\\n"
02065                                 "}\n";
02066         }else
02067         {
02068                 //loop unroll
02069                 out <<  "#define REPEAT4(FUNCTION)\\\n"
02070                                 "FUNCTION(0);\\\n"
02071                                 "FUNCTION(1);\\\n"
02072                                 "FUNCTION(2);\\\n"
02073                                 "FUNCTION(3);\n";
02074         }
02075         //tex(X)(Y)
02076         //X: (CLR) (CENTER 0, LEFT -1, RIGHT +1)  
02077         //Y: (CDU) (CENTER 0, DOWN -1, UP    +1) 
02078         out <<  "#define THRESHOLD0 " << threshold0 << "\n"
02079                         "#define THRESHOLD1 " << threshold1 << "\n"
02080                         "#define THRESHOLD2 " << threshold2 << "\n";
02081 
02082         out<<
02083         "uniform sampler2DRect tex; uniform sampler2DRect texU;\n"
02084         "uniform sampler2DRect texD; void main ()\n"
02085         "{\n"
02086         "       vec2 TexRU = vec2(gl_TexCoord[2].x, gl_TexCoord[4].y); \n"
02087         "       vec4 ccc = texture2DRect(tex, gl_TexCoord[0].xy);\n"
02088         "       vec4 clc = texture2DRect(tex, gl_TexCoord[1].xy);\n"
02089         "       vec4 crc = texture2DRect(tex, gl_TexCoord[2].xy);\n"
02090         "       vec4 ccd = texture2DRect(tex, gl_TexCoord[3].xy);\n"
02091         "       vec4 ccu = texture2DRect(tex, gl_TexCoord[4].xy);\n"
02092         "       vec4 cld = texture2DRect(tex, gl_TexCoord[5].xy);\n"
02093         "       vec4 clu = texture2DRect(tex, gl_TexCoord[6].xy);\n"
02094         "       vec4 crd = texture2DRect(tex, gl_TexCoord[7].xy);\n"
02095         "       vec4 cru = texture2DRect(tex, TexRU.xy);\n"
02096         "       vec4  cc = ccc;\n"
02097         "       vec4  v1[4], v2[4];\n"
02098         "       v1[0] = vec4(clc.g, ccc.g, ccd.b, ccc.b);\n"
02099         "       v1[1] = vec4(ccc.r, crc.r, ccd.a, ccc.a);\n"
02100         "       v1[2] = vec4(clc.a, ccc.a, ccc.r, ccu.r);\n"
02101         "       v1[3] = vec4(ccc.b, crc.b, ccc.g, ccu.g);\n"
02102         "       v2[0] = vec4(cld.a, clc.a, ccd.a, ccc.a);\n"
02103         "       v2[1] = vec4(ccd.b, ccc.b, crd.b, crc.b);\n"
02104         "       v2[2] = vec4(clc.g, clu.g, ccc.g, ccu.g);\n"
02105         "       v2[3] = vec4(ccc.r, ccu.r, crc.r, cru.r);\n"
02106 
02107         //test against 8 neighbours
02108         //use variable to identify type of extremum
02109         //1.0 for local maximum and -1.0 for minimum
02110         <<
02111         "       vec4 key = vec4(0.0); \n"
02112         "       #define KEYTEST_STEP0(i) \\\n"
02113         "       {\\\n"
02114         "               bvec4 test1 = greaterThan(vec4(cc[i]), max(v1[i], v2[i])), test2 = lessThan(vec4(cc[i]), min(v1[i], v2[i]));\\\n"
02115         "               key[i] = cc[i] > THRESHOLD0 && all(test1)?1.0: 0.0;\\\n"
02116         "               key[i] = cc[i] < -THRESHOLD0 && all(test2)? -1.0: key[i];\\\n"
02117         "       }\n"
02118         "       REPEAT4(KEYTEST_STEP0);\n"
02119         "       if(gl_TexCoord[0].x < 1.0) {key.rb = vec2(0.0);}\n"
02120         "       if(gl_TexCoord[0].y < 1.0) {key.rg = vec2(0.0);}\n"
02121         "       gl_FragColor = vec4(0.0);\n"
02122         "       if(any(notEqual(key, vec4(0.0)))) {\n";
02123 
02124         //do edge supression first.. 
02125         //vector v1 is < (-1, 0), (1, 0), (0,-1), (0, 1)>
02126         //vector v2 is < (-1,-1), (-1,1), (1,-1), (1, 1)>
02127 
02128         out<<
02129         "       float fxx[4], fyy[4], fxy[4], fx[4], fy[4];\n"
02130         "       #define EDGE_SUPPRESION(i) \\\n"
02131         "       if(key[i] != 0)\\\n"
02132         "       {\\\n"
02133         "               vec4 D2 = v1[i].xyzw - cc[i];\\\n"
02134         "               vec2 D4 = v2[i].xw - v2[i].yz;\\\n"
02135         "               vec2 D5 = 0.5*(v1[i].yw-v1[i].xz); \\\n"
02136         "               fx[i] = D5.x;   fy[i] = D5.y ;\\\n"
02137         "               fxx[i] = D2.x + D2.y;\\\n"
02138         "               fyy[i] = D2.z + D2.w;\\\n"
02139         "               fxy[i] = 0.25*(D4.x + D4.y);\\\n"
02140         "               float fxx_plus_fyy = fxx[i] + fyy[i];\\\n"
02141         "               float score_up = fxx_plus_fyy*fxx_plus_fyy; \\\n"
02142         "               float score_down = (fxx[i]*fyy[i] - fxy[i]*fxy[i]);\\\n"
02143         "               if( score_down <= 0 || score_up > THRESHOLD2 * score_down)key[i] = 0;\\\n"
02144         "       }\n"
02145         "       REPEAT4(EDGE_SUPPRESION);\n"
02146         "       if(any(notEqual(key, vec4(0.0)))) {\n";
02147         
02149         //read 9 pixels of upper/lower level
02150         out<<
02151         "       vec4  v4[4], v5[4], v6[4];\n"
02152         "       ccc = texture2DRect(texU, gl_TexCoord[0].xy);\n"
02153         "       clc = texture2DRect(texU, gl_TexCoord[1].xy);\n"
02154         "       crc = texture2DRect(texU, gl_TexCoord[2].xy);\n"
02155         "       ccd = texture2DRect(texU, gl_TexCoord[3].xy);\n"
02156         "       ccu = texture2DRect(texU, gl_TexCoord[4].xy);\n"
02157         "       cld = texture2DRect(texU, gl_TexCoord[5].xy);\n"
02158         "       clu = texture2DRect(texU, gl_TexCoord[6].xy);\n"
02159         "       crd = texture2DRect(texU, gl_TexCoord[7].xy);\n"
02160         "       cru = texture2DRect(texU, TexRU.xy);\n"
02161         "       vec4 cu = ccc;\n"
02162         "       v4[0] = vec4(clc.g, ccc.g, ccd.b, ccc.b);\n"
02163         "       v4[1] = vec4(ccc.r, crc.r, ccd.a, ccc.a);\n"
02164         "       v4[2] = vec4(clc.a, ccc.a, ccc.r, ccu.r);\n"
02165         "       v4[3] = vec4(ccc.b, crc.b, ccc.g, ccu.g);\n"
02166         "       v6[0] = vec4(cld.a, clc.a, ccd.a, ccc.a);\n"
02167         "       v6[1] = vec4(ccd.b, ccc.b, crd.b, crc.b);\n"
02168         "       v6[2] = vec4(clc.g, clu.g, ccc.g, ccu.g);\n"
02169         "       v6[3] = vec4(ccc.r, ccu.r, crc.r, cru.r);\n"
02170         <<
02171         "       #define KEYTEST_STEP1(i)\\\n"
02172         "       if(key[i] == 1.0)\\\n"
02173         "       {\\\n"
02174         "               bvec4 test = lessThan(vec4(cc[i]), max(v4[i], v6[i])); \\\n"
02175         "               if(cc[i] < cu[i] || any(test))key[i] = 0.0; \\\n"
02176         "       }else if(key[i] == -1.0)\\\n"
02177         "       {\\\n"
02178         "               bvec4 test = greaterThan(vec4(cc[i]), min(v4[i], v6[i])); \\\n"
02179         "               if(cc[i] > cu[i] || any(test) )key[i] = 0.0; \\\n"
02180         "       }\n"
02181         "       REPEAT4(KEYTEST_STEP1);\n"
02182         "       if(any(notEqual(key, vec4(0.0)))) { \n"
02183         <<
02184         "       ccc = texture2DRect(texD, gl_TexCoord[0].xy);\n"
02185         "       clc = texture2DRect(texD, gl_TexCoord[1].xy);\n"
02186         "       crc = texture2DRect(texD, gl_TexCoord[2].xy);\n"
02187         "       ccd = texture2DRect(texD, gl_TexCoord[3].xy);\n"
02188         "       ccu = texture2DRect(texD, gl_TexCoord[4].xy);\n"
02189         "       cld = texture2DRect(texD, gl_TexCoord[5].xy);\n"
02190         "       clu = texture2DRect(texD, gl_TexCoord[6].xy);\n"
02191         "       crd = texture2DRect(texD, gl_TexCoord[7].xy);\n"
02192         "       cru = texture2DRect(texD, TexRU.xy);\n"
02193         "       vec4 cd = ccc;\n"
02194         "       v5[0] = vec4(clc.g, ccc.g, ccd.b, ccc.b);\n"
02195         "       v5[1] = vec4(ccc.r, crc.r, ccd.a, ccc.a);\n"
02196         "       v5[2] = vec4(clc.a, ccc.a, ccc.r, ccu.r);\n"
02197         "       v5[3] = vec4(ccc.b, crc.b, ccc.g, ccu.g);\n"
02198         "       v6[0] = vec4(cld.a, clc.a, ccd.a, ccc.a);\n"
02199         "       v6[1] = vec4(ccd.b, ccc.b, crd.b, crc.b);\n"
02200         "       v6[2] = vec4(clc.g, clu.g, ccc.g, ccu.g);\n"
02201         "       v6[3] = vec4(ccc.r, ccu.r, crc.r, cru.r);\n"
02202         <<
02203         "       #define KEYTEST_STEP2(i)\\\n"
02204         "       if(key[i] == 1.0)\\\n"
02205         "       {\\\n"
02206         "               bvec4 test = lessThan(vec4(cc[i]), max(v5[i], v6[i]));\\\n"
02207         "               if(cc[i] < cd[i] || any(test))key[i] = 0.0; \\\n"
02208         "       }else if(key[i] == -1.0)\\\n"
02209         "       {\\\n"
02210         "               bvec4 test = greaterThan(vec4(cc[i]), min(v5[i], v6[i]));\\\n"
02211         "               if(cc[i] > cd[i] || any(test))key[i] = 0.0; \\\n"
02212         "       }\n"
02213         "       REPEAT4(KEYTEST_STEP2);\n"
02214         "       float keysum = dot(abs(key), vec4(1, 1, 1, 1)) ;\n"
02215         "       //assume there is only one keypoint in the four. \n"
02216         "       if(keysum==1.0) {\n";
02217 
02219         if(GlobalUtil::_SubpixelLocalization)
02220 
02221         out <<
02222         "       vec3 offset = vec3(0, 0, 0); \n"
02223         "       #define TESTMOVE_KEYPOINT(idx) \\\n"
02224         "       if(key[idx] != 0) \\\n"
02225         "       {\\\n"
02226         "               cu[0] = cu[idx];        cd[0] = cd[idx];        cc[0] = cc[idx];        \\\n"
02227         "               v4[0] = v4[idx];        v5[0] = v5[idx];                                                \\\n"
02228         "               fxy[0] = fxy[idx];      fxx[0] = fxx[idx];      fyy[0] = fyy[idx];      \\\n"
02229         "               fx[0] = fx[idx];        fy[0] = fy[idx];                                                \\\n"
02230         "       }\n"
02231         "       TESTMOVE_KEYPOINT(1);\n"
02232         "       TESTMOVE_KEYPOINT(2);\n"
02233         "       TESTMOVE_KEYPOINT(3);\n"
02234         <<
02235                 
02236         "       float fs = 0.5*( cu[0] - cd[0] );                               \n"
02237         "       float fss = cu[0] + cd[0] - cc[0] - cc[0];\n"
02238         "       float fxs = 0.25 * (v4[0].y + v5[0].x - v4[0].x - v5[0].y);\n"
02239         "       float fys = 0.25 * (v4[0].w + v5[0].z - v4[0].z - v5[0].w);\n"
02240         "       vec4 A0, A1, A2 ;                       \n"
02241         "       A0 = vec4(fxx[0], fxy[0], fxs, -fx[0]); \n"
02242         "       A1 = vec4(fxy[0], fyy[0], fys, -fy[0]); \n"
02243         "       A2 = vec4(fxs, fys, fss, -fs);  \n"
02244         "       vec3 x3 = abs(vec3(fxx[0], fxy[0], fxs));               \n"
02245         "       float maxa = max(max(x3.x, x3.y), x3.z);        \n"
02246         "       if(maxa >= 1e-10 ) \n"
02247         "       {                                                                                               \n"
02248         "               if(x3.y ==maxa )                                                        \n"
02249         "               {                                                                                       \n"
02250         "                       vec4 TEMP = A1; A1 = A0; A0 = TEMP;     \n"
02251         "               }else if( x3.z == maxa )                                        \n"
02252         "               {                                                                                       \n"
02253         "                       vec4 TEMP = A2; A2 = A0; A0 = TEMP;     \n"
02254         "               }                                                                                       \n"
02255         "               A0 /= A0.x;                                                                     \n"
02256         "               A1 -= A1.x * A0;                                                        \n"
02257         "               A2 -= A2.x * A0;                                                        \n"
02258         "               vec2 x2 = abs(vec2(A1.y, A2.y));                \n"
02259         "               if( x2.y > x2.x )                                                       \n"
02260         "               {                                                                                       \n"
02261         "                       vec3 TEMP = A2.yzw;                                     \n"
02262         "                       A2.yzw = A1.yzw;                                                \n"
02263         "                       A1.yzw = TEMP;                                                  \n"
02264         "                       x2.x = x2.y;                                                    \n"
02265         "               }                                                                                       \n"
02266         "               if(x2.x >= 1e-10) {                                                             \n"
02267         "                       A1.yzw /= A1.y;                                                         \n"
02268         "                       A2.yzw -= A2.y * A1.yzw;                                        \n"
02269         "                       if(abs(A2.z) >= 1e-10) {\n"
02270         "                               offset.z = A2.w /A2.z;                              \n"
02271         "                               offset.y = A1.w - offset.z*A1.z;                            \n"
02272         "                               offset.x = A0.w - offset.z*A0.z - offset.y*A0.y;        \n"
02273         "                               bool test = (abs(cc[0] + 0.5*dot(vec3(fx[0], fy[0], fs), offset ))>THRESHOLD1) ;\n"
02274         "                               if(!test || any( greaterThan(abs(offset), vec3(1.0)))) key = vec4(0.0);\n"
02275         "                       }\n"
02276         "               }\n"
02277         "       }\n"
02278         <<"\n"
02279         "       float keyv = dot(key, vec4(1.0, 2.0, 3.0, 4.0));\n"
02280         "       gl_FragColor = vec4(keyv,  offset);\n"
02281         "       }}}}\n"
02282         "}\n"   <<'\0';
02283 
02284         else out << "\n"
02285         "       float keyv = dot(key, vec4(1.0, 2.0, 3.0, 4.0));\n"
02286         "       gl_FragColor =  vec4(keyv, 0, 0, 0);\n"
02287         "       }}}}\n"
02288         "}\n"   <<'\0';
02289 
02290         ProgramGLSL * program = new ProgramGLSL(buffer); 
02291         s_keypoint = program ;
02292 
02293         //parameter
02294         _param_dog_texu = glGetUniformLocation(*program, "texU");
02295         _param_dog_texd = glGetUniformLocation(*program, "texD");
02296 }
02297 void ShaderBagPKSL::SetDogTexParam(int texU, int texD)
02298 {
02299         glUniform1i(_param_dog_texu, 1);
02300         glUniform1i(_param_dog_texd, 2);
02301 }
02302 void ShaderBagPKSL::SetGenListStepParam(int tex, int tex0)
02303 {
02304         glUniform1i(_param_genlist_step_tex0, 1);       
02305 }
02306 
02307 void ShaderBagPKSL::SetGenVBOParam(float width, float fwidth,float size)
02308 {
02309         float sizes[4] = {size*3.0f, fwidth, width, 1.0f/width};
02310         glUniform4fv(_param_genvbo_size, 1, sizes);
02311 }
02312 void ShaderBagPKSL::SetGradPassParam(int texP)
02313 {
02314         glUniform1i(_param_grad_pass_texp, 1);
02315 }
02316 
02317 void ShaderBagPKSL::LoadDescriptorShader()
02318 {
02319         GlobalUtil::_DescriptorPPT = 16;
02320         LoadDescriptorShaderF2();
02321     s_rect_description = LoadDescriptorProgramRECT();
02322 }
02323 
02324 ProgramGLSL* ShaderBagPKSL::LoadDescriptorProgramRECT()
02325 {
02326         //one shader outpout 128/8 = 16 , each fragout encodes 4
02327         //const double twopi = 2.0*3.14159265358979323846;
02328         //const double rpi  = 8.0/twopi;
02329         char buffer[10240];
02330         ostrstream out(buffer, 10240);
02331 
02332         out<<setprecision(8);
02333 
02334         out<<"\n"
02335         "#define M_PI 3.14159265358979323846\n"
02336         "#define TWO_PI (2.0*M_PI)\n"
02337         "#define RPI 1.2732395447351626861510701069801\n"
02338         "#define WF size.z\n"
02339         "uniform sampler2DRect tex;                     \n"
02340         "uniform sampler2DRect gtex;                    \n"
02341         "uniform sampler2DRect otex;                    \n"
02342         "uniform vec4           dsize;                          \n"
02343         "uniform vec3           size;                           \n"
02344         "void main()                    \n"
02345         "{\n"
02346         "       vec2 dim        = size.xy;      //image size                    \n"
02347         "       float index = dsize.x*floor(gl_TexCoord[0].y * 0.5) + gl_TexCoord[0].x;\n"
02348         "       float idx = 8.0* fract(index * 0.125) + 8.0 * floor(2.0* fract(gl_TexCoord[0].y * 0.5));                \n"
02349         "       index = floor(index*0.125)+ 0.49;  \n"
02350         "       vec2 coord = floor( vec2( mod(index, dsize.z), index*dsize.w)) + 0.5 ;\n"
02351         "       vec2 pos = texture2DRect(tex, coord).xy;                \n"
02352         "       vec2 wsz = texture2DRect(tex, coord).zw;\n"
02353     "   float aspect_ratio = wsz.y / wsz.x;\n"
02354     "   float aspect_sq = aspect_ratio * aspect_ratio; \n"
02355         "       vec2 spt  = wsz * 0.25; vec2 ispt = 1.0 / spt; \n";
02356 
02357         //here cscs is actually (cos, sin, -cos, -sin) * (factor: 3)*sigma
02358         //and rots is  (cos, sin, -cos, -sin ) /(factor*sigma)
02359         //devide the 4x4 sift grid into 16 1x1 block, and each corresponds to a shader thread
02360         //To use linear interoplation, 1x1 is increased to 2x2, by adding 0.5 to each side
02361         out<<
02362         "       vec4 temp; vec2 pt;                             \n"
02363     "   pt.x = pos.x + fract(idx*0.25) * wsz.x;                         \n"
02364         "       pt.y = pos.y + (floor(idx*0.25) + 0.5) * spt.y;                 \n";
02365         
02366         //get a horizontal bounding box of the rotated rectangle
02367         out<<
02368     "   vec4 sz;                                        \n"
02369         "       sz.xy = max(pt - spt, vec2(2,2));\n"
02370         "       sz.zw = min(pt + spt, dim - 3);         \n"
02371         "       sz = floor(sz * 0.5)+0.5;"; //move sample point to pixel center
02372         //get voting for two box
02373 
02374         out<<"\n"
02375         "       vec4 DA, DB;   vec2 spos;                       \n"
02376         "       DA = DB  = vec4(0, 0, 0, 0);            \n"
02377         "       vec4 nox = vec4(0, 1.0, 0.0, 1.0);                                      \n"
02378         "       vec4 noy = vec4(0, 0.0, 1.0, 1.0);                                      \n"
02379         "       for(spos.y = sz.y; spos.y <= sz.w;      spos.y+=1.0)                            \n"
02380         "       {                                                                                                                               \n"
02381         "               for(spos.x = sz.x; spos.x <= sz.z;      spos.x+=1.0)                    \n"
02382         "               {                                                                                                                       \n"
02383         "                       vec2 tpt = spos * 2.0 - pt - 0.5;                                       \n"
02384     "                   vec4 nx = (tpt.x + nox) * ispt.x;                                                               \n"
02385     "                   vec4 ny = (tpt.y + noy) * ispt.y;                       \n"
02386         "                       vec4 nxn = abs(nx), nyn = abs(ny);                                              \n"
02387     "                   bvec4 inside = lessThan(max(nxn, nyn) , vec4(1.0));     \n"
02388         "                       if(any(inside))\n"
02389         "                       {\n"
02390         "                               vec4 gg = texture2DRect(gtex, spos);\n"
02391         "                               vec4 oo = texture2DRect(otex, spos);\n"
02392     //"               vec4 cc = cos(oo), ss = sin(oo); \n"
02393     //"               oo = atan(ss* aspect_ratio, cc); \n"
02394     //"               gg = gg * sqrt(ss * ss * aspect_sq + cc * cc); \n "
02395         "                               vec4 theta0 = (- oo)*RPI;\n"
02396         "                               vec4 theta = 8.0 * fract(1.0 + 0.125 * theta0);                 \n"
02397         "                               vec4 theta1 = floor(theta);                                                             \n"
02398         "                               vec4 weight = (1 - nxn) * (1 - nyn) * gg; \n"
02399         "                               vec4 weight2 = (theta - theta1) * weight;                               \n"
02400         "                               vec4 weight1 = weight - weight2;                                                \n"
02401         "                               for(int i = 0;i < 4; i++)\n"
02402         "                               {\n"
02403         "                                       if(inside[i])\n"
02404         "                                       {\n"
02405         "                                               DA += vec4(equal(vec4(theta1[i]), vec4(0, 1, 2, 3)))*weight1[i]; \n"
02406         "                                               DA += vec4(equal(vec4(theta1[i]), vec4(7, 0, 1, 2)))*weight2[i]; \n"
02407         "                                               DB += vec4(equal(vec4(theta1[i]), vec4(4, 5, 6, 7)))*weight1[i]; \n"
02408         "                                               DB += vec4(equal(vec4(theta1[i]), vec4(3, 4, 5, 6)))*weight2[i]; \n"
02409         "                                       }\n"
02410         "                               }\n"
02411         "                       }\n"
02412         "               }\n"
02413         "       }\n";
02414         out<<
02415         "        gl_FragData[0] = DA; gl_FragData[1] = DB;\n"
02416         "}\n"<<'\0';
02417 
02418         ProgramGLSL * program =  new ProgramGLSL(buffer); 
02419         if(program->IsNative()) 
02420         {
02421                 return program;
02422         }
02423         else
02424         {
02425                 delete program;
02426                 return NULL;
02427         }
02428 }
02429 
02430 ProgramGLSL* ShaderBagPKSL::LoadDescriptorProgramPKSL()
02431 {
02432         //one shader outpout 128/8 = 16 , each fragout encodes 4
02433         //const double twopi = 2.0*3.14159265358979323846;
02434         //const double rpi  = 8.0/twopi;
02435         char buffer[10240];
02436         ostrstream out(buffer, 10240);
02437 
02438         out<<setprecision(8);
02439 
02440         out<<"\n"
02441         "#define M_PI 3.14159265358979323846\n"
02442         "#define TWO_PI (2.0*M_PI)\n"
02443         "#define RPI 1.2732395447351626861510701069801\n"
02444         "#define WF size.z\n"
02445         "uniform sampler2DRect tex;                     \n"
02446         "uniform sampler2DRect gtex;                    \n"
02447         "uniform sampler2DRect otex;                    \n"
02448         "uniform vec4           dsize;                          \n"
02449         "uniform vec3           size;                           \n"
02450         "void main()                    \n"
02451         "{\n"
02452         "       vec2 dim        = size.xy;      //image size                    \n"
02453         "       float index = dsize.x*floor(gl_TexCoord[0].y * 0.5) + gl_TexCoord[0].x;\n"
02454         "       float idx = 8.0* fract(index * 0.125) + 8.0 * floor(2.0* fract(gl_TexCoord[0].y * 0.5));                \n"
02455         "       index = floor(index*0.125)+ 0.49;  \n"
02456         "       vec2 coord = floor( vec2( mod(index, dsize.z), index*dsize.w)) + 0.5 ;\n"
02457         "       vec2 pos = texture2DRect(tex, coord).xy;                \n"
02458         "       if(any(lessThan(pos.xy, vec2(1.0))) || any(greaterThan(pos.xy, dim-1.0))) "
02459         "       //discard;      \n" 
02460         "       { gl_FragData[0] = gl_FragData[1] = vec4(0.0); return; }\n"
02461         "       float anglef = texture2DRect(tex, coord).z;\n"
02462         "       if(anglef > M_PI) anglef -= TWO_PI;\n"
02463         "       float sigma = texture2DRect(tex, coord).w; \n"
02464         "       float spt  = abs(sigma * WF);   //default to be 3*sigma \n";
02465         //rotation
02466         out<<
02467         "       vec4 cscs, rots;                                                \n"
02468         "       cscs.x = cos(anglef); cscs.y = sin(anglef);     \n"
02469         "       cscs.zw = - cscs.xy;                                                    \n"
02470         "       rots = cscs /spt;                                                               \n"
02471         "       cscs *= spt; \n";
02472 
02473         //here cscs is actually (cos, sin, -cos, -sin) * (factor: 3)*sigma
02474         //and rots is  (cos, sin, -cos, -sin ) /(factor*sigma)
02475         //devide the 4x4 sift grid into 16 1x1 block, and each corresponds to a shader thread
02476         //To use linear interoplation, 1x1 is increased to 2x2, by adding 0.5 to each side
02477         out<<
02478         "       vec4 temp; vec2 pt, offsetpt;                           \n"
02479         "       /*the fraction part of idx is .5*/                      \n"
02480         "       offsetpt.x = 4.0* fract(idx*0.25) - 2.0;                                \n"
02481         "       offsetpt.y = floor(idx*0.25) - 1.5;                     \n"
02482         "       temp = cscs.xwyx*offsetpt.xyxy;                         \n"
02483         "       pt = pos + temp.xz + temp.yw;                           \n";
02484         
02485         //get a horizontal bounding box of the rotated rectangle
02486         out<<
02487         "       vec2 bwin = abs(cscs.xy);                                       \n"
02488         "       float bsz = bwin.x + bwin.y;                                    \n"
02489         "       vec4 sz;                                        \n"
02490         "       sz.xy = max(pt - vec2(bsz), vec2(2,2));\n"
02491         "       sz.zw = min(pt + vec2(bsz), dim - 3);           \n"
02492         "       sz = floor(sz * 0.5)+0.5;"; //move sample point to pixel center
02493         //get voting for two box
02494 
02495         out<<"\n"
02496         "       vec4 DA, DB;   vec2 spos;                       \n"
02497         "       DA = DB  = vec4(0, 0, 0, 0);            \n"
02498         "       vec4 nox = vec4(0, rots.xy, rots.x + rots.y);                                   \n"
02499         "       vec4 noy = vec4(0, rots.wx, rots.w + rots.x);                                   \n"
02500         "       for(spos.y = sz.y; spos.y <= sz.w;      spos.y+=1.0)                            \n"
02501         "       {                                                                                                                               \n"
02502         "               for(spos.x = sz.x; spos.x <= sz.z;      spos.x+=1.0)                    \n"
02503         "               {                                                                                                                       \n"
02504         "                       vec2 tpt = spos * 2.0 - pt - 0.5;                                       \n"
02505         "                       vec4 temp = rots.xywx * tpt.xyxy;                                               \n"
02506         "                       vec2 temp2 = temp.xz + temp.yw;                                         \n"
02507         "                       vec4 nx = temp2.x + nox;                                                                \n"
02508         "                       vec4 ny = temp2.y + noy;                        \n"
02509         "                       vec4 nxn = abs(nx), nyn = abs(ny);                                              \n"
02510         "                       bvec4 inside = lessThan(max(nxn, nyn) , vec4(1.0));     \n"
02511         "                       if(any(inside))\n"
02512         "                       {\n"
02513         "                               vec4 gg = texture2DRect(gtex, spos);\n"
02514         "                               vec4 oo = texture2DRect(otex, spos);\n"
02515         "                               vec4 theta0 = (anglef - oo)*RPI;\n"
02516         "                               vec4 theta = 8.0 * fract(1.0 + 0.125 * theta0);                 \n"
02517         "                               vec4 theta1 = floor(theta);                                                             \n"
02518         "                               vec4 diffx = nx + offsetpt.x, diffy = ny + offsetpt.y;  \n"
02519         "                               vec4 ww = exp(-0.125 * (diffx * diffx + diffy * diffy ));       \n"
02520         "                               vec4 weight = (1 - nxn) * (1 - nyn) * gg * ww; \n"
02521         "                               vec4 weight2 = (theta - theta1) * weight;                               \n"
02522         "                               vec4 weight1 = weight - weight2;                                                \n"
02523         "                               for(int i = 0;i < 4; i++)\n"
02524         "                               {\n"
02525         "                                       if(inside[i])\n"
02526         "                                       {\n"
02527         "                                               DA += vec4(equal(vec4(theta1[i]), vec4(0, 1, 2, 3)))*weight1[i]; \n"
02528         "                                               DA += vec4(equal(vec4(theta1[i]), vec4(7, 0, 1, 2)))*weight2[i]; \n"
02529         "                                               DB += vec4(equal(vec4(theta1[i]), vec4(4, 5, 6, 7)))*weight1[i]; \n"
02530         "                                               DB += vec4(equal(vec4(theta1[i]), vec4(3, 4, 5, 6)))*weight2[i]; \n"
02531         "                                       }\n"
02532         "                               }\n"
02533         "                       }\n"
02534         "               }\n"
02535         "       }\n";
02536         out<<
02537         "        gl_FragData[0] = DA; gl_FragData[1] = DB;\n"
02538         "}\n"<<'\0';
02539 
02540         ProgramGLSL * program =  new ProgramGLSL(buffer); 
02541         if(program->IsNative()) 
02542         {
02543                 return program;
02544         }
02545         else
02546         {
02547                 delete program;
02548                 return NULL;
02549         }
02550 }
02551 
02552 void ShaderBagPKSL::LoadDescriptorShaderF2()
02553 {
02554 
02555         ProgramGLSL * program = LoadDescriptorProgramPKSL();
02556         if( program )
02557         {
02558                 s_descriptor_fp = program;
02559                 _param_descriptor_gtex = glGetUniformLocation(*program, "gtex");
02560                 _param_descriptor_otex = glGetUniformLocation(*program, "otex");
02561                 _param_descriptor_size = glGetUniformLocation(*program, "size");
02562                 _param_descriptor_dsize = glGetUniformLocation(*program, "dsize");
02563         }
02564 }
02565 
02566 
02567 
02568 void ShaderBagPKSL::SetSimpleOrientationInput(int oTex, float sigma, float sigma_step)
02569 {
02570         glUniform1i(_param_orientation_gtex, 1);
02571         glUniform2f(_param_orientation_size, sigma, sigma_step);
02572 }
02573 
02574 
02575 void ShaderBagPKSL::SetFeatureOrientationParam(int gtex, int width, int height, float sigma, int otex, float step)
02576 {
02578         glUniform1i(_param_orientation_gtex, 1);        
02579         glUniform1i(_param_orientation_otex, 2);        
02580 
02581         float size[4];
02582         size[0] = (float)width;
02583         size[1] = (float)height;
02584         size[2] = sigma;
02585         size[3] = step;
02586         glUniform4fv(_param_orientation_size, 1, size);
02587 
02588 }
02589 
02590 void ShaderBagPKSL::SetFeatureDescirptorParam(int gtex, int otex, float dwidth, float fwidth,  float width, float height, float sigma)
02591 {
02592     if(sigma == 0 && s_rect_description)
02593     {
02594         //rectangle description mode
02595         s_rect_description->UseProgram();       
02596         GLint param_descriptor_gtex = glGetUniformLocation(*s_rect_description, "gtex");
02597                 GLint param_descriptor_otex = glGetUniformLocation(*s_rect_description, "otex");
02598                 GLint param_descriptor_size = glGetUniformLocation(*s_rect_description, "size");
02599                 GLint param_descriptor_dsize = glGetUniformLocation(*s_rect_description, "dsize");
02601             glUniform1i(param_descriptor_gtex, 1);      
02602             glUniform1i(param_descriptor_otex, 2);      
02603 
02604             float dsize[4] ={dwidth, 1.0f/dwidth, fwidth, 1.0f/fwidth};
02605             glUniform4fv(param_descriptor_dsize, 1, dsize);
02606             float size[3];
02607             size[0] = width;
02608             size[1] = height;
02609             size[2] = GlobalUtil::_DescriptorWindowFactor;
02610             glUniform3fv(param_descriptor_size, 1, size);
02611     }else
02612     {
02614             glUniform1i(_param_descriptor_gtex, 1);     
02615             glUniform1i(_param_descriptor_otex, 2);     
02616 
02617 
02618             float dsize[4] ={dwidth, 1.0f/dwidth, fwidth, 1.0f/fwidth};
02619             glUniform4fv(_param_descriptor_dsize, 1, dsize);
02620             float size[3];
02621             size[0] = width;
02622             size[1] = height;
02623             size[2] = GlobalUtil::_DescriptorWindowFactor;
02624             glUniform3fv(_param_descriptor_size, 1, size);
02625     }
02626 
02627 }
02628 
02629 
02630 void ShaderBagPKSL::SetGenListEndParam(int ktex)
02631 {
02632         glUniform1i(_param_genlist_end_ktex, 1);
02633 }
02634 void ShaderBagPKSL::SetGenListInitParam(int w, int h)
02635 {
02636         float bbox[4] = {(w -1.0f) * 0.5f +0.25f, (w-1.0f) * 0.5f - 0.25f,  (h - 1.0f) * 0.5f + 0.25f, (h-1.0f) * 0.5f - 0.25f};
02637         glUniform4fv(_param_genlist_init_bbox, 1, bbox);
02638 }
02639 
02640 void ShaderBagPKSL::SetMarginCopyParam(int xmax, int ymax)
02641 {
02642         float truncate[4];
02643         truncate[0] = (xmax - 0.5f) * 0.5f; //((xmax + 1)  >> 1) - 0.5f;
02644         truncate[1] = (ymax - 0.5f) * 0.5f; //((ymax + 1)  >> 1) - 0.5f;
02645         truncate[2] = (xmax %2 == 1)? 0.0f: 1.0f;
02646         truncate[3] = truncate[2] +  (((ymax % 2) == 1)? 0.0f : 2.0f);
02647         glUniform4fv(_param_margin_copy_truncate, 1,  truncate);
02648 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Defines


rgbd_registration
Author(s): Ross Kidson
autogenerated on Thu May 23 2013 15:36:53