Implement Cross-Polarization Antenna Factors - Now clean

Detector.h

@@ -217,7 +217,11 @@ enum EAntennaType {
   eVPol, // (bottom) Vpol
   eVPolTop, // top Vpol
   eHPol, // Hpol
-  eTx // transmitter 
+  eTx, // transmitter
+  eTxCross, // transmitter
+  eVPolCross,     // Cross-pol VPol
+  eVPolTopCross,  // Cross-pol VPol Top
+  eHPolCross      // Cross-pol HPol
 };
 
 class Detector {
@@ -236,6 +240,15 @@ class Detector {
         vector<vector<double> > Hgain;
         vector<vector<double> > Hphase;
         vector<double> Freq;
+        vector<double> FreqCross;
+
+        // Cross-pol vectors
+        vector<vector<double> > VgainCross;
+        vector<vector<double> > VphaseCross;
+        vector<vector<double> > VgainTopCross;
+        vector<vector<double> > VphaseTopCross;
+        vector<vector<double> > HgainCross;
+        vector<vector<double> > HphaseCross;
         double antenna_source_medium_n;    
 
         //Define impedance and gain for receiving antenna
@@ -250,11 +263,14 @@ class Detector {
         //Define impedance and gain for transmitting antenna
         vector<double> RealImpedanceTx;
         vector<double> ImagImpedanceTx;
-        int Tx_freq_init;
-        int Tx_freq_width;
+        double Tx_freq_init;
+        double Tx_freq_width;
         vector<double> TxFreq;
+        vector<double> TxFreqCross;
         vector<vector<double> > Txgain;
         vector<vector<double> > Txphase;
+        vector<vector<double> > TxgainCross;
+        vector<vector<double> > TxphaseCross;
         double Txantenna_source_medium_n;    
         void ReadImpedance(string filename, vector<double> *TempRealImpedance, vector<double> *TempImagImpedance);
         void ReadAllAntennaImpedance(Settings *settings1);
@@ -351,6 +367,11 @@ class Detector {
         vector <double> transVTop_databin;
         vector <double> transH_databin;
 
+        // Cross-pol data bins
+        vector <double> transVCross_databin;
+        vector <double> transVTopCross_databin;
+        vector <double> transHCross_databin;
+
 
         void ReadAmplifierNoiseFigure(Settings *settings1);
         vector< vector<double> > amplifierNoiseFig_ch;
@@ -390,7 +411,7 @@ class Detector {
         double GetGain(double freq, double theta, double phi, int ant_m, int ant_o, double antenna_target_medium_n);    //read antenna gain at certain angle, certain type, and certain orientation
         double GetGain(double freq, double theta, double phi, int ant_m, double antenna_target_medium_n);   //read antenna gain at certain angle, certain type. (orientation : default)
 
-        double GetGain_1D_OutZero(double freq, double theta, double phi, int ant_m, double antenna_target_medium_n, int string_number=0, int ant_number=0, bool useInTransmitterMode=false);   //read antenna gain at certain angle, certain type. (orientation : default) and use 1-D interpolation to get gain, if freq bigger than freq range, return 0 gain
+        double GetGain_1D_OutZero(double freq, double theta, double phi, int ant_m, double antenna_target_medium_n, int string_number=0, int ant_number=0, bool useInTransmitterMode=false, bool useCrossPol=false);   //read antenna gain at certain angle, certain type. (orientation : default) and use 1-D interpolation to get gain, if freq bigger than freq range, return 0 gain
 
         //Creating function to interpolate antenna impedance to frequency binning.
         double GetImpedance(double freq, int ant_m=0, int ant_number=0, bool useInTransmitterMode=false);
@@ -400,7 +421,7 @@ class Detector {
 
         double GetAntPhase(double freq, double theta, double phi, int ant_m, double antenna_target_medium_n); // return antenna phase with 2-D interpolation
 
-        double GetAntPhase_1D(double freq, double theta, double phi, int ant_m, double antenna_target_medium_n, bool useInTransmitterMode=false); // return antenna phase with 1-D interpolation
+        double GetAntPhase_1D(double freq, double theta, double phi, int ant_m, double antenna_target_medium_n, bool useInTransmitterMode=false, bool useCrossPol=false); // return antenna phase with 1-D interpolation
 
 
         double GetFilterGain(int bin) { return FilterGain[bin]; }   // same bin with Vgain, Hgain

Primaries.cc

@@ -2426,11 +2426,30 @@ void Interaction::PickExactGlobal(IceModel *antarctica, Detector *detector, Sett
 
   //Convert source vector into array coordinates
   TVector3 sourceStationVector = AraGeomTool::Instance()->convertArrayToStationCoords(stationId, sourceArrayVector);
-
+
+  // Calculate average antenna depth
+  double avgDepth = 0.;
+  int count = 0;
+  for (int chan = 0; chan < 16; chan++) {
+      double antDepth = AraGeomTool::Instance()->getStationInfo(settings1->DETECTOR_STATION_ARAROOT,settings1->DETECTOR_YEAR)->fAntInfo[chan].antLocation[2];
+      avgDepth += antDepth;
+      count++;
+  }
+
+  avgDepth /= double(count);
+
+
+  // SOURCE_DEPTH is referenced to the surface,
+  // while PickExact() later adds avgZ.
+  // Therefore shift only the Z coordinate.
+  double sourceX = sourceStationVector[0];
+  double sourceY = sourceStationVector[1];
+  double sourceZ = sourceStationVector[2] - avgDepth;
+
   //Calculate posnu
-  double R = sqrt(pow(sourceStationVector[0],2) + pow(sourceStationVector[1],2) + pow(sourceStationVector[2],2));
-  double phi = (360+(atan2(sourceStationVector[1], sourceStationVector[0]))*180/PI)*PI/180;
-  double theta = acos((sourceStationVector[2])/R);
+  double R = sqrt(pow(sourceX,2) + pow(sourceY,2) + pow(sourceZ,2));
+  double phi = (360+(atan2(sourceY, sourceX))*180/PI)*PI/180;
+  double theta = acos((sourceZ)/R);
   settings1->POSNU_THETA = theta;
   settings1->POSNU_PHI = phi;
   settings1->POSNU_R = R;