text: /*#******************************************************************************
  ** IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
  **
  ** By downloading, copying, installing or using the software you agree to this license.
  ** If you do not agree to this license, do not download, install,
  ** copy or use the software.
  **
  **
  ** bioinspired : interfaces allowing OpenCV users to integrate Human Vision System models. Presented models originate from Jeanny Herault's original research and have been reused and adapted by the author&collaborators for computed vision applications since his thesis with Alice Caplier at Gipsa-Lab.
  ** Use: extract still images & image sequences features, from contours details to motion spatio-temporal features, etc. for high level visual scene analysis. Also contribute to image enhancement/compression such as tone mapping.
  **
  ** Maintainers : Listic lab (code author current affiliation & applications) and Gipsa Lab (original research origins & applications)
  **
  **  Creation - enhancement process 2007-2011
  **      Author: Alexandre Benoit (benoit.alexandre.vision@gmail.com), LISTIC lab, Annecy le vieux, France
  **
  ** Theses algorithm have been developped by Alexandre BENOIT since his thesis with Alice Caplier at Gipsa-Lab (www.gipsa-lab.inpg.fr) and the research he pursues at LISTIC Lab (www.listic.univ-savoie.fr).
  ** Refer to the following research paper for more information:
  ** Benoit A., Caplier A., Durette B., Herault, J., "USING HUMAN VISUAL SYSTEM MODELING FOR BIO-INSPIRED LOW LEVEL IMAGE PROCESSING", Elsevier, Computer Vision and Image Understanding 114 (2010), pp. 758-773, DOI: http://dx.doi.org/10.1016/j.cviu.2010.01.011
  ** This work have been carried out thanks to Jeanny Herault who's research and great discussions are the basis of all this work, please take a look at his book:
  ** Vision: Images, Signals and Neural Networks: Models of Neural Processing in Visual Perception (Progress in Neural Processing),By: Jeanny Herault, ISBN: 9814273686. WAPI (Tower ID): 113266891.
  **
  ** The retina filter includes the research contributions of phd/research collegues from which code has been redrawn by the author :
  ** _take a look at the retinacolor.hpp module to discover Brice Chaix de Lavarene color mosaicing/demosaicing and the reference paper:
  ** ====> B. Chaix de Lavarene, D. Alleysson, B. Durette, J. Herault (2007). "Efficient demosaicing through recursive filtering", IEEE International Conference on Image Processing ICIP 2007
  ** _take a look at imagelogpolprojection.hpp to discover retina spatial log sampling which originates from Barthelemy Durette phd with Jeanny Herault. A Retina / V1 cortex projection is also proposed and originates from Jeanny's discussions.
  ** ====> more informations in the above cited Jeanny Heraults's book.
  **
  **                          License Agreement
  **               For Open Source Computer Vision Library
  **
  ** Copyright (C) 2000-2008, Intel Corporation, all rights reserved.
  ** Copyright (C) 2008-2011, Willow Garage Inc., all rights reserved.
  **
  **               For Human Visual System tools (bioinspired)
  ** Copyright (C) 2007-2011, LISTIC Lab, Annecy le Vieux and GIPSA Lab, Grenoble, France, all rights reserved.
  **
  ** Third party copyrights are property of their respective owners.
  **
  ** Redistribution and use in source and binary forms, with or without modification,
  ** are permitted provided that the following conditions are met:
  **
  ** * Redistributions of source code must retain the above copyright notice,
  **    this list of conditions and the following disclaimer.
  **
  ** * Redistributions in binary form must reproduce the above copyright notice,
  **    this list of conditions and the following disclaimer in the documentation
  **    and/or other materials provided with the distribution.
  **
  ** * The name of the copyright holders may not be used to endorse or promote products
  **    derived from this software without specific prior written permission.
  **
  ** This software is provided by the copyright holders and contributors "as is" and
  ** any express or implied warranties, including, but not limited to, the implied
  ** warranties of merchantability and fitness for a particular purpose are disclaimed.
  ** In no event shall the Intel Corporation or contributors be liable for any direct,
  ** indirect, incidental, special, exemplary, or consequential damages
  ** (including, but not limited to, procurement of substitute goods or services;
  ** loss of use, data, or profits; or business interruption) however caused
  ** and on any theory of liability, whether in contract, strict liability,
  ** or tort (including negligence or otherwise) arising in any way out of
  ** the use of this software, even if advised of the possibility of such damage.
  *******************************************************************************/
  
  /**
  * @class RetinaColor a color multilexing/demultiplexing (demosaicing) based on a human vision inspiration. Different mosaicing strategies can be used, included random sampling !
  * => please take a look at the nice and efficient demosaicing strategy introduced by B.Chaix de Lavarene, take a look at the cited paper for more mathematical details
  * @brief Retina color sampling model which allows classical bayer sampling, random and potentially several other method ! Low color errors on corners !
  * -> Based on the research of:
  *		.Brice Chaix Lavarene (chaix@lis.inpg.fr)
  *		.Jeanny Herault (herault@lis.inpg.fr)
  *		.David Alleyson (david.alleyson@upmf-grenoble.fr)
  *      .collaboration: alexandre benoit (benoit.alexandre.vision@gmail.com or benoit@lis.inpg.fr)
  * Please cite: B. Chaix de Lavarene, D. Alleysson, B. Durette, J. Herault (2007). "Efficient demosaicing through recursive filtering", IEEE International Conference on Image Processing ICIP 2007
  * @author Alexandre BENOIT, benoit.alexandre.vision@gmail.com, LISTIC / Gipsa-Lab, France: www.gipsa-lab.inpg.fr/
  * Creation date 2007
  */
  
  #ifndef RETINACOLOR_HPP_
  #define RETINACOLOR_HPP_
  
  #include "basicretinafilter.hpp"
  
  //#define __RETINACOLORDEBUG //define RETINACOLORDEBUG in order to display debug data
  
  namespace cv
  {
  namespace bioinspired
  {
      class RetinaColor: public BasicRetinaFilter
      {
      public:
          /**
          * @typedef which allows to select the type of photoreceptors color sampling
          */
  
          /**
          * constructor of the retina color processing model
          * @param NBrows: number of rows of the input image
          * @param NBcolumns: number of columns of the input image
          * @param samplingMethod: the chosen color sampling method
          */
          RetinaColor(const unsigned int NBrows, const unsigned int NBcolumns, const int samplingMethod=RETINA_COLOR_BAYER);
  
          /**
          * standard destructor
          */
          virtual ~RetinaColor();
  
          /**
          * function that clears all buffers of the object
          */
          void clearAllBuffers();
  
          /**
          * resize retina color filter object (resize all allocated buffers)
          * @param NBrows: the new height size
          * @param NBcolumns: the new width size
          */
          void resize(const unsigned int NBrows, const unsigned int NBcolumns);
  
  
          /**
          * color multiplexing function: a demultiplexed RGB frame of size M*N*3 is transformed into a multiplexed M*N*1 pixels frame where each pixel is either Red, or Green or Blue
          * @param inputRGBFrame: the input RGB frame to be processed
          * @return, nothing but the multiplexed frame is available by the use of the getMultiplexedFrame() function
          */
          inline void runColorMultiplexing(const std::valarray<float> &inputRGBFrame) { runColorMultiplexing(inputRGBFrame, *_multiplexedFrame); }
  
          /**
          * color multiplexing function: a demultipleed RGB frame of size M*N*3 is transformed into a multiplexed M*N*1 pixels frame where each pixel is either Red, or Green or Blue if using RGB images
          * @param demultiplexedInputFrame: the demultiplexed input frame to be processed of size M*N*3
          * @param multiplexedFrame: the resulting multiplexed frame
          */
          void runColorMultiplexing(const std::valarray<float> &demultiplexedInputFrame, std::valarray<float> &multiplexedFrame);
  
          /**
          * color demultiplexing function: a multiplexed frame of size M*N*1 pixels is transformed into a RGB demultiplexed M*N*3 pixels frame
          * @param multiplexedColorFrame: the input multiplexed frame to be processed
          * @param adaptiveFiltering: specifies if an adaptive filtering has to be perform rather than standard filtering (adaptive filtering allows a better rendering)
          * @param maxInputValue: the maximum input data value (should be 255 for 8 bits images but it can change in the case of High Dynamic Range Images (HDRI)
          * @return, nothing but the output demultiplexed frame is available by the use of the getDemultiplexedColorFrame() function, also use getLuminance() and getChrominance() in order to retreive either luminance or chrominance
          */
          void runColorDemultiplexing(const std::valarray<float> &multiplexedColorFrame, const bool adaptiveFiltering=false, const float maxInputValue=255.0);
  
          /**
          * activate color saturation as the final step of the color demultiplexing process
          * -> this saturation is a sigmoide function applied to each channel of the demultiplexed image.
          * @param saturateColors: boolean that activates color saturation (if true) or desactivate (if false)
          * @param colorSaturationValue: the saturation factor
          * */
          void setColorSaturation(const bool saturateColors=true, const float colorSaturationValue=4.0) { _saturateColors=saturateColors; _colorSaturationValue=colorSaturationValue; }
  
          /**
          * set parameters of the low pass spatio-temporal filter used to retreive the low chrominance
          * @param beta: gain of the filter (generally set to zero)
          * @param tau: time constant of the filter (unit is frame for video processing), typically 0 when considering static processing, 1 or more if a temporal smoothing effect is required
          * @param k: spatial constant of the filter (unit is pixels), typical value is 2.5
          */
          void setChrominanceLPfilterParameters(const float beta, const float tau, const float k) { setLPfilterParameters(beta, tau, k); }
  
          /**
          * apply to the retina color output the Krauskopf transformation which leads to an opponent color system: output colorspace if Acr1cr2 if input of the retina was LMS color space
          * @param result: the input buffer to fill with the transformed colorspace retina output
          * @return true if process ended successfully
          */
          bool applyKrauskopfLMS2Acr1cr2Transform(std::valarray<float> &result);
  
          /**
          * apply to the retina color output the CIE Lab color transformation
          * @param result: the input buffer to fill with the transformed colorspace retina output
          * @return true if process ended successfully
          */
          bool applyLMS2LabTransform(std::valarray<float> &result);
  
          /**
          * @return the multiplexed frame result (use this after function runColorMultiplexing)
          */
          inline const std::valarray<float> &getMultiplexedFrame() const { return *_multiplexedFrame; }
  
          /**
          * @return the demultiplexed frame result (use this after function runColorDemultiplexing)
          */
          inline const std::valarray<float> &getDemultiplexedColorFrame() const { return _demultiplexedColorFrame; }
  
          /**
          * @return the luminance of the processed frame (use this after function runColorDemultiplexing)
          */
          inline const std::valarray<float> &getLuminance() const { return *_luminance; }
  
          /**
          * @return the chrominance of the processed frame (use this after function runColorDemultiplexing)
          */
          inline const std::valarray<float> &getChrominance() const { return _chrominance; }
  
          /**
          * standard 0 to 255 image clipping function appled to RGB images (of size M*N*3 pixels)
          * @param inputOutputBuffer: the image to be normalized (rewrites the input), if no parameter, then, the built in buffer reachable by getOutput() function is normalized
          * @param maxOutputValue: the maximum value allowed at the output (values superior to it would be clipped
          */
          void clipRGBOutput_0_maxInputValue(float *inputOutputBuffer, const float maxOutputValue=255.0);
  
          /**
          * standard 0 to 255 image normalization function appled to RGB images (of size M*N*3 pixels)
          * @param maxOutputValue: the maximum value allowed at the output (values superior to it would be clipped
          */
          void normalizeRGBOutput_0_maxOutputValue(const float maxOutputValue=255.0);
  
          /**
          * return the color sampling map: a Nrows*Mcolumns image in which each pixel value is the ofsset adress which gives the adress of the sampled pixel on an Nrows*Mcolumns*3 color image ordered by layers: layer1, layer2, layer3
          */
          inline const std::valarray<unsigned int> &getSamplingMap() const { return _colorSampling; }
  
          /**
          * function used (to bypass processing) to manually set the color output
          * @param demultiplexedImage: the color image (luminance+chrominance) which has to be written in the object buffer
          */
          inline void setDemultiplexedColorFrame(const std::valarray<float> &demultiplexedImage) { _demultiplexedColorFrame=demultiplexedImage; }
  
      protected:
  
          // private functions
          int _samplingMethod;
          bool _saturateColors;
          float _colorSaturationValue;
          // links to parent buffers (more convienient names
          TemplateBuffer<float> *_luminance;
          std::valarray<float> *_multiplexedFrame;
          // instance buffers
          std::valarray<unsigned int> _colorSampling; // table (size (_nbRows*_nbColumns) which specifies the color of each pixel
          std::valarray<float> _RGBmosaic;
          std::valarray<float> _tempMultiplexedFrame;
          std::valarray<float> _demultiplexedTempBuffer;
          std::valarray<float> _demultiplexedColorFrame;
          std::valarray<float> _chrominance;
          std::valarray<float> _colorLocalDensity;// buffer which contains the local density of the R, G and B photoreceptors for a normalization use
          std::valarray<float> _imageGradient;
  
          // variables
          float _pR, _pG, _pB; // probabilities of color R, G and B
          bool _objectInit;
  
          // protected functions
          void _initColorSampling();
          void _interpolateImageDemultiplexedImage(float *inputOutputBuffer);
          void _interpolateSingleChannelImage111(float *inputOutputBuffer);
          void _interpolateBayerRGBchannels(float *inputOutputBuffer);
          void _applyRIFfilter(const float *sourceBuffer, float *destinationBuffer);
          void _getNormalizedContoursImage(const float *inputFrame, float *outputFrame);
          // -> special adaptive filters dedicated to low pass filtering on the chrominance (skeeps filtering on the edges)
          void _adaptiveSpatialLPfilter(const float *inputFrame,  float *outputFrame);
          void _adaptiveHorizontalCausalFilter_addInput(const float *inputFrame, float *outputFrame, const unsigned int IDrowStart, const unsigned int IDrowEnd); // TBB parallelized
          void _adaptiveVerticalAnticausalFilter_multGain(float *outputFrame, const unsigned int IDcolumnStart, const unsigned int IDcolumnEnd);
          void _computeGradient(const float *luminance);
          void _normalizeOutputs_0_maxOutputValue(void);
  
          // color space transform
          void _applyImageColorSpaceConversion(const std::valarray<float> &inputFrame, std::valarray<float> &outputFrame, const float *transformTable);
  
  #ifdef MAKE_PARALLEL
          /******************************************************
          ** IF some parallelizing thread methods are available, then, main loops are parallelized using these functors
          ** ==> main idea paralellise main filters loops, then, only the most used methods are parallelized... TODO : increase the number of parallelised methods as necessary
          ** ==> functors names = Parallel_$$$ where $$$= the name of the serial method that is parallelised
          ** ==> functors constructors can differ from the parameters used with their related serial functions
          */
  
          /* Template :
          class Parallel_ : public cv::ParallelLoopBody
          {
          private:
  
          public:
          Parallel_()
          : {}
  
          virtual void operator()( const cv::Range& r ) const {
  
          }
          }:
          */
          class Parallel_adaptiveHorizontalCausalFilter_addInput: public cv::ParallelLoopBody
          {
          private:
              float *outputFrame;
              const float *inputFrame, *imageGradient;
              unsigned int nbColumns;
          public:
              Parallel_adaptiveHorizontalCausalFilter_addInput(const float *inputImg, float *bufferToProcess, const float *imageGrad, const unsigned int nbCols)
                  :outputFrame(bufferToProcess), inputFrame(inputImg), imageGradient(imageGrad), nbColumns(nbCols) { }
  
              virtual void operator()( const Range& r ) const CV_OVERRIDE
              {
                  float* outputPTR=outputFrame+r.start*nbColumns;
                  const float* inputPTR=inputFrame+r.start*nbColumns;
                  const float *imageGradientPTR= imageGradient+r.start*nbColumns;
                  for (int IDrow=r.start; IDrow!=r.end; ++IDrow)
                  {
                      float result=0;
                      for (unsigned int index=0; index<nbColumns; ++index)
                      {
                          result = *(inputPTR++) + (*imageGradientPTR++)* result;
                          *(outputPTR++) = result;
                      }
                  }
              }
          };
  
          class Parallel_adaptiveVerticalAnticausalFilter_multGain: public cv::ParallelLoopBody
          {
          private:
              float *outputFrame;
              const float *imageGradient;
              unsigned int nbRows, nbColumns;
              float filterParam_gain;
          public:
              Parallel_adaptiveVerticalAnticausalFilter_multGain(float *bufferToProcess, const float *imageGrad, const unsigned int nbRws, const unsigned int nbCols, const float  gain)
                  :outputFrame(bufferToProcess), imageGradient(imageGrad), nbRows(nbRws), nbColumns(nbCols), filterParam_gain(gain) { }
  
              virtual void operator()( const Range& r ) const CV_OVERRIDE {
                  float* offset=outputFrame+nbColumns*nbRows-nbColumns;
                  const float* gradOffset= imageGradient+nbColumns*nbRows-nbColumns;
                  for (int IDcolumn=r.start; IDcolumn!=r.end; ++IDcolumn)
                  {
                      float result=0;
                      float *outputPTR=offset+IDcolumn;
                      const float *imageGradientPTR=gradOffset+IDcolumn;
                      for (unsigned int index=0; index<nbRows; ++index)
                      {
                          result = *(outputPTR) + *(imageGradientPTR) * result;
                          *(outputPTR) = filterParam_gain*result;
                          outputPTR-=nbColumns;
                          imageGradientPTR-=nbColumns;
                      }
                  }
              }
          };
  
          class Parallel_computeGradient: public cv::ParallelLoopBody
          {
          protected:
              float *imageGradient;
              const float *luminance;
              unsigned int nbColumns, doubleNbColumns, nbRows, nbPixels;
          public:
              Parallel_computeGradient(const unsigned int nbCols, const unsigned int nbRws, const float *lum, float *imageGrad)
              :imageGradient(imageGrad), luminance(lum), nbColumns(nbCols), doubleNbColumns(2*nbCols), nbRows(nbRws), nbPixels(nbRws*nbCols) { }
  
              virtual void operator()( const Range& r ) const CV_OVERRIDE {
                  for (int idLine=r.start;idLine!=r.end;++idLine)
                  {
                      for (unsigned int idColumn=2;idColumn<nbColumns-2;++idColumn)
                      {
                          const unsigned int pixelIndex=idColumn+nbColumns*idLine;
  
                          // horizontal and vertical local gradients
                          const float verticalGrad=fabs(luminance[pixelIndex+nbColumns]-luminance[pixelIndex-nbColumns]);
                          const float horizontalGrad=fabs(luminance[pixelIndex+1]-luminance[pixelIndex-1]);
  
                          // neighborhood horizontal and vertical gradients
                          const float verticalGrad_p=fabs(luminance[pixelIndex]-luminance[pixelIndex-doubleNbColumns]);
                          const float horizontalGrad_p=fabs(luminance[pixelIndex]-luminance[pixelIndex-2]);
                          const float verticalGrad_n=fabs(luminance[pixelIndex+doubleNbColumns]-luminance[pixelIndex]);
                          const float horizontalGrad_n=fabs(luminance[pixelIndex+2]-luminance[pixelIndex]);
  
                          const float horizontalGradient=0.5f*horizontalGrad+0.25f*(horizontalGrad_p+horizontalGrad_n);
                          const float verticalGradient=0.5f*verticalGrad+0.25f*(verticalGrad_p+verticalGrad_n);
  
                          // compare local gradient means and fill the appropriate filtering coefficient value that will be used in adaptative filters
                          if (horizontalGradient<verticalGradient)
                          {
                              imageGradient[pixelIndex+nbPixels]=0.06f;
                              imageGradient[pixelIndex]=0.57f;
                          }
                          else
                          {
                              imageGradient[pixelIndex+nbPixels]=0.57f;
                              imageGradient[pixelIndex]=0.06f;
                          }
                      }
                  }
              }
          };
  
  #endif
      };
  }// end of namespace bioinspired
  }// end of namespace cv
  
  #endif /*RETINACOLOR_HPP_*/
  