/[PAMELA software]/DarthVader/OrbitalInfo/src/InclinationInfo.cpp
ViewVC logotype

Annotation of /DarthVader/OrbitalInfo/src/InclinationInfo.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.8 - (hide annotations) (download)
Tue May 15 14:31:12 2012 UTC (12 years, 6 months ago) by mocchiut
Branch: MAIN
Changes since 1.7: +6 -0 lines
Reprocessing bugs fixed, ToF bugs fixed, new software versions, new quaternions, IGRF bug fixed, code cleanup

1 pamelaprod 1.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    
22     using namespace std;
23    
24     // InclinationInfoI()::InclinationInfoI() {
25     // // memset(time,0,6*sizeof(double));
26     // // memset(quad,0,6*4*sizeof(double));
27     // };
28    
29     void InclinationInfoI::fill(TArrayC* data){
30     short extIndex = 0;
31     short innIndex = 0;
32     long tempData = 0;
33     for (int i = 0; i < 6; i++){
34     extIndex = 20*i;
35     time[i] = (((data->At(extIndex) << 24) & 0xFF000000) +
36     ((data->At(extIndex + 1) << 16) & 0x00FF0000) + ((data->At(extIndex + 2) << 8) & 0x0000FF00) +
37     (data->At(extIndex + 3) & 0x000000FF))/128.0;
38     for (int j = 0; j < 4; j++){
39     innIndex = extIndex + 4*j;
40     tempData = ((data->At(innIndex + 4) << 24) & 0xFF000000) + ((data->At(innIndex + 5) << 16) & 0x00FF0000) + ((data->At(innIndex + 6) << 8) & 0x0000FF00) + (data->At(innIndex + 7) & 0x000000FF);
41     if (data->At(innIndex + 4) >> 8) {
42     quat[i][j] = (~tempData * -1.0)/1073741824.0;
43     } else {
44     quat[i][j] = tempData / 1073741824.0;
45     }
46     }
47     }
48     }
49    
50 mocchiut 1.4 void InclinationInfoI::clear() {
51     for(UInt_t i = 0; i < 6; i++){
52     time[i]=0;
53     for(UInt_t j = 0; j < 4; j++) quat[i][j]=0;
54     }
55     return ;
56     }
57 pamelaprod 1.1
58    
59     Quaternions::Quaternions()
60     : InclinationInfoI()
61     {
62     }
63    
64    
65     Quaternions::~Quaternions()
66     {
67     }
68    
69     InclinationInfo::InclinationInfo()
70     : TObject()
71     {
72     }
73    
74     InclinationInfo::~InclinationInfo()
75     {
76     }
77    
78 mocchiut 1.7 /*Sine::Sine()
79     : TObject()
80     {
81     }
82    
83     Sine::~Sine()
84     {
85     }*/
86    
87 mocchiut 1.4 short int Sign_1(double_t a, Int_t b){
88 pamelaprod 1.1 if(a>0){b=1;}
89     if(a<0){b=-1;}
90     else{b=0;}
91     return b;
92     }
93    
94    
95     /******************************************************************************************************************/
96     /******************************************************************************************************************/
97     //********************* ***************************************************************/
98     //********************* COORDINATE SYSTEMS ***************************************************************/
99     //********************* ***************************************************************/
100     //*****************************************************************************************************************/
101     //*****************************************************************************************************************/
102     //
103     // ZISK
104     // +
105     // / \ YOSK ZOSK (Directed by Radius)
106     // | _ _.
107     // | |\ /|
108     // | \ /
109     // | \ /
110     // |.__..__ \ /
111     // Orbit _._.***| **.\/_ XOSK (Directed by velocity)
112 mocchiut 1.4 // .* | (X0,Y0,Z0) **--.___|
113 pamelaprod 1.1 // _** | / *. /
114     // .* | * *
115     // * ..****|***.. / R *
116     // .* | .*.
117     // .* | / *.
118     // * EARTH | / * YISK
119 mocchiut 1.4 // * | /_ _ _*_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _|
120 pamelaprod 1.1 // * / * /
121     // * / .*
122     // *. / .*
123     // **/*******
124     // /
125     // /
126     // /
127     // /
128     // /
129     // /
130     // |/
131     // *--
132     // XISK
133     //
134     //****************************************************************************************************/
135     //****************************************************************************************************/
136    
137    
138 mocchiut 1.4 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){
139 pamelaprod 1.1
140 mocchiut 1.4 double_t a = 360/(2*TMath::Pi());
141 pamelaprod 1.1
142 mocchiut 1.4 TMatrixD Xij(3,3);
143     Xij(0,0) = 1; Xij(0,1) = 0; Xij(0,2) = 0;
144     Xij(1,0) = 0; Xij(1,1) = 0; Xij(1,2) = 1;
145     Xij(2,0) = 0; Xij(2,1) = -1; Xij(2,2) = 0;
146    
147     TMatrixD Zij(3,3);
148     Zij(0,0) = 0; Zij(0,1) = 0; Zij(0,2) = -1;
149     Zij(1,0) = -1; Zij(1,1) = 0; Zij(1,2) = 0;
150     Zij(2,0) = 0; Zij(2,1) = 1; Zij(2,2) = 0;
151 pamelaprod 1.1
152     TMatrixD Pij(3,3);
153     Pij(0,0) = pow(q0,2)+pow(q1,2)-pow(q2,2)-pow(q3,2);
154 mocchiut 1.4 Pij(0,1) = /*2*(q1*q2+q0*q3);/*/ 2*(q1*q2-q0*q3);
155     Pij(0,2) = /*2*(q1*q3-q0*q2);/*/ 2*(q1*q3+q0*q2);
156     Pij(1,0) = /*2*(q1*q2-q0*q3);/*/ 2*(q1*q2+q0*q3);
157 pamelaprod 1.1 Pij(1,1) = pow(q0,2)-pow(q1,2)+pow(q2,2)-pow(q3,2);
158 mocchiut 1.4 Pij(1,2) = /*2*(q2*q3+q0*q1);/*/ 2*(q2*q3-q0*q1);
159     Pij(2,0) = /*2*(q1*q3+q0*q2);/*/ 2*(q1*q3-q0*q2);
160     Pij(2,1) = /*2*(q2*q3-q0*q1);/*/ 2*(q2*q3+q0*q1);
161 pamelaprod 1.1 Pij(2,2) = pow(q0,2)-pow(q1,2)-pow(q2,2)+pow(q3,2);
162    
163     TMatrixD Aij(3,3);
164 mocchiut 1.4
165 pamelaprod 1.1 Double_t C1 = y0*Vz0 - z0*Vy0;
166     Double_t C2 = z0*Vx0 - x0*Vz0;
167     Double_t C3 = x0*Vy0 - y0*Vx0;
168     Double_t C = sqrt(pow(C1,2) + pow(C2,2) + pow(C3,2));
169     Double_t V0 = sqrt(pow(Vx0,2)+pow(Vy0,2) + pow(Vz0,2));
170 mocchiut 1.4 // Double_t R0 = sqrt(pow(x0,2)+pow(y0,2) + pow(z0,2));
171 pamelaprod 1.1 Aij(0,0) = /*(C2*z0-C3*y0)/(C*R0);/*/Vx0/V0;
172     Aij(0,1) = C1/C;
173     Aij(0,2) = /*x0/R0;/*/(Vy0*C3-Vz0*C2)/(V0*C);
174     Aij(1,0) = /*(C3*x0-C1*z0)/(C*R0);/*/Vy0/V0;
175     Aij(1,1) = C2/C;
176     Aij(1,2) = /*y0/R0;/*/(Vz0*C1-Vx0*C3)/(V0*C);
177     Aij(2,0) = /*(C1*y0-C2*x0)/(C*R0);/*/Vz0/V0;
178     Aij(2,1) = C3/C;
179     Aij(2,2) = /*x0/R0;/*/(Vx0*C2-Vy0*C1)/(V0*C);
180     Aij.Invert();
181    
182 mocchiut 1.4 TMatrixD Full_(3,3);
183 pamelaprod 1.1
184 mocchiut 1.4 Full_ = Aij*(Pij*Zij);
185    
186     //Double_t u13 = Full_(0,2);
187     //Double_t u23 = Full_(1,2);
188     //Double_t u22 = Full_(1,1);
189     //Double_t u33 = Full_(2,2);
190     //Double_t u21 = Full_(1,0);
191    
192     Double_t u13 = Full_(0,0);
193     Double_t u23 = -Full_(1,0);
194     Double_t u22 = Full_(1,1);
195     Double_t u33 = Full_(2,0);
196     Double_t u21 = Full_(1,2);
197 pamelaprod 1.1
198     Tangazh = a*atan(-u13/u33);
199     Kren = a*atan(-u23/sqrt(1 - pow(u23,2)));
200 mocchiut 1.4 Ryskanie = a*atan(u21/u22);
201 pamelaprod 1.1
202 mocchiut 1.8 Full_.Delete();
203     Aij.Delete();
204     Pij.Delete();
205     Zij.Delete();
206     Xij.Delete();
207    
208 pamelaprod 1.1 return ;
209     }
210    
211    
212 mocchiut 1.5 void InclinationInfo::Clear(Option_t *t){
213 mocchiut 1.4 //Int_t gyh = 0;
214     }
215    
216    
217 pamelaprod 1.1 //ClassImp(McmdItem)
218     ClassImp(InclinationInfoI)
219     ClassImp(Quaternions)
220     ClassImp(InclinationInfo)
221 mocchiut 1.7 //ClassImp(Sine)

  ViewVC Help
Powered by ViewVC 1.1.23