Commit 59191a42 authored by Austin Sanders's avatar Austin Sanders
Browse files

Clean history

parents
Loading
Loading
Loading
Loading

.gitignore

0 → 100644
+11 −0
Original line number Diff line number Diff line
*.make
*.o
*.so
*.so.*
*.a
.DS_Store
CMakeCache.txt
CMakeFiles/
CTestTestfile.cmake
Makefile
cmake_install.cmake

camera_models/CAHV.cpp

0 → 100644
+410 −0
Original line number Diff line number Diff line
#include <iostream>
#include <vector>
#include <limits>
#include <stdlib.h>
#include <iomanip>

#include "CAHV.h"
#include "../common/GeometricTransf.h"

using namespace std;

namespace at
{
    CAHV::CAHV()
    {

    }

    //not sure if this constructor is needed.
    CAHV::CAHV(vector<float> c, vector<float> a, vector<float> h, vector<float> v,
	     vector<float> offset, long int width, long int height)
    {
      m_height = height;
      m_width = width;
      m_offset.resize(3);
      m_c.resize(3);
      m_a.resize(3);
      m_h.resize(3);
      m_v.resize(3);
      for (int i=0; i < 3; i++){
	m_c[i]=c[i]+offset[i];
	m_a[i]=a[i];
	m_h[i]=h[i];
	m_v[i]=v[i];
	m_offset[i]=offset[i];
      }
    }

    CAHV::CAHV(vector<float> c, vector<float> a, vector<float> h, vector<float> v,
               vector<float> offset, vector<float> quaternion,
               long int width, long int height)
    {
        m_height = height;
        m_width = width;
        m_offset.resize(3);
        m_quaternion.resize(4);
        m_c.resize(3);
        m_a.resize(3);
        m_h.resize(3);
        m_v.resize(3);

        double c_in[3];
        double a_in[3];
        double h_in[3];
        double v_in[3];
        double temp_quaternion[4];

        for (int i=0; i < 3; i++){
	  c_in[i]=c[i];
	  a_in[i]=a[i];
	  h_in[i]=h[i];
	  v_in[i]=v[i];
	  m_offset[i]=offset[i];
        }
	
        for (int i=0; i < 4; i++){
            m_quaternion[i] = quaternion[i];
            temp_quaternion[i] = quaternion[i];
        }
  
        double c_out[3]; 
        rotateWithQuaternion( c_in, c_out, temp_quaternion);
        double a_out[3]; 
        rotateWithQuaternion( a_in, a_out, temp_quaternion);
        double h_out[3];
        rotateWithQuaternion( h_in, h_out, temp_quaternion);
        double v_out[3];
        rotateWithQuaternion( v_in, v_out, temp_quaternion);

        for (int i=0; i < 3; i++){
	  m_c[i]=c_out[i] + m_offset[i];//Sept 02, 2014 with Larry
	  m_a[i]=a_out[i];
	  m_h[i]=h_out[i];
	  m_v[i]=v_out[i];	
        }
    }

    CAHV::CAHV(string cahv_calibration)
    {
	std::ifstream fin (cahv_calibration.c_str());
	std::string line;

	bool found_offset = false;
	bool found_quaternion = false;
	int found = 0;

	double c[3],  a[3], h[3], v[3];
	m_quaternion.resize(4);
	m_offset.resize(3);

	m_height = -1;
	m_width = -1;

	// Check for proper loading
	if (fin.is_open()){// File loaded correctly
	  // Continue loading until the end of the file
	  while (!fin.eof()){
	    line.clear();
	    std::getline(fin, line);
	    
	    if (!line.empty() && line.at(0) != '#'){ // skip comment & empty lines

	      std::stringstream lineStream(line);
	      std::string keyword("");
	      lineStream >> keyword;
	      
	      if (keyword.compare(std::string("WIDTH_HEIGHT"))==0){
		found++;
	       
		// Check number of parameters
		vector<string> contents = split(line);
		if (contents.size() < 3) {// Not enough parameters
		  // Throw failed load
		  throw 1;
		}
		
		// Get the parameters
		m_width = atof(contents[1].c_str());
		m_height = atof(contents[2].c_str());
	      }
	      else if (keyword.compare(std::string("C"))==0){
		  found++;
		  
		  // Check number of parameters
		  vector<string> contents = split(line);
		  if (contents.size() < 4){// Not enough parameters
		    // Throw failed load
		    throw 1;
		  }
		  
		  // Get the parameters
		  c[0] = atof(contents[1].c_str());
		  c[1] = atof(contents[2].c_str());
		  c[2] = atof(contents[3].c_str());	
	      }
	      else if (keyword.compare(std::string("A"))==0){
		found++;
		lineStream>> a[0] >> a[1] >> a[2];
		
		// Check number of parameters
		vector<string> contents = split(line);
		if (contents.size() < 4){ // Not enough parameters
		    // Throw failed load
		    throw 1;
		  }
		
		// Get the parameters
		a[0] = atof(contents[1].c_str());
		a[1] = atof(contents[2].c_str());
		a[2] = atof(contents[3].c_str());
	      }
	      else if (keyword.compare(std::string("H"))==0){
		found++;
		lineStream >> h[0] >> h[1] >> h[2];
		
		// Check number of parameters
		vector<string> contents = split(line);
		if (contents.size() < 4){// Not enough parameters
		  // Throw failed load
		  throw 1;
		}
		
		// Get the parameters
		h[0] = atof(contents[1].c_str());
		h[1] = atof(contents[2].c_str());
		h[2] = atof(contents[3].c_str());
	      }
	      else if (keyword.compare(std::string("V"))==0){
		found++;		  
		// Check number of parameters
		vector<string> contents = split(line);
		if (contents.size() < 4){// Not enough parameters 
		  // Throw failed load
		  throw 1;
		}
		
		// Get the parameters
		v[0] = atof(contents[1].c_str());
		v[1] = atof(contents[2].c_str());
		v[2] = atof(contents[3].c_str());
	      }
	      else if (keyword.compare(std::string("QUATERNION"))==0){
		found_quaternion = true;	  
		// Check number of parameters
		vector<string> contents = split(line);
		if (contents.size() < 5){ // Not enough parameters
		  // Throw failed load
		  throw 1;
		}
		
		// Get the parameters
		m_quaternion[0] = atof(contents[1].c_str());
		m_quaternion[1] = atof(contents[2].c_str());
		m_quaternion[2] = atof(contents[3].c_str());
		m_quaternion[3] = atof(contents[4].c_str());
	      }
	      else if (keyword.compare(std::string("OFFSET"))==0){

		found_offset = true;
		
		// Check number of parameters
		vector<string> contents = split(line);
		if (contents.size() < 4){// Not enough parameters
		  // Throw failed load
		  throw 1;
		}
		
		// Get the parameters
		m_offset[0] = atof(contents[1].c_str());
		m_offset[1] = atof(contents[2].c_str());
		m_offset[2] = atof(contents[3].c_str());
	      }
	    }
	  }
	  
	  // Check for correct number of parameters
	  if (found < 5){// Not enough parameters
	    // Throw failed load exception
	    throw 1;
	  }
	  
	  // Add offset to the current position if an offset was loaded
	  if (found_offset){// Offset was found
	    // Add the offset to the current position
	    for (int i = 0; i < 3; i++){
	      m_c[i] = c[i] + m_offset[i];
	    }
	  }
	  
	  else{ // No offset found
	    // Use c without modification
	    for (int i = 0; i < 3; i++){
	      m_c[i] = c[i];
	    }
	  }
	  
	  // Apply quaternion if one was loaded
	  if (found_quaternion){ // Quaternion was found
	    // Copy the quaternion
	    double temp_quaternion[4];
	    for (int i = 0; i < 4; i++){
	      temp_quaternion[i] = m_quaternion[i];
	    } 
	    
	    // Apply the quaternion
	    double c_out[3]; 
	    rotateWithQuaternion( c, c_out, temp_quaternion);
	    double a_out[3]; 
	    rotateWithQuaternion( a, a_out, temp_quaternion);
	    double h_out[3];
	    rotateWithQuaternion( h, h_out, temp_quaternion);
	    double v_out[3];
	    rotateWithQuaternion( v, v_out, temp_quaternion);
	    
	    // Save the results as the new cavh parameters
	    for(int i = 0; i < 3; i++){
	      m_c[i] = c_out[i];
	      m_a[i] = a_out[i];
	      m_h[i] = h_out[i];
	      m_v[i] = v_out[i];
	    }
	  }
	  
	  else{ // No quaternion was found
	    // Use parameters without modification
	    for(int i = 0; i < 3; i++){
	      m_c[i] = c[i];
	      m_a[i] = a[i];
	      m_h[i] = h[i];
	      m_v[i] = v[i];
	    }
	  }
	}
	
	else{// File did not load correctly
	  // Throw failed open exception
	  throw 0;
	}
    }
  
    //build cahv (in pixels) from pinhole
    //focalLength, opticalCenter, imgWidth and imgHeight in pixels
    CAHV::CAHV(Eigen::Vector2d focalLength, Eigen::Vector2d opticalCenter, int imgWidth, int imgHeight, 
	       Eigen::Matrix3d rotation, Eigen::Vector3d translation)
    {
      
      //m_c = rotation*centerProj + translation;
      for (int i=0; i < 3; i++){  
	m_c[i]=translation[i];  
      }

      for (int i = 0; i < 3; i++){
	m_a[i] = rotation(2,i);
      }

      //vr_real32 fH = f/pixelSize[0], fV = f/pixelSize[1];

      for (int i = 0; i < 3; i++){
	m_h[i] = focalLength[0] * rotation(0,i) + opticalCenter[0] * m_a[i];
      }
      for (int i = 0; i < 3; i++){
	m_v[i] = focalLength[1] * rotation(1,i) + opticalCenter[1] * m_a[i];
      }

      m_height = (long int)imgHeight;
      m_width = (long int)imgWidth;
    }

    vector<float> CAHV::pixelToVector(vector<int> thisPixel)
    {
        //vertical component
        Eigen::Vector3f va = m_v + (-thisPixel[1])*m_a;
        //horizontal component
        Eigen::Vector3f vb = m_h + (-thisPixel[0])*m_a;
        //cross product
        Eigen:: Vector3f pVec = va.cross(vb);
        //normalization
        //cout<<"pVec="<<pVec<<endl;
        //cout<<"sqnorm="<<sqrt(pVec.squaredNorm())<<endl;
        pVec = pVec / pVec.norm();
        //cout<<"nVec="<<pVec<<endl;

        // The vector should be pointing in the same directions as A, if it
        // isn't, flip it:
        if (pVec.dot(m_a) < 0.0)
            pVec *= -1;

        vector<float> thisVector(3);

        thisVector[0] = pVec[0];
        thisVector[1] = pVec[1];
        thisVector[2] = pVec[2];


        return thisVector;
    }


  vector<int>
  CAHV::vectorToPixel(vector<float> thisVector)
  {
    vector<int> thisPixel(2);
    Eigen::Vector3f pVec(thisVector.data());
    double dDot = pVec.dot(m_a);

    // Ara probably wants to make this condition > 0.0... just sayin'
    // -- LJE
    if (dDot != 0.0)
    {
      thisPixel[0] = pVec.dot(m_h) / dDot;
      thisPixel[1] = pVec.dot(m_v) / dDot;
    }
    else
    {
      thisPixel[0] = numeric_limits<int>::max();
      thisPixel[1] = numeric_limits<int>::max();
    }

    return thisPixel;
  }

  void CAHV::getImagePlane(double a_unit[3], double h_unit[3], double v_unit[3])
  {
    // Extract image plane vectors from H and V
    Eigen::Vector3f h_img_plane = (m_h - m_a.dot(m_h) * m_a) / (m_a.cross(m_h).norm());
    Eigen::Vector3f v_img_plane = (m_v - m_a.dot(m_v) * m_a) / (m_a.cross(m_v).norm());

    // Copy to the output vectors
    for(int i = 0; i < 3; i++)
    {
      a_unit[i] = m_a(i);
      h_unit[i] = h_img_plane(i);
      v_unit[i] = v_img_plane(i);
    }
  }

  // Writes the current CAHV model to a file so that it can be loaded later
  void CAHV::writeFile(string outputFilename)
  {
    // Open the file
    ofstream output(outputFilename.c_str());

    // Write all of the parameters
    output << std::setprecision(8) << "C " << m_c[0] << " " << m_c[1] << " " << m_c[2] << endl;
    output << std::setprecision(8) << "A " << m_a[0] << " " << m_a[1] << " " << m_a[2] << endl;
    output << std::setprecision(8) << "H " << m_h[0] << " " << m_h[1] << " " << m_h[2] << endl;
    output << std::setprecision(8) << "V " << m_v[0] << " " << m_v[1] << " " << m_v[2] << endl;

    output << "WIDTH_HEIGHT " << m_width << " " << m_height << endl;
  }
  void CAHV::print()
  {
    cout << "CAHV model params: " << endl;
    cout << "c = " << m_c[0] << ", " << m_c[1] << ", " << m_c[2] << endl;
    cout << "a = " << m_a[0] << ", " << m_a[1] << ", " << m_a[2] << endl;
    cout << "h = " << m_h[0] << ", " << m_h[1] << ", " << m_h[2] << endl;
    cout << "v = " << m_v[0] << ", " << m_v[1] << ", " << m_v[2] << endl;
  }

}

camera_models/CAHV.h

0 → 100644
+75 −0
Original line number Diff line number Diff line
#ifndef _CAHV_H_
#define _CAHV_H_

#include <iostream>
#include <cmath>
#include <Eigen/Geometry>
using namespace std;

namespace at 
{

    class CAHV
    {
    public:

        long int m_height;
        long int m_width;
        Eigen::Vector3f m_c;
        Eigen::Vector3f m_a;
        Eigen::Vector3f m_h;
        Eigen::Vector3f m_v;
        vector<float> m_offset;
        vector<float> m_quaternion;

        CAHV();
        CAHV(vector<float> c, vector<float> a, vector<float> h, vector<float> v,
             vector<float>offset, long int width, long int height);
        CAHV(vector<float> c, vector<float> a, vector<float> h, vector<float> v,
             vector<float> offset, vector<float> quaternion,
             long int width, long int height);
	CAHV(Eigen::Vector2d focalLength, Eigen::Vector2d opticalCenter, int imgWidth, int imgHeight, 
	     Eigen::Matrix3d rotation, Eigen::Vector3d translation);

	// This constructor will throw exceptions when loading fails
        // This will delete the object
        // Will apply a quaternion if one is found
        // 0 means that the file was not found
        // 1 means that the file was misconfigured
	CAHV(string cahv_calibration);
        ~CAHV() { };
    
        vector<float> pixelToVector(vector<int> thisPixel);
        vector<int>   vectorToPixel(vector<float> thisVector);

	void getImagePlane(double a_unit[3], double h_unit[3], double v_unit[3]);

	// Writes the current CAHV model to a file so that it can be loaded later
        // Saves final CAHV, will not save a quaternion if one was used to create this object
	void writeFile(std::string outputFilename);
	void print();

    private:

	// From: http://stackoverflow.com/questions/236129/how-to-split-a-string-in-c
	// These functions split a string given a delimiter
	// Should these go to common/string_util
	std::vector<std::string> &split(const std::string &s, char delim, std::vector<std::string> &elems) {
	    std::stringstream ss(s);
	    std::string item;
	    while (std::getline(ss, item, delim)) {
		elems.push_back(item);
	    }
	    return elems;
	}


	std::vector<std::string> split(const std::string &s, char delim = ' ') {
	    std::vector<std::string> elems;
	    split(s, delim, elems);
	    return elems;
	}
    };
}

#endif
+32 −0
Original line number Diff line number Diff line
project(CAMERA_MODELS)
cmake_minimum_required(VERSION 2.6 FATAL_ERROR)

set( CMAKE_RUNTIME_OUTPUT_DIRECTORY ${PROJECT_SOURCE_DIR}/tests )
set( CMAKE_MODULE_PATH ${PROJECT_SOURCE_DIR}/../cmake/Modules/ ${CMAKE_MODULE_PATH} )

include(SetupInstall)
include(SetupRPATH)

find_package(Eigen3)

include_directories(${EIGEN3_INCLUDE_DIR} ../common)

setup_rpath()

#set(CMAKE_VERBOSE_MAKEFILE ON)
add_executable (cahv_to_pinhole cahv_to_pinhole.cc 
		CAHV.cc  pinhole.cc
		../common/GeometricTransf.cc
		../common/StringUtils.cc)
target_link_libraries (cahv_to_pinhole)

add_executable (pinhole_to_cahv pinhole_to_cahv.cc 
		CAHV.cc  pinhole.cc
		../common/GeometricTransf.cc
		../common/StringUtils.cc)
target_link_libraries (pinhole_to_cahv)




+286 −0
Original line number Diff line number Diff line
//Eigen includes
#include <Eigen/Dense>
#include <Eigen/Geometry>

#include "IsisInterfaceATK.h"


using namespace std;
using namespace Isis;

namespace at
{
// A program demonstrating how to call point_to_pixel(), pixel_to_vector(),
// and camera_center() for an ISIS frame camera.

// Very important note: This code assumes that pixel indices start from 0,
// not from 1 like in ISIS.

// Utilities
string print_vec(vector<double> const& vec){
  ostringstream oss;
  oss.precision(18);
  oss << "(";
  for (int s = 0; s < (int)vec.size()-1; s++){
    oss << vec[s] << ", ";
  }
  if (!vec.empty()) oss << vec[vec.size()-1];
  oss << ")";
  return oss.str();
}

vector<double> diff(vector<double> const& vec1, vector<double> const& vec2){
  vector<double> ans = vec1;
  for (int i = 0; i < (int)vec1.size(); i++){
    ans[i] -= vec2[i];
  }
  return ans;
}

vector<double> prod(double val, vector<double> const& vec){
  vector<double> ans = vec;
  for (int i = 0; i < (int)vec.size(); i++){
    ans[i] *= val;
  }
  return ans;
}

vector<double> div(vector<double> const& vec, double val){
  vector<double> ans = vec;
  for (int i = 0; i < (int)vec.size(); i++){
    ans[i] /= val;
  }
  return ans;
}

vector<double> normalize(vector<double> const& vec){
  double len = 0;
  for (int i = 0; i < (int)vec.size(); i++){
    len += vec[i]*vec[i];
  }
  len = sqrt(len);
  vector<double> ans = vec;
  for (int i = 0; i < (int)vec.size(); i++){
    ans[i] /= len;
  }
  return ans;
}

  IsisInterface::IsisInterface( string const& file ) {
    // Opening labels and camera
    Isis::FileName ifilename( QString::fromStdString(file) );
    m_label = new Isis::Pvl();
    m_label->read( ifilename.expanded() );
    
    // Opening Isis::Camera
    m_cube = new Isis::Cube(QString::fromStdString(file));
    m_camera = Isis::CameraFactory::Create( *m_cube );
  }
  
  IsisInterface::~IsisInterface(){
    delete m_label;
    delete m_camera;
    delete m_cube;
  }
  
  IsisInterface* IsisInterface::open( string const& filename ) {
    IsisInterface* result;
    
    // Opening Labels (This should be done somehow though labels)
    Isis::FileName ifilename( QString::fromStdString(filename) );
    Isis::Pvl label;
    label.read( ifilename.expanded() );
    
    Isis::Cube tempCube(QString::fromStdString(filename));
    Isis::Camera* camera = Isis::CameraFactory::Create( tempCube );
    
    switch ( camera->GetCameraType() ) {
    case 0:
      // Framing Camera
      if ( camera->HasProjection() ){
	throw string("Don't support Isis Camera Type with projection at this moment");
	//result = new IsisInterfaceMapFrame( filename );
      }
      else{
	result = new IsisInterfaceFrame( filename );
      }
      break;
    case 2:
      // Linescan Camera
      if ( camera->HasProjection() ){
	throw string("Don't support Projected Linescan Isis Camera Type at this moment");
	//result = new IsisInterfaceMapLineScan( filename );
      }
      else{
	result = new IsisInterfaceLineScan( filename );
      }
      break;
    default:
      throw string("Don't support unknown Isis Camera Type at this moment");
    }
    return result;
  }
  
  IsisInterfaceFrame::IsisInterfaceFrame( string const& filename ) :
    IsisInterface(filename), m_alphacube( *m_cube ) {
    
    // Gutting Isis::Camera
    m_distortmap = m_camera->DistortionMap();
    m_focalmap   = m_camera->FocalPlaneMap();
    m_detectmap  = m_camera->DetectorMap();
    
    // Calculating Center (just once)
    m_center.resize(3);
    m_camera->instrumentPosition(&m_center[0]);
    for (int i = 0; i < (int)m_center.size(); i++){
      m_center[i] *= 1000; // convert to meters
    }
    //Calculating Rotation (just once)
    std::vector<double> rot_inst = m_camera->instrumentRotation()->Matrix();
    std::vector<double> rot_body = m_camera->bodyRotation()->Matrix();
    Eigen::Matrix3d RotInst;
    RotInst(0,0)=rot_inst[0];   RotInst(0,1)=rot_inst[1];   RotInst(0,2)=rot_inst[2];
    RotInst(1,0)=rot_inst[3];   RotInst(1,1)=rot_inst[4];   RotInst(1,2)=rot_inst[5];
    RotInst(2,0)=rot_inst[6];   RotInst(2,1)=rot_inst[7];   RotInst(2,2)=rot_inst[8];
    Eigen::Matrix3d RotBody;
    RotBody(0,0)=rot_body[0];   RotBody(0,1)=rot_body[1];   RotBody(0,2)=rot_body[2];
    RotBody(1,0)=rot_body[3];   RotBody(1,1)=rot_body[4];   RotBody(1,2)=rot_body[5];
    RotBody(2,0)=rot_body[6];   RotBody(2,1)=rot_body[7];   RotBody(2,2)=rot_body[8];
    Eigen::Matrix3d Rot = RotBody*(RotInst.transpose());
    m_rotation.resize(9);
    m_rotation[0]=Rot(0,0);    m_rotation[1]=Rot(0,1);   m_rotation[2]=Rot(0,2);
    m_rotation[3]=Rot(1,0);    m_rotation[4]=Rot(1,1);   m_rotation[5]=Rot(1,2);
    m_rotation[6]=Rot(2,0);    m_rotation[7]=Rot(2,1);   m_rotation[8]=Rot(2,2);
    //= m_camera->bodyRotation()->Matrix();
    //R_body*transpose(R_inst)
  }
  
  vector<double>
  IsisInterfaceFrame::point_to_pixel( vector<double> const& point ) const{
    
    // Note: The smallest pixel is (0, 0), not (1, 1)!
    
    vector<double> look = normalize( diff(point, m_center) );
    
    vector<double> lookB_copy(3);
    std::copy( look.begin(), look.end(), lookB_copy.begin() );
    lookB_copy = m_camera->bodyRotation()->J2000Vector(lookB_copy);
    lookB_copy = m_camera->instrumentRotation()->ReferenceVector(lookB_copy);
    std::copy( lookB_copy.begin(), lookB_copy.end(), look.begin() );
    look = div(prod(m_camera->FocalLength(), look), fabs(look[2]) );
    
    // Back Projecting
    m_distortmap->SetUndistortedFocalPlane( look[0], look[1] );
    m_focalmap->SetFocalPlane( m_distortmap->FocalPlaneX(),
			       m_distortmap->FocalPlaneY() );
    m_detectmap->SetDetector( m_focalmap->DetectorSample(),
			      m_focalmap->DetectorLine() );
    
    vector<double> ans(2);
    ans[0] = m_alphacube.BetaSample( m_detectmap->ParentSample() ) - 1;
    ans[1] = m_alphacube.BetaLine( m_detectmap->ParentLine() ) - 1;
    return ans;
  }
  
  vector<double> 
  IsisInterfaceFrame::pixel_to_vector( vector<double> const& px ) const {
    
    // Note: The smallest pixel is (0, 0), not (1, 1)!
    
    m_detectmap->SetParent( m_alphacube.AlphaSample(px[0]+1),
			    m_alphacube.AlphaLine(px[1]+1));
    m_focalmap->SetDetector( m_detectmap->DetectorSample(),
			     m_detectmap->DetectorLine() );
    m_distortmap->SetFocalPlane( m_focalmap->FocalPlaneX(),
				 m_focalmap->FocalPlaneY() );
    vector<double> look(3);
    look[0] = m_distortmap->UndistortedFocalPlaneX();
    look[1] = m_distortmap->UndistortedFocalPlaneY();
    look[2] = m_distortmap->UndistortedFocalPlaneZ();
    
    look = normalize( look );
    look = m_camera->instrumentRotation()->J2000Vector(look);
    look = m_camera->bodyRotation()->ReferenceVector(look);
    return look;
  }
  
  vector<double> IsisInterfaceFrame::camera_center() const{
    return m_center;
  }
  vector<double> IsisInterfaceFrame::camera_rotation() const{
    return m_rotation;
  }

  // IsisInterfaceLineScan class implementation
  IsisInterfaceLineScan::IsisInterfaceLineScan( std::string const& filename ):
    IsisInterface(filename), m_alphacube( *m_cube ) {
    
    // Gutting Isis::Camera
    m_distortmap = m_camera->DistortionMap();
    m_focalmap   = m_camera->FocalPlaneMap();
    m_detectmap  = m_camera->DetectorMap();  
  }
  
  vector<double>
  IsisInterfaceLineScan::point_to_pixel( vector<double> const& point ) const
  {
    using namespace Isis;
    vector<double> pix(2);
    
  // Note: The smallest pixel is (0, 0), not (1, 1)!
  SurfacePoint surfP(Displacement(point[0], Displacement::Meters),
    Displacement(point[1], Displacement::Meters),
    Displacement(point[2], Displacement::Meters));
  
  if(m_camera->SetUniversalGround(surfP.GetLatitude().degrees(),
				  surfP.GetLongitude().degrees(),
				  surfP.GetLocalRadius().meters())){
    pix[0] = m_camera->Sample()-1;
    pix[1] = m_camera->Line()-1;
    //cout<<"pix="<<pix[0]<<", "<<pix[1]<<endl;
   }
  else{
    cout<<"Could not project point into the camera."<<endl;
    pix[0]=0;
    pix[1]=0;
    //cout<<"lat="<<surfP.GetLatitude().degrees()<<endl;
    //cout<<"lon="<<surfP.GetLongitude().degrees()<<endl;
    //cout<<"rad="<<surfP.GetLocalRadius().meters()<<endl;
  }
  return pix;
 }
  
vector<double> 
IsisInterfaceLineScan::pixel_to_vector( vector<double> const& px ) const {

  // Note: The smallest pixel is (0, 0), not (1, 1)!

  m_detectmap->SetParent( m_alphacube.AlphaSample(px[0]+1),
                          m_alphacube.AlphaLine(px[1]+1));
  m_focalmap->SetDetector( m_detectmap->DetectorSample(),
                           m_detectmap->DetectorLine() );
  m_distortmap->SetFocalPlane( m_focalmap->FocalPlaneX(),
                               m_focalmap->FocalPlaneY() );

  vector<double> look(3);
  look[0] = m_distortmap->UndistortedFocalPlaneX();
  look[1] = m_distortmap->UndistortedFocalPlaneY();
  look[2] = m_distortmap->UndistortedFocalPlaneZ();

  look = normalize( look );
  look = m_camera->instrumentRotation()->J2000Vector(look);
  look = m_camera->bodyRotation()->ReferenceVector(look);

  return look;
}

vector<double> IsisInterfaceLineScan::camera_center() const{
  throw string("Computation of camera center is not implemented")
    + " for linescan cameras.";
}
vector<double> IsisInterfaceLineScan::camera_rotation() const{
  throw string("Computation of camera rotation is not implemented")
    + " for linescan cameras.";
}
}