21 #define ader_ref(a_1,a_2) ader[(a_2)*(vkalNTrkM*3+3) + (a_1) - (vkalNTrkM*3+4)]
22 #define dcv_ref(a_1,a_2) dcv[(a_2)*6 + (a_1) - 7]
23 #define useWeightScheme 1
30 long int ic, jc,
it, jt;
36 long int NVar = (NTRK + 1) * 3;
47 for (
it = 1;
it <= NTRK; ++
it) {
54 for (jt = 1; jt <= NTRK; ++jt) {
55 for (
int k = 0;
k < 3; ++
k) {
56 for (
int l = 0;
l < 3; ++
l) {
58 for (
int j = 0; j < 2; ++j) {
59 dpipj[
k][
l] += vk->
tmpArr[jt-1]->drdp[j][
k] * drdpy[j][
l];
63 for (
int k = 1;
k <= 3; ++
k) {
64 for (
int l = 1;
l <= 3; ++
l) {
84 for (
it = 1;
it <= NTRK; ++
it) {
102 for (
it = 1;
it <= NTRK; ++
it) {
103 for (jt =
it; jt <= NTRK; ++jt) {
120 for (
i=1;
i<=NVar-1; ++
i) {
121 for (j =
i+1; j<=NVar; ++j) {
144 for(
i=0;
i<NVar+1;
i++){ ta[
i] = tab.data() +
i*(NVar+1);}
145 for (
i=1;
i<=NVar; ++
i)
for (j =
i; j<=NVar; ++j) ta[
i][j] = ta[j][
i] =
ader_ref(
i,j);
153 for(
i=0;
i<NVar+1;
i++){ tv[
i] = tvb.data() +
i*(NVar+1); tr[
i] = trb.data() +
i*(NVar+1);}
155 vkSVDCmp( ta.data(), NVar, NVar, tw.data(), tv.data());
158 for(
i=1;
i<NVar+1;
i++)
if(fabs(tw[
i])>tmax)tmax=fabs(tw[
i]);
159 for(
i=1;
i<NVar+1;
i++)
if(fabs(tw[
i])/tmax < 1.
e-18) tw[
i]=0.;
160 for(
i=1;
i<=NVar;
i++){
for(j=1; j<=NVar; j++){
161 tr[
i][j]=0.;
for(
int k=1;
k<=NVar;
k++)
if(tw[
k]!=0.) tr[
i][j] += ta[
i][
k]*tv[j][
k]/tw[
k];
164 for (
i=1;
i<=NVar; ++
i)
for (j=1; j<=NVar; ++j)
ader_ref(
i,j)=tr[
i][j];
170 for(
i=1;
i<=NVar;
i++){
for(j=
i; j<=NVar; j++){
171 double mcheck=0.;
for(
int k=1;
k<=NVar;
k++) mcheck+=ta[
i][
k]*
ader_ref(
k,j);
172 if(
i!=j) maxDiff = (maxDiff > std::abs(mcheck)) ? maxDiff : std::abs(mcheck);
173 if(
i==j) maxDiff = (maxDiff > std::abs(1.-mcheck)) ? maxDiff : std::abs(1.-mcheck);
176 if(maxDiff>0.1)
return -1;
181 for (
i=1;
i<=NVar;
i++) {
182 for (j=1; j<=NVar; j++) {
197 for (
it=1;
it<=NTRK;
it++) {
200 - vcov[1] * t_trk->
wbci[1]
201 - vcov[3] * t_trk->
wbci[2];
203 - vcov[2] * t_trk->
wbci[1]
204 - vcov[4] * t_trk->
wbci[2];
206 - vcov[4] * t_trk->
wbci[1]
207 - vcov[5] * t_trk->
wbci[2];
209 - vcov[1] * t_trk->
wbci[4]
210 - vcov[3] * t_trk->
wbci[5];
212 - vcov[2] * t_trk->
wbci[4]
213 - vcov[4] * t_trk->
wbci[5];
215 - vcov[4] * t_trk->
wbci[4]
216 - vcov[5] * t_trk->
wbci[5];
218 - vcov[1] * t_trk->
wbci[7]
219 - vcov[3] * t_trk->
wbci[8];
221 - vcov[2] * t_trk->
wbci[7]
222 - vcov[4] * t_trk->
wbci[8];
224 - vcov[4] * t_trk->
wbci[7]
225 - vcov[5] * t_trk->
wbci[8];
238 for (
it = 1;
it<=NTRK; ++
it) {
240 for (jt=1; jt<=NTRK; ++jt) {
272 std::vector<std::vector< Vect3DF> > tf0t;
273 std::vector< Vect3DF > th0t;
274 std::vector< double > taa;
275 std::vector< Vect3DF > tmpVec;
285 tf0t.push_back( tmpVec );
289 double **R =
new double*[totNC];
for(
ic=0;
ic<totNC;
ic++) R[
ic]=
new double[NVar];
290 double **RC=
new double*[totNC];
for(
ic=0;
ic<totNC;
ic++)RC[
ic]=
new double[NVar];
291 double *RCRt=
new double[totNC*totNC];
292 for(
ic=0;
ic<totNC;
ic++){
296 for(
it=1;
it<=NTRK;
it++){
303 for(
ic=0;
ic<totNC;
ic++){
304 for(j=0; j<NVar; j++){ RC[
ic][j]=0;
309 for(
ic=0;
ic<totNC;
ic++){
310 for(jc=0; jc<totNC; jc++){ RCRt[
ic*totNC + jc] =0.;
311 for(
i=0;
i<NVar;
i++) RCRt[
ic*totNC + jc] += RC[
ic][
i]*R[jc][
i];
314 dsinv(totNC, RCRt, totNC, &IERR);
315 if ( IERR != 0)
return IERR;
317 for(
i=0;
i<NVar;
i++){
318 for(j=0; j<NVar; j++){
double COR=0.;
319 for(
ic=0;
ic<totNC;
ic++){
320 for(jc=0; jc<totNC; jc++){
321 COR += RC[
ic][
i]*RC[jc][j]*RCRt[
ic*totNC +jc];
328 for(
ic=0;
ic<totNC;
ic++)
delete[] R[
ic];
330 for(
ic=0;
ic<totNC;
ic++)
delete[] RC[
ic];
346 for (
i = 1;
i <= 6; ++
i) {
347 for (j = 1; j <= 6; ++j) {
349 for (
ic=1;
ic<=NVar; ++
ic) {
351 for (jc=1; jc<=NVar; ++jc) {
352 if(
dcv_ref(j, jc)==0.)
continue;
365 #undef useWeightScheme