#include <G4Tubs.hh>
#include <G4UnionSolid.hh>

#include <RAT/AVPipe.hh>
#include <RAT/Log.hh>
#include <RAT/GeoFactory.hh>
using namespace RAT;
using namespace RAT::geo;

#include <CLHEP/Units/SystemOfUnits.h>
#include <CLHEP/Units/PhysicalConstants.h>
#include <string>
using namespace CLHEP;
using namespace std;

G4VSolid*
AVPipe::Construct( const std::string& name,
                   DBLinkPtr table ) const
{
  const vector<double> xPositions = table->GetDArray( "x" );
  const vector<double> yPositions = table->GetDArray( "y" );
  const vector<double> zPositions = table->GetDArray( "z" );
  const double rMin = table->GetD( "r_min" ) * mm;
  const double rMax = table->GetD( "r_max" ) * mm;

  Log::Assert( xPositions.size() == yPositions.size() && xPositions.size() == zPositions.size(), "AVPipe::Construct: Position array size mismatch," );

  G4VSolid* pipe = ConstructPipe( name,
                                  G4ThreeVector( xPositions[0] * mm, yPositions[0] * mm, zPositions[0] * mm ),
                                  G4ThreeVector( xPositions[1] * mm, yPositions[1] * mm, zPositions[1] * mm ),
                                  rMin, rMax );

  for( size_t iPos = 2; iPos < xPositions.size(); iPos++ )
    {
      G4ThreeVector start( xPositions[iPos - 1] * mm, yPositions[iPos - 1] * mm, zPositions[iPos - 1] * mm );
      G4ThreeVector end( xPositions[iPos] * mm, yPositions[iPos] * mm, zPositions[iPos] * mm );
      pipe = new G4UnionSolid( name, pipe, ConstructPipe( name, start, end, rMin, rMax ), NULL, G4ThreeVector( 0.0, 0.0, 0.0 ) );
    }
  return pipe;
}

G4VSolid*
AVPipe::ConstructPipe( const std::string& name,
                       const G4ThreeVector start,
                       const G4ThreeVector end,
                       const double rMin,
                       const double rMax ) const
{
  const G4ThreeVector direction = (end - start);
  const double length = direction.mag();
  G4Tubs* pipe = new G4Tubs( name + "_solid", rMin, rMax, length / 2.0, 0.0, twopi );
  const G4ThreeVector mid = ( start + end ) / 2.0;
  G4RotationMatrix* rotation = GeoFactory::DirectionToRotation( direction );
  return new G4DisplacedSolid( name + "_displaced_solid", pipe, rotation, mid );
}