| 1 | /*************************************************************************** | 
| 2 | *   Copyright (C) 2006 by pamelaprod                                      * | 
| 3 | *   pamelaprod@P1.pamela                                                  * | 
| 4 | *                                                                         * | 
| 5 | *   This program is free software; you can redistribute it and/or modify  * | 
| 6 | *   it under the terms of the GNU General Public License as published by  * | 
| 7 | *   the Free Software Foundation; either version 2 of the License, or     * | 
| 8 | *   (at your option) any later version.                                   * | 
| 9 | *                                                                         * | 
| 10 | *   This program is distributed in the hope that it will be useful,       * | 
| 11 | *   but WITHOUT ANY WARRANTY; without even the implied warranty of        * | 
| 12 | *   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the         * | 
| 13 | *   GNU General Public License for more details.                          * | 
| 14 | *                                                                         * | 
| 15 | *   You should have received a copy of the GNU General Public License     * | 
| 16 | *   along with this program; if not, write to the                         * | 
| 17 | *   Free Software Foundation, Inc.,                                       * | 
| 18 | *   59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.             * | 
| 19 | ***************************************************************************/ | 
| 20 | #include <InclinationInfo.h> | 
| 21 | #include <TMath.h> | 
| 22 | #include <TMatrixD.h> | 
| 23 |  | 
| 24 | using namespace std; | 
| 25 |  | 
| 26 | // InclinationInfoI()::InclinationInfoI() { | 
| 27 | //   //  memset(time,0,6*sizeof(double)); | 
| 28 | //   // memset(quad,0,6*4*sizeof(double)); | 
| 29 | // }; | 
| 30 |  | 
| 31 | void InclinationInfoI::fill(TArrayC* data){ | 
| 32 | short extIndex = 0; | 
| 33 | short innIndex = 0; | 
| 34 | long tempData = 0; | 
| 35 | for (int i = 0; i < 6; i++){ | 
| 36 | extIndex = 20*i; | 
| 37 | time[i] = (((data->At(extIndex) << 24) & 0xFF000000) + | 
| 38 | ((data->At(extIndex + 1) << 16) & 0x00FF0000) + ((data->At(extIndex + 2) << 8) & 0x0000FF00) + | 
| 39 | (data->At(extIndex + 3) & 0x000000FF))/128.0; | 
| 40 | for (int j = 0; j < 4; j++){ | 
| 41 | innIndex = extIndex + 4*j; | 
| 42 | tempData = ((data->At(innIndex + 4) << 24) & 0xFF000000) + ((data->At(innIndex + 5) << 16) & 0x00FF0000) + ((data->At(innIndex + 6) << 8) & 0x0000FF00) + (data->At(innIndex + 7) & 0x000000FF); | 
| 43 | if (data->At(innIndex + 4) >> 8) { | 
| 44 | quat[i][j] = (~tempData * -1.0)/1073741824.0; | 
| 45 | } else { | 
| 46 | quat[i][j] = tempData / 1073741824.0; | 
| 47 | } | 
| 48 | } | 
| 49 | } | 
| 50 | } | 
| 51 |  | 
| 52 | void InclinationInfoI::clear() { | 
| 53 | for(UInt_t i = 0; i < 6; i++){ | 
| 54 | time[i]=0; | 
| 55 | for(UInt_t j = 0; j < 4; j++) quat[i][j]=0; | 
| 56 | } | 
| 57 | return ; | 
| 58 | } | 
| 59 |  | 
| 60 |  | 
| 61 | Quaternions::Quaternions() | 
| 62 | : InclinationInfoI() | 
| 63 | { | 
| 64 | } | 
| 65 |  | 
| 66 |  | 
| 67 | Quaternions::~Quaternions() | 
| 68 | { | 
| 69 | } | 
| 70 |  | 
| 71 | InclinationInfo::InclinationInfo() | 
| 72 | : TObject() | 
| 73 | { | 
| 74 | } | 
| 75 |  | 
| 76 | InclinationInfo::~InclinationInfo() | 
| 77 | { | 
| 78 | } | 
| 79 |  | 
| 80 | short int Sign_1(double_t a, Int_t b){ | 
| 81 | if(a>0){b=1;} | 
| 82 | if(a<0){b=-1;} | 
| 83 | else{b=0;} | 
| 84 | return b; | 
| 85 | } | 
| 86 |  | 
| 87 |  | 
| 88 | /******************************************************************************************************************/ | 
| 89 | /******************************************************************************************************************/ | 
| 90 | //*********************                             ***************************************************************/ | 
| 91 | //*********************     COORDINATE SYSTEMS      ***************************************************************/ | 
| 92 | //*********************                             ***************************************************************/ | 
| 93 | //*****************************************************************************************************************/ | 
| 94 | //*****************************************************************************************************************/ | 
| 95 | // | 
| 96 | //                                  ZISK | 
| 97 | //                                 + | 
| 98 | //                                / \       YOSK      ZOSK (Directed by Radius) | 
| 99 | //                                 |       _        _. | 
| 100 | //                                 |      |\        /| | 
| 101 | //                                 |        \      / | 
| 102 | //                                 |         \    / | 
| 103 | //                                 |.__..__   \  / | 
| 104 | //                Orbit     _._.***|        **.\/_        XOSK (Directed by velocity) | 
| 105 | //                        .*       | (X0,Y0,Z0) **--.___| | 
| 106 | //                     _**         |        /     *.    / | 
| 107 | //                   .*            |       *        * | 
| 108 | //                  *        ..****|***.. /  R       * | 
| 109 | //                         .*      |    .*. | 
| 110 | //                        .*       |   /  *. | 
| 111 | //                        * EARTH  |  /    *                                    YISK | 
| 112 | //                        *        | /_ _  _*_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _| | 
| 113 | //                        *        /       *                                   / | 
| 114 | //                         *      /      .* | 
| 115 | //                          *.   /      .* | 
| 116 | //                            **/******* | 
| 117 | //                             / | 
| 118 | //                            / | 
| 119 | //                           / | 
| 120 | //                          / | 
| 121 | //                         / | 
| 122 | //                        / | 
| 123 | //                      |/ | 
| 124 | //                      *-- | 
| 125 | //                  XISK | 
| 126 | // | 
| 127 | //****************************************************************************************************/ | 
| 128 | //****************************************************************************************************/ | 
| 129 |  | 
| 130 |  | 
| 131 | void InclinationInfo::TransAngle(Double_t x0, Double_t y0, Double_t z0, Double_t Vx0, Double_t Vy0, Double_t Vz0, Double_t q0, Double_t q1, Double_t q2, Double_t q3){ | 
| 132 |  | 
| 133 | double_t a = 360/(2*TMath::Pi()); | 
| 134 |  | 
| 135 | TMatrixD Xij(3,3); | 
| 136 | Xij(0,0) = 1; Xij(0,1) = 0; Xij(0,2) = 0; | 
| 137 | Xij(1,0) = 0; Xij(1,1) = 0; Xij(1,2) = 1; | 
| 138 | Xij(2,0) = 0; Xij(2,1) = -1; Xij(2,2) = 0; | 
| 139 |  | 
| 140 | TMatrixD Zij(3,3); | 
| 141 | Zij(0,0) = 0; Zij(0,1) = 0; Zij(0,2) = -1; | 
| 142 | Zij(1,0) = -1; Zij(1,1) = 0; Zij(1,2) = 0; | 
| 143 | Zij(2,0) = 0; Zij(2,1) = 1; Zij(2,2) = 0; | 
| 144 |  | 
| 145 | TMatrixD Pij(3,3); | 
| 146 | Pij(0,0) = pow(q0,2)+pow(q1,2)-pow(q2,2)-pow(q3,2); | 
| 147 | Pij(0,1) = /*2*(q1*q2+q0*q3);/*/ 2*(q1*q2-q0*q3); | 
| 148 | Pij(0,2) = /*2*(q1*q3-q0*q2);/*/ 2*(q1*q3+q0*q2); | 
| 149 | Pij(1,0) = /*2*(q1*q2-q0*q3);/*/ 2*(q1*q2+q0*q3); | 
| 150 | Pij(1,1) = pow(q0,2)-pow(q1,2)+pow(q2,2)-pow(q3,2); | 
| 151 | Pij(1,2) = /*2*(q2*q3+q0*q1);/*/ 2*(q2*q3-q0*q1); | 
| 152 | Pij(2,0) = /*2*(q1*q3+q0*q2);/*/ 2*(q1*q3-q0*q2); | 
| 153 | Pij(2,1) = /*2*(q2*q3-q0*q1);/*/ 2*(q2*q3+q0*q1); | 
| 154 | Pij(2,2) = pow(q0,2)-pow(q1,2)-pow(q2,2)+pow(q3,2); | 
| 155 |  | 
| 156 | TMatrixD Aij(3,3); | 
| 157 |  | 
| 158 | Double_t C1 = y0*Vz0 - z0*Vy0; | 
| 159 | Double_t C2 = z0*Vx0 - x0*Vz0; | 
| 160 | Double_t C3 = x0*Vy0 - y0*Vx0; | 
| 161 | Double_t C  = sqrt(pow(C1,2) + pow(C2,2) + pow(C3,2)); | 
| 162 | Double_t V0 = sqrt(pow(Vx0,2)+pow(Vy0,2) + pow(Vz0,2)); | 
| 163 | //    Double_t R0 = sqrt(pow(x0,2)+pow(y0,2) + pow(z0,2)); | 
| 164 | Aij(0,0) = /*(C2*z0-C3*y0)/(C*R0);/*/Vx0/V0; | 
| 165 | Aij(0,1) = C1/C; | 
| 166 | Aij(0,2) = /*x0/R0;/*/(Vy0*C3-Vz0*C2)/(V0*C); | 
| 167 | Aij(1,0) = /*(C3*x0-C1*z0)/(C*R0);/*/Vy0/V0; | 
| 168 | Aij(1,1) = C2/C; | 
| 169 | Aij(1,2) = /*y0/R0;/*/(Vz0*C1-Vx0*C3)/(V0*C); | 
| 170 | Aij(2,0) = /*(C1*y0-C2*x0)/(C*R0);/*/Vz0/V0; | 
| 171 | Aij(2,1) = C3/C; | 
| 172 | Aij(2,2) = /*x0/R0;/*/(Vx0*C2-Vy0*C1)/(V0*C); | 
| 173 | Aij.Invert(); | 
| 174 |  | 
| 175 | TMatrixD Full_(3,3); | 
| 176 |  | 
| 177 | Full_ = Aij*(Pij*Zij); | 
| 178 |  | 
| 179 | //Double_t u13 = Full_(0,2); | 
| 180 | //Double_t u23 = Full_(1,2); | 
| 181 | //Double_t u22 = Full_(1,1); | 
| 182 | //Double_t u33 = Full_(2,2); | 
| 183 | //Double_t u21 = Full_(1,0); | 
| 184 |  | 
| 185 | Double_t u13 = Full_(0,0); | 
| 186 | Double_t u23 = -Full_(1,0); | 
| 187 | Double_t u22 = Full_(1,1); | 
| 188 | Double_t u33 = Full_(2,0); | 
| 189 | Double_t u21 = Full_(1,2); | 
| 190 |  | 
| 191 | Tangazh = a*atan(-u13/u33); | 
| 192 | Kren = a*atan(-u23/sqrt(1 - pow(u23,2))); | 
| 193 | Ryskanie = a*atan(u21/u22); | 
| 194 |  | 
| 195 | return ; | 
| 196 | } | 
| 197 |  | 
| 198 |  | 
| 199 | void InclinationInfo::Clear(){ | 
| 200 | //Int_t gyh = 0; | 
| 201 | } | 
| 202 |  | 
| 203 |  | 
| 204 | //ClassImp(McmdItem) | 
| 205 | ClassImp(InclinationInfoI) | 
| 206 | ClassImp(Quaternions) | 
| 207 | ClassImp(InclinationInfo) |