//*****************************************************************************
//       FillMissingCotEntry <runnumber> <dbid_write> <dbid_read>
//    this executable will make an DB entry using the beamline entries in the
//    same store.
//
//      06/20/02  Hartmut Stadie  IEKP Karlsruhe
//*****************************************************************************
//=============================================================================
// RCS Current Revision Record
//-----------------------------------------------------------------------------
// $Source: /cdf/code/cdfcvs/run2/Alignment/test/FillMissingCotEntry.cc,v $
// $Revision: 1.8 $
// $Date: 2003/02/26 16:57:35 $
// $Author: stadie $
// $State: Exp $
// $Locker:  $
//*****************************************************************************

#include <cmath>
#include <iostream>
#include <strstream>
#include <vector>

#include "CalibDB/RunListKey.hh"      
#include "Alignment/CotBeam.hh"
#include "AlignmentDB/CotBeamPosition.Defs.hh"
#include "AlignmentDB/SiAlignFrame.Defs.hh"
#include "CalibDB/Helpers.hh"
#include "CalibDB/RunListDefs.hh" 
#include "RunConfigDB/BEAM_CONDITIONS.Defs.hh"
#include "DBManager/ConControl.hh"
#include "DBManager/IoEntry.hh"
#include "CLHEP/Matrix/Vector.h"
#include "CLHEP/Matrix/SymMatrix.h"

//following the example in RunConfigDB/Test/BeamConditions.cc
int storeNumber(int run)
{
  static bool initialized = false;
  static PackageManager m;
  if( ! initialized)
    {
      ConControl::instance()->addIoEntry(new IoEntry("offline","OTL","cdf_reader/reader@cdfofprd"));
      m.setDefaultID("offline");
      initialized = true;
    }
  // Define the manager & container/iterator for returned values.
  static BEAM_CONDITIONS_mgr dmgr1("BEAM_CONDITIONS");
  static BEAM_CONDITIONSContainer_ptr rvals1;
  static BEAM_CONDITIONSContainer::iterator it1;
  int store = -1;
  // Set the key -- just simple run selection + cardname selection 
  std::ostringstream srt;
  //srt << run << ends;
  srt << run;
  std::string sql_string = "WHERE RC_RUNNUMBER = " ; 
  sql_string += srt.str();
  ROKey dbk(sql_string);
  
  // Do the fetch from the DB and print out returned rows.
  Result rc = dmgr1.get(dbk,rvals1);
  if (rc == Result::success)
    {
      //std::cout << "Rows got:" << std::endl;
      it1 = rvals1->begin();
      if(it1 != rvals1->end()) store = it1->tevStore();
    }
  else std::cout << "could not read RunConfigDB(error = "<< rc << ")" << std::endl;
  return store;
}

std::string keyAlgorithm(const std::string& dbid,
			 const RunListKey& rlk)
{
  RunList_mgr runIom(dbid,"RunList");
  RunList_var row;

  if( !runIom.isValid() ) 
    {
      std::cerr<< "Error opening RunList database" << std::endl;
      return rlk.algorithm();
    }
  
  if(( runIom.get(rlk,row)!=Result::success ) ||
     ( row->size() != 1 ))
    {
      std::cerr << "Failure of runIom.get() to retrieve db contents" << std::endl;  
      return rlk.algorithm();
    }
  //std::cout << "new rlk:" << *(row->begin()) << std::endl;
  return row->begin()->algorithm();
}

bool meanBeamLine(CotBeamPosition& newcot, std::vector<CotBeamPositionContainer> input)
{
  HepVector x(4,0);
  HepSymMatrix G(4,0);
  int ntracks = 0;
  double sumz = 0;
  double sumsigmaz = 0;
  
  for(std::vector<CotBeamPositionContainer>::const_iterator bpiter = input.begin() ;  bpiter != input.end() ; ++bpiter)
    {
      int ntracks2 = bpiter->begin()->Statistics0();
      HepVector x2(4,0);
      HepSymMatrix G2(4,0);

      x2[0] = bpiter->begin()->BeamX();
      x2[1] = bpiter->begin()->BeamY();
      x2[2] = bpiter->begin()->SlopeX();
      x2[3] = bpiter->begin()->SlopeY();
      
      G2[0][0] = bpiter->begin()->FitCov00();
      G2[0][1] = bpiter->begin()->FitCov01();
      G2[0][2] = bpiter->begin()->FitCov02();
      G2[0][3] = bpiter->begin()->FitCov03(); 
      G2[1][1] = bpiter->begin()->FitCov11();
      G2[1][2] = bpiter->begin()->FitCov12();
      G2[1][3] = bpiter->begin()->FitCov13();
      G2[2][2] = bpiter->begin()->FitCov22();
      G2[2][3] = bpiter->begin()->FitCov23();
      G2[3][3] = bpiter->begin()->FitCov33();
 
      sumz += ntracks2 * bpiter->begin()->BeamZ();
      sumsigmaz += ntracks2 * bpiter->begin()->WidthZ();     
      //       std::cout << "input:  n tracks:" << ntracks2 << std::endl;
//       std::cout << "x:" << x2;
//       std::cout << "G:" << G2;
      int ierr = 0;
      ntracks += ntracks2;
      x += x2 * ntracks2;
      G += G2.inverse(ierr) * ntracks2;
     
      if( ierr != 0 ) 
	{
	  std::cout << "FillMissingCotEntry: matrix inversion failed." 
		    << std::endl;
	  return false;
	}
    }
  int ierr = 0;
  x /= ntracks;
  G /= ntracks;
  G.invert(ierr);
  if( ierr != 0 ) 
    {
      std::cout << "FillMissingCotEntry: matrix inversion failed." 
		    << std::endl;
      return false;
    }
//   std::cout << "result:    n tracks:" << ntracks << std::endl;
//   std::cout << "x:" << x;
//   std::cout << "G:" << G;
  //fill the table
  std::ostringstream os;
  
 os << (long)1 << " ";
  os << x[0] << " " << x[1] << " " << x[2] << " " 
     << x[3] << " ";
  os << (double)0 << " " << (double)0 << " " << (double)0 << " ";
  os << sumz/ntracks << " " << sumsigmaz/ntracks << " ";
  os << 0 << " " << 0 << " " << 0 << " " << 0 << " ";
  os << 0 << " " << 0 << " " << 0 << " ";
  os << 0 << " " << 0 << " ";
  os << 0 << " " << 0 << " " << 0 << " " << 0 << " ";
  os << 0 << " " << 0 << " " << 0 << " " << 0 << " ";
  //cov
  os << G[0][0] << " " << G[0][1] << " " << G[0][2] << " " << G[0][3] << " ";
  os << G[1][1] << " " << G[1][2] << " " << G[1][3] << " ";
  os << G[2][2] << " " << G[2][3] << " ";
  os << G[3][3] << " ";
 
  //statistics
  os << ntracks << " " << 0 << " ";
  //flags
  std::vector<CotBeamPositionContainer>::const_iterator bpiter = input.begin();
  
  long flag0,flag1;
  long hist = (bpiter->begin()->Flag0() & 0xFFFF);
  long algo = (((bpiter->begin()->Flag0()>>16) & 0xFFFF) >= 20) ? 15 : 5;
  long srce = 10;// 5 online, 10 offline
  flag0 = (algo << 16) + hist;
  flag1 =  srce;
  os << flag0 << " " << flag1 << " ";
  //alignment
  os << 0 << " " << 0 << " " << 0 << " "<< 0 << " ";
  os << std::ends;
  
  std::istringstream is(os.str());    
  is >> newcot;
  return true;
}

int main(int argc, char* argv[])
{
  if(argc != 6)
    {
      std::cout << "USAGE: FillMissingCotEntry <runnumber> <history> <dbid_write> <dbid_read> <algorithm>"
		<< std::endl;
      return 1;
    }
  int runnumber = atoi(argv[1]);
  int history = atoi(argv[2]);
  std::string dbidwrite = argv[3];
  std::string dbidread  = argv[4];
  int algorithm = atoi(argv[5]);
  
  //get list with all runs that have a COT beamline
  RunList_mgr runIom(dbidread,"RunList");
  if( !runIom.isValid() ) 
    {
      std::cout << "FillMissingCotEntry: Error opening RunList database" 
	       << std::endl;
      return 1;
    }
  RunListKey rlk("GOOD","CotBeamPosition",Run::undefined,
		 Version::undefined);
  RunList_var runlist;
  if(runIom.get(rlk,runlist)!=Result::success) 
    {
      std::cout << "FillMissingCotEntry: Failure of runIom.get() to retrieve db contents" << std::endl;  
      return 1;
    }

  //find neighbors. 
  int diffdown = 1000000;
  int diffup = 1000000;
  int rundown = 0,runup = 0;
  int store = storeNumber(runnumber); 

  CotBeam cotbeam(dbidread);
  cotbeam.setHistoryCode(history);
  cotbeam.setAlgorithm(algorithm);
  for(RunList::iterator iter = runlist->begin(); iter != runlist->end();
      ++iter)
    {
      int newrun = iter->run();
      int diff = runnumber - newrun;
      
      //They should be within 50 runs.
      if(abs(diff) > 50) continue;
 
      if( diff == 0 )
	{  
	  //cotbeam.setAlgorithm(algorithm-5);
	  //if(cotbeam.loadRun(newrun) == 0)
	  //  {
	  //    std::cout << "FillMissingCotEntry: run " << runnumber << " has an entry already." << std::endl;
	  //    return 1;
	  //   }
	  // cotbeam.setAlgorithm(algorithm);
	  if(cotbeam.loadRun(newrun) == 0)
	    {
	      std::cout << "FillMissingCotEntry: run " << runnumber << " has an entry already." << std::endl;
	      return 1;
	    }
	}
      else if( diff < 0 )
	{
	  if(( -diff < diffdown ) && (storeNumber(newrun) == store))
	    {
	      //can we load the entry
	      if( cotbeam.loadRun(newrun) == 0 )
		{
		  diffdown = -diff;
		  rundown = newrun;
		}
	    }
	}
      else
	if(( diff < diffup ) && (storeNumber(newrun) == store))
	  {
	    //can we load the entry
	    if( cotbeam.loadRun(newrun) == 0 )
	      {
		diffup = diff;
		runup = newrun;
	      }
	  }
    }
 
  if( (runup == 0) && (rundown == 0)) 
    {
      std::cout << "FillMissingCotEntry: no other run with COT beam line information found for store " << store << "." << std::endl;
       return 2;
     }
  std::cout << "FillMissingCotEntry: found: store: " << store << " run= " << rundown << "," 
	    << runup << std::endl;
  
  //make vector for good entries
  std::vector<CotBeamPositionContainer> beampositions;
  RunListKey oldkey;
 
  if(runup != 0)
    {
      if( cotbeam.loadRun(runup) == 0 )
	{
	  CotBeamPositionContainer sbpc = cotbeam.getContainer();

	  beampositions.push_back(sbpc);
	  oldkey = cotbeam.getKey();
	}
    }
  if(rundown != 0)
    {
      if( cotbeam.loadRun(rundown) == 0 )
	{
	  CotBeamPositionContainer sbpc = cotbeam.getContainer();
	 
	  beampositions.push_back(sbpc);
	  oldkey = cotbeam.getKey();
	}
    }
  if(beampositions.empty())
    {
      std::cout << "FillMissingCotEntry: no good entry for the COT beam line found." << std::endl;
      return 3;
    }
  
  CotBeamPosition newbp;
  CotBeamPositionContainer oldcontainer;
  meanBeamLine(newbp, beampositions);

  std::cout << "FillMissingCotEntry:new beampos entry:" << newbp << std::endl;
   
//   //fill DB
  CotBeamPositionContainer contDB;
  CotBeamPosition_mgr iomCotBeam(dbidwrite,"CotBeamPosition");
  if(iomCotBeam.isValid()==false)
    {
      std::cout <<"FillMissingCotEntry: CotBeamPosition_mgr initialization failure" << std::endl;
      return  false;
    }                                                                
  RunListKey qk(runnumber,Version::generate); //just add a new entry
  qk.setDataStatus("GOOD");
  qk.newVersion();
  qk.setAlgorithm(keyAlgorithm(dbidread,oldkey));
  contDB.clear();
  contDB.push_back(newbp);
  //add parents 
  for(std::vector<CotBeamPositionContainer>::const_iterator bpiter = beampositions.begin() ;  
      bpiter != beampositions.end() ; ++bpiter)
    contDB.addParent(bpiter->key());
  
  if(iomCotBeam.put(qk,contDB) != Result::success)
    {
      std::cout << "FillMissingCotEntry: Put failed for CotBeamPosition Run = " << runnumber
		<< " in database " << dbidwrite << "." << std::endl; 
      return false;
    }
  std::cout << "FillMissingCotEntry: making COT DB entry...   done" << std::endl;
  
}