#ifndef _NOMAGICK_
#include <Magick++.h>
#endif

#include <iostream>
#include <string>
#include <fstream>
#include <math.h>
#include "imagealign.h"
#include "image2d.h"
#include <fftw3.h>
#include <tclap/CmdLine.h>
#include "bdbenv.h"
#include "dbimage.h"
#include "pso.h"
#include "histequalize.h"
#include "tile.h"
#include "tileimage.h"
#include "dbalnimage.h"
#include "alignimage.h"

int main( int argc, char *argv[] )
{
  try
    {
      ///set up the command line parsing
      TCLAP::CmdLine cmd("Image alignment program", ' ', "0.001");

      TCLAP::ValueArg<std::string> dbdir("d", "db", "database dirctory", true, "", "databaase");

      TCLAP::ValueArg<int> todo("a", "aid", "aligned image id", true, 0, "aligned image id" );

      TCLAP::ValueArg<int> todo1("b", "aid1", "aligned image id", true, 0, "aligned image id to align" );

      TCLAP::ValueArg<std::string> name("n", "name", "job name", true, "", "name");

      TCLAP::ValueArg<std::string> track("t", "trackdb", "tracker sqlite file", false, "", "tracker output");

      //      clustermanager cm;
#ifndef _NOMAGICK_
      Magick::Image mim;
#endif

      cmd.add(dbdir);
      cmd.add(todo);
      cmd.add(todo1);
      cmd.add( name );
      cmd.add( track );

      cmd.parse(argc, argv);

      bdbenv bdb;
      dbimage imstore;
      dbalnimage dba;

      bdb.open( dbdir.getValue() );

      imstore.init( bdb.envptr() );
      dba.init( bdb.envptr() );

      bdb.adddb( &imstore );
      bdb.adddb( &dba );

      bdb.cleanlog();

      alignimage aim;
      aim.init( &dba );
      aim.load( );

      std::vector< alignimage::alignjob > jobs;

      aim.getrowalnjobs( jobs, 0 );

      std::vector< image2d* > imsa, imsb;
      std::vector< int > alnids;
      //      for( int i = 5; i < 21; i++ )
      for( int i = 5; i < 15; i++ )
        {
          alnids.push_back( i );
          imsa.push_back( new image2d );
          imsb.push_back( new image2d );
          
          imsa.push_back( new image2d );
          imsb.push_back( new image2d );

          imsa.push_back( new image2d );
          imsb.push_back( new image2d );          
        }
      image2d imsamap, imsbmap;
      pso mypso;
      psotracker psostat;
      std::string trackdbname;
      trackdbname = track.getValue();
      if( trackdbname.length() > 0 )
        psostat.init( trackdbname );
      for( unsigned i = 0; i < jobs.size(); i++ )
        {
          if( jobs[i].aid2 == todo1.getValue() && jobs[i].aid1 == todo.getValue() )
            {
              std::cout << "aligning: " << jobs[i].aid1 << "-"
                        << jobs[i].aid2 << " " << jobs[i].lowx << " "
                        << jobs[i].lowy << " " <<  jobs[i].highx
                        << " " << jobs[i].highy << "\n";
              
              dba.get_aligned_mask( jobs[i].aid1, &imsamap );
              dba.get_aligned_mask( jobs[i].aid2, &imsbmap );
              
              histequalize heq;
              heq.init( QUANTIZEFACTOR );
              std::cout << "getting channel 0 images\n";
              std::vector< psotrack_slice > psoslicedata;
              psoslicedata.clear();

              for( unsigned j = 0; j < alnids.size(); j++ )
                {
                  psotrack_slice pts;
                  pts.psoslice = j;
                  pts.dbslice = alnids[j];
                  pts.dbtilex1 = pts.dbtiley1 = -1;
                  pts.dbalignedtile1 = jobs[i].aid1;
                  pts.dbtilex2 = pts.dbtiley2 = -1;
                  pts.dbalignedtile2 = jobs[i].aid2;
                  pts.dbchannel = 1;
                  psoslicedata.push_back( pts );
                  std::cout << "getting channel " << j << " 1 images\n";
                  dba.get_aligned_image( jobs[i].aid1, alnids[j], 1, imsa[j] );
                  dba.get_aligned_image( jobs[i].aid2, alnids[j], 1, imsb[j] );

                  pts.dbchannel = 0;
                  pts.psoslice = j+alnids.size();
                  psoslicedata.push_back( pts );
                  std::cout << "getting channel " << j << " 0 images\n";
                  dba.get_aligned_image( jobs[i].aid1, alnids[j], 0, imsa[j+alnids.size()] );
                  dba.get_aligned_image( jobs[i].aid2, alnids[j], 0, imsb[j+alnids.size()] );

                  pts.dbchannel = 2;
                  pts.psoslice = j+(2*alnids.size());
                  psoslicedata.push_back( pts );
                  std::cout << "getting channel " << j << " 2 images\n";
                  dba.get_aligned_image( jobs[i].aid1, alnids[j], 2, imsa[j+(2*alnids.size())] );
                  dba.get_aligned_image( jobs[i].aid2, alnids[j], 2, imsb[j+(2*alnids.size())] );
                }
              std::cout << "equalizing\n";
              for( unsigned j = 0; j < imsa.size(); j++ )
                {
                  heq.init( QUANTIZEFACTOR );                  
                  heq.addimage( imsa[j] );
                  heq.addimage( imsb[j] );

                  heq.equalize( imsa[j] );
                  heq.equalize( imsb[j] );
                }
              box mybox;
              mybox.lowx = jobs[i].lowx;
              mybox.lowy = jobs[i].lowy*0.75;
              mybox.highx = jobs[i].highx;
              mybox.highy = jobs[i].highy;
              if( trackdbname.length() > 0 )
                psostat.set_runconfig( psoslicedata, 
                                       800, mybox);
              std::cout << "optimizing...\n";  
              mypso.init( imsa, imsb, imsamap, imsbmap );

              if( trackdbname.length() > 0)
                mypso.mpoptimize( 800, jobs[i].lowx, 0.75*jobs[i].lowy, 
                                  jobs[i].highx, jobs[i].highy, name.getValue(),
                                  &psostat );
              else
                mypso.mpoptimize( 800, jobs[i].lowx, 0.75*jobs[i].lowy, 
                                  jobs[i].highx, jobs[i].highy, name.getValue(),
                                  NULL );

              std::vector< transformation > trans;
              trans.clear();
              mypso.getresult( trans );
              std::cout << "RESULTS " << jobs[i].aid1 << "-"
                        << jobs[i].aid2 << " " << jobs[i].lowx << " "
                        << jobs[i].lowy << " " <<  jobs[i].highx
                        << " " << jobs[i].highy << "\n";
              for( unsigned j = 0; j < trans.size(); j++ )
                std::cout << trans[j].transx << " " << trans[j].transy 
                          << " : " << trans[j].score << "\n";
            }
        }
      
      
      for( unsigned k = 0; k < imsa.size(); k++ )
        {
          if( imsa[k] != NULL )
            delete imsa[k];
          else
            std::cout << "ERROR: a " << alnids[k] << "\n";
        }
      for( unsigned k = 0; k < imsb.size(); k++ )
        {
          if( imsb[k] != NULL )
            delete imsb[k];
          else
            std::cout << "ERROR: b " << alnids[k] << "\n";
        }

      //      cm.end();      
      dba.close();
      imstore.close();
      bdb.close();
    }
  catch (TCLAP::ArgException &e)
    {
      std::cout << "ERROR: " << e.error() 
                << " for arg " << e.argId() << "\n";
      return 0;
    }
  return 0;
}
