ExampleProjectionParametersPerspective.cpp

/* 
This file is part of the BIAS library (Basic ImageAlgorithmS).

Copyright (C) 2003-2009    (see file CONTACT for details)
Multimediale Systeme der Informationsverarbeitung
Institut fuer Informatik
Christian-Albrechts-Universitaet Kiel


BIAS is free software; you can redistribute it and/or modify
it under the terms of the GNU Lesser General Public License as published by
the Free Software Foundation; either version 2.1 of the License, or
(at your option) any later version.

BIAS is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
GNU Lesser General Public License for more details.

You should have received a copy of the GNU Lesser General Public License
along with BIAS; if not, write to the Free Software
Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
*/


/**
   @example ExampleProjectionParametersPerspective.cpp
   @relates ProjectionParametersPerspective
   @brief small example demonstrating the ProjectionParametersPerspective
   by calculating the lens distortion after brown (?)
   @ingroup g_examples
   @author grauel 9/2007
*/

#include <Geometry/ProjectionParametersPerspective.hh>
#include <Base/Image/Image.hh>
#include <Base/Image/ImageIO.hh>

using namespace BIAS;
using namespace std;


int main(int argc, char *argv[])
{
  ProjectionParametersPerspective ppp;

// used .IOR - parameters
// 1                     
//-999        
//-20.30377    cam const in mm (?)
//-0.06444     princpx    normalized or mm (?)
//-0.15281     princpy
//-2.87698e-004  k1
//6.73824e-007   k2
// .85           r0
//-6.23736e-010  
// 1.97915e-005  k3
//-2.11616e-005  k4
//-1.38247e-004  ascpect
// 8.89027e-006  skew
// 23.60000      sensorsize in mm
// 15.80000      sensorsize in mm
// 3872          pixel x
// 2592          pixel y 

  double CamConst = 20.30377; //in mm - perpendicular distance of projection center to image plane
  double SensorSizeX = 23.6;  //in mm - 
  double SensorSizeY = 15.8;  //in mm

  unsigned int ImgWidth  = 3872;
  unsigned int ImgHeight = 2592;

  double principalX = (double) ImgWidth*0.5  *  (1-0.06444/SensorSizeX);
  double principalY = (double) ImgHeight*0.5 *  (1-0.15281/SensorSizeY);
  double focalX = (double)ImgWidth*0.5*(CamConst/SensorSizeX );//in pixel
  double aspectrat = (1-1.38247e-004)*(double)ImgWidth/(double)ImgHeight;
  
  //ppp.SetDistortionType(DISTYPE_BROWN);  //not neede anymore
                                           //also set by calling
                                           //SetUndistortionBrown()
  double kc1 = -2.87698e-004;
  double kc2 = 6.73824e-007 ;
  double kc3 = 1.97915e-005;
  double kc4 = -2.11616e-005;

  double r0  = .85;  // second constant root of polynomial

  ppp.SetUndistortionBrown(kc1,kc2,kc3,kc4,r0);
  ppp.SetFocalLengthAndAspect(focalX,aspectrat);
  ppp.SetPrincipal(principalX,principalY);    
  ppp.SetIdealImageSize(ImgWidth,ImgHeight);

#ifdef BIAS_HAVE_XML2
  string fname = "ProjParPersExample01.xml";
  int res = ppp.XMLWrite(fname);

  if (res!=0)
  {
    cerr << "error writing "<< fname << endl;
    return -3;
  }
  cout << "XMLWrite wrote "<< fname << endl;
  
  kc1 = 0;
  kc2 = 0;
  kc3 = 0;
  kc4 = 0;
  r0  = 0;

  res = ppp.XMLRead(fname);
  if (res!=0)
  {
    cerr << "error reading "<< fname << endl;
    return -3;
  }
  cout << "XMLWrite read "<< fname << endl;  

  BIAS_ProjParaPersp_DISTORTION_TYPE distype = ppp.GetDistortionType();
  switch(distype)
  {
    case DISTYPE_BROWN:
      cout<<"the file contained brown distortion parameters - as expected. everything went fine."<<endl;
      ppp.GetUndistortionBrown(kc1,kc2,kc3,kc4,r0);
      cout<<"parameters are: kc1="<<kc1<<" "<< 
                            "kc2="<<kc2<<" "<<
                            "kc3="<<kc3<<" "<<
                            "kc4="<<kc4<<" "<<
                            "r0="<<r0<<" "<<endl;
    break;
    default:
    cerr<<"unexpected parameter type. check io functions!"<<endl;
  }

#else
  cout<<"did not XMLWrite because BIAS_HAVE_XML2 not available."<<endl;
#endif

  //generate distortion images
  ImageBase DistImgVert;
  ImageBase DistImgHorz; 
  ImageBase DistImgMagn;

  DistImgHorz.Init(ImgWidth, ImgHeight, 1, ImageBase::ST_float);
  DistImgVert.Init(ImgWidth, ImgHeight, 1, ImageBase::ST_float);
  DistImgMagn.Init(ImgWidth, ImgHeight, 1, ImageBase::ST_float);

  float** pImgHorz=((float **)(DistImgHorz.GetImageDataArray()));
  float** pImgVert=((float **)(DistImgVert.GetImageDataArray()));
  float** pImgMagn=((float **)(DistImgMagn.GetImageDataArray()));

  for(unsigned x = 0; x<ImgWidth; x++)
    for(unsigned y = 0; y<ImgHeight; y++)
    {
      HomgPoint2D point2d(x,y);
      
      ppp.Undistort(point2d);
      //ppp.Distort(point2d);
      
      //calculate magnitude of horizontal displacement 
      float disx = (float) point2d[0] - (float) x;      
      pImgHorz[y][x*DistImgHorz.GetChannelCount()] = disx;
      float disy = (float) point2d[1] - (float) y;      
      pImgVert[y][x*DistImgVert.GetChannelCount()] = disy;
      float magn = sqrt((disx*disx)+(disy*disy));
      pImgMagn[y][x*DistImgMagn.GetChannelCount()] = magn;
    }
  
  if(ImageIO::Save("DisplacementX.mip",DistImgHorz)==-1)
  {
    cerr<<" can't save displacement image." << endl;
    return -3;
  }
  if(ImageIO::Save("DisplacementY.mip",DistImgVert)==-1)
  {
    cerr<<" can't save displacement image." << endl;
    return -3;
  }
  if(ImageIO::Save("DisplacementMagn.mip",DistImgMagn)==-1)
  {
    cerr<<" can't save displacement image." << endl;
    return -3;
  }

  cout<< " displacement images saved." << endl;

  return 0;
}