UsefulAdu5Pat.h

#ifndef USEFULADU5PAT_H
#define USEFULADU5PAT_H

//Includes
#include <TObject.h>
#include <TGraph.h>
#include <TVector3.h>
#include "Adu5Pat.h"
#include "AnitaConventions.h"
#include "RampdemReader.h"
#include "AnitaGeomTool.h"

namespace AnitaStaticAdu5Offsets{
  const Double_t heading = 0;
  const Double_t pitch = 0;
  const Double_t roll = 0;
}

class UsefulAdu5Pat: public Adu5Pat
{

 public:
  UsefulAdu5Pat();          
  UsefulAdu5Pat(const Adu5Pat *patPtr); 
  ~UsefulAdu5Pat();         

  int getSourceLonAndLatAtDesiredAlt(Double_t phiWave, Double_t thetaWave, Double_t &sourceLon, Double_t &sourceLat, Double_t desiredAlt);

  
  int getSourceLonAndLatAltZero(Double_t phiWave, Double_t thetaWave, Double_t &sourceLon, Double_t &sourceLat);
  
  int getSourceLonAndLatAtAlt(Double_t phiWave, Double_t thetaWave, Double_t &sourceLon, Double_t &sourceLat,Double_t &sourceAltitude);

  int getSourceLonAndLatAtAlt2(Double_t phiWave, Double_t thetaWave, Double_t &sourceLon, Double_t &sourceLat,Double_t &sourceAltitude, Int_t maxLoopIter = -1, TVector3* sourcePosOptional = NULL) const;



  TVector3 getUnitVectorAlongThetaWavePhiWave(double thetaWave, double phiWave) const;
  
  int getSourceLonAndLatAtAlt3(Double_t phiWave, Double_t thetaWave, Double_t &sourceLon, Double_t &sourceLat,Double_t &sourceAltitude,
                   double* deltaAltIfNoIntersectection = NULL, bool returnBestPositionIfNoIntersection = false) const;
  int traceBackToContinent3(Double_t phiWave, Double_t thetaWave, 
                Double_t * lon, Double_t * lat, Double_t *alt, Double_t * theta_adjustment_required) const;


  
  

  int traceBackToContinent(Double_t phiWave, Double_t thetaWave, 
                           Double_t * lon, Double_t * lat, Double_t *alt, Double_t * theta_adjustment_required, 
                           Double_t max_theta_adjustment = TMath::Pi()/180, Int_t max_iter = 10) const; 




  void getThetaAndPhiWave(Double_t sourceLon, Double_t sourceLat, Double_t sourceAlt, Double_t &thetaWave, Double_t &phiWave);


  void getThetaAndPhiWave2(Double_t sourceLon, Double_t sourceLat, Double_t sourceAlt, Double_t &thetaWave, Double_t &phiWave, TVector3* sourcePos = NULL) const;

  void getThetaAndPhiWaveCart(TVector3 * sourcePos, Double_t & thetaWave, Double_t & phiWave) const; 

  void getThetaAndPhiWaveOfRayAtInfinity(const TVector3 & p0, const TVector3 & v0 , Double_t & thetaWave, Double_t & phiWave,
                                          Bool_t plus_infinity = true, Double_t eps = 0.00001 * TMath::DegToRad(), Double_t step = 10e7, TVector3 * testPosition = 0) const; 

  void getThetaAndPhiWaveWaisDivide(Double_t &thetaWave, Double_t &phiWave);


  
  void getThetaAndPhiWaveLDB(Double_t &thetaWave, Double_t &phiWave);

  void getThetaWaveAtBase(Double_t baseLon, Double_t baseLat, Double_t baseAlt, Double_t &thetaWave);


  
  void getThetaAndPhiWaveWillySeavey(Double_t &thetaWave, Double_t &phiWave);

  
  /*  For a the Williams Field borehole antenna calculates the payload centric azimuthal and elevation angles of the plane wave incident at the payload.
   * @param phiWave Reference to a Double_t in which to store the azimuthal angle of plane wave (in payload centric coordinates with phi equals zero lying in the direction the ADU5 fore antenna)
   * @param thetaWave Reference to a Double_t in which to store the elevation angle of plane wave (in payload centric coordinates with phi equals zero lying in the direction the ADU5 fore antenna)
   */
  void getThetaAndPhiWaveWillyBorehole(Double_t &thetaWave,Double_t &phiWave);

  
  void getThetaAndPhiWaveTaylorDome(Double_t &thetaWave, Double_t &phiWave);


  
  Double_t getDeltaTExpected(Int_t ant1, Int_t ant2,Double_t sourceLon, Double_t sourceLat, Double_t sourceAlt); 


  
  Double_t getDeltaTExpected(Int_t ant1, Int_t ant2,Double_t phiWave, Double_t thetaWave);

  Double_t getDeltaTExpected(Int_t ant1,Int_t ant2,Double_t cosPhi,Double_t sinPhi,Double_t cosTheta,Double_t sinTheta);


  
  Double_t getDeltaTWillySeavey(Int_t ant1, Int_t ant2);


  
  Double_t getDeltaTWillyBorehole(Int_t ant1, Int_t ant2);


  
  Double_t getDeltaTTaylor(Int_t ant1, Int_t ant2);
  


  Double_t getPhiWave() const {
    return fPhiWave;
  }

  Double_t getThetaWave() const {
    return fThetaWave;
  }

  Double_t getSourceLongitude() const {
    return fSourceLongitude;
  }
  Double_t getSourceLatitude() const{
    return fSourceLatitude;
  }

  Double_t getSourceAltitude() const {
    return fSourceAltitude;
  }

  UInt_t getTaylorDomeTriggerTimeNs() const;                                
  UInt_t getWaisDivideTriggerTimeNs() const;                                
  UInt_t getSipleTriggerTimeNs() const;                                 
  UInt_t getLDBTriggerTimeNs() const;                                   
  UInt_t getTriggerTimeNsFromSource(Double_t sourceLat, Double_t sourceLong, Double_t sourceAlt) const; 
  Double_t getDistanceFromSource(Double_t sourceLat, Double_t sourceLong, Double_t sourceAlt) const;    
  void setIncludeGroupDelay(Int_t flag){fIncludeGroupDelay=flag;}                   

  RampdemReader *fRampdemReader;
  AnitaGeomTool * fUPGeomTool; 

  Double_t getAngleBetweenPayloadAndSource(Double_t sourceLon, Double_t sourceLat, Double_t sourceAlt); //ACG additional 
  void getSunPosition(Double_t& phiDeg, Double_t& thetaDeg) const;
  Double_t getAzimuthOfSunRelativeToNorth() const;
  Double_t getAzimuthOfSun() const;
  Double_t getDifferencePointingToSun(Double_t phiAngle, Bool_t inputInDegrees=true) const;

  void getThetaAndPhiWaveHiCal(Double_t& thetaWave, Double_t& phiWave); 


  static Double_t getReflectionAngle(Double_t plAlt, Double_t imageEl, Double_t imageAlt);

  Double_t getReflectionAngle(Double_t imageEl, Double_t imageAlt) const;


  int astronomicalCoordinates(Double_t phiWave, Double_t thetaWave, Double_t * RA = 0, Double_t * dec = 0, Double_t * l = 0, Double_t * b = 0) const; 

  int fromRADec (Double_t RA,  Double_t dec, Double_t *phi, Double_t *theta) const; 
  int fromRADec (int N, const Double_t * RA,  const Double_t * dec, Double_t *phi, Double_t *theta) const; 
  

  inline void setDebug(bool db){fDebug = db;}

  inline bool getDebug() const {return fDebug;}

  inline void setInterpSurfaceAboveGeoid(bool i){fInterpSurfaceAboveGeoid = i;} // *TOGGLE *GETTER=getInterpSurfaceAboveGeoid
  inline bool getInterpSurfaceAboveGeoid() const {return fInterpSurfaceAboveGeoid;}

  inline void setSurfaceCloseEnoughInter(double closeEnough = 1.0){fSurfaceCloseEnoughIter = closeEnough;}

  inline double getSurfaceCloseEnoughInter() const {return fSurfaceCloseEnoughIter;}


  inline void setMaxLoopIterations(Int_t n = 50){fMaxLoopIterations = n;}

  inline int getMaxLoopIterations() const {return fMaxLoopIterations;}
  

  inline double surfaceAboveGeoid(Double_t lon, Double_t lat) const {
    if(fInterpSurfaceAboveGeoid){
      return RampdemReader::BilinearInterpolatedSurfaceAboveGeoid(lon, lat);
    }
    else{
      return RampdemReader::SurfaceAboveGeoid(lon, lat);
    }
  }


  void updateCartesianBalloonInfo();

  // mutable std::vector<TGraph*> grTests;

 private:
  Int_t     fIncludeGroupDelay;     
  TVector3  fSourcePos;         
  Double_t  fSourceLongitude;       
  Double_t  fSourceLatitude;        
  Double_t  fSourceAltitude;        
  Double_t  fThetaWave;         
  Double_t  fPhiWave;           
  Double_t  fBalloonCoords[3];      
  TVector3  fBalloonPos;            
  Double_t  fBalloonTheta;          
  Double_t  fBalloonPhi;            
  Double_t  fBalloonHeight;         
  Float_t       fBalloonLonCache;       
  Float_t       fBalloonLatCache;       
  Float_t       fBalloonAltCache;       
  Bool_t    fDebug;             
  Bool_t        fInterpSurfaceAboveGeoid;   
  Double_t      fSurfaceCloseEnoughIter;    
  Int_t         fMaxLoopIterations;             
  ClassDef(UsefulAdu5Pat,0);            

  //optimisation stuff
  Double_t getDeltaTTaylorOpt(Int_t ant1, Int_t ant2, Double_t *deltaR, Double_t *deltaZ, Double_t *deltaPhi);
  Double_t getDeltaTExpectedOpt(Int_t ant1, Int_t ant2,Double_t sourceLon, Double_t sourceLat, Double_t sourceAlt, Double_t *deltaR, Double_t *deltaZ, Double_t *deltaPhi);

  Double_t getDeltaTSeaveyOpt(Int_t ant1, Int_t ant2, Double_t *deltaR, Double_t *deltaZ, Double_t *deltaPhi);
  Double_t getDeltaTExpectedSeaveyOpt(Int_t ant1, Int_t ant2,Double_t sourceLon, Double_t sourceLat, Double_t sourceAlt, Double_t *deltaR, Double_t *deltaZ, Double_t *deltaPhi);
  Double_t getGroupDelay(Double_t phiToAntBoresight);  
  

  void accountForPitchAndRollInPhiWaveThetaWave(Double_t& phiWave, Double_t& thetaWave) const;
  
};


#endif //USEFULADU5PAT_H