23 #include "TMatrixDSym.h"
24 #include "TMatrixDSymEigen.h"
72 {
for(
size_t i=0;
i<ifoList.size();
i++) ifoList[
i]->TFmap.Forward(k); }
74 {
for(
size_t i=0;
i<ifoList.size();
i++) ifoList[
i]->TFmap.Inverse(k); }
83 int setTimeShifts(
size_t=1,
double=1.,
size_t=0,
size_t=0,
84 const char* = NULL,
const char* =
"w",
size_t* = NULL);
96 virtual size_t initwc(
double,
double);
104 void setSkyMaps(
double,
double=0.,
double=180.,
double=0.,
double=360.);
108 void setSkyMaps(
int);
116 bool setndm(
size_t,
size_t,
bool=
true,
int=1);
117 bool SETNDM(
size_t,
size_t,
bool=
true,
int=1);
122 inline double getNDM(
size_t i,
size_t j) {
return NDM[
i][
j]; }
131 void setDelayFilters(
detector* = NULL);
132 void setDelayFilters(
char*,
char* = NULL);
136 void writeFilter(
const char *
fname);
139 void readFilter(
const char*);
143 double setVeto(
double=5.);
152 size_t readMDClog(
char*,
double=0.,
int=11,
int=12);
158 size_t readSEGlist(
char*,
int=1);
172 void setDelayIndex(
double rate);
177 void setDelayIndex(
int=0);
184 size_t setIndexMode(
size_t=0);
190 return getifo(0)->tau.getSkyIndex(theta,phi);
202 void delay(
double theta,
double phi);
218 long coherence(
double,
double=0.,
double=0.);
225 long getNetworkPixels(
int LAG,
double Eo,
double DD=1., TH1F*
hist=NULL);
239 size_t netcut(
double,
char=
'L',
size_t=0,
int=1);
248 long subNetCut(
int lag,
float subnet=0.6,
float subcut=0.33, TH2F*
hist=NULL);
266 long likelihoodB(
char=
'E',
double=sqrt(2.),
int=0,
size_t=0,
int=-1,
bool=
false);
273 long likelihoodI(
char=
'P',
double=sqrt(2.),
int=0,
size_t=0,
int=-1,
bool=
false);
283 long likelihood2G(
char mode,
int lag,
int ID, TH2F*
hist=NULL);
284 long likelihoodWP(
char mode,
int lag,
int ID, TH2F*
hist=NULL);
287 long likelihood(
char=
'E',
double=sqrt(2.),
int=0,
size_t=0,
int=-1,
bool=
false);
292 size_t setRank(
double,
double=0.);
304 if(!wc_List.size())
return 0;
306 for(
size_t n=0;
n<nLag;
n++) m += wc_List[
n].
cluster(kt,kf);
312 if(!wc_List.size())
return 0;
314 for(
size_t n=0;
n<nLag;
n++) {
315 m += wc_List[
n].esize();
324 if(!wc_List.size())
return 0;
325 if(lag>=(
int)nLag) lag=nLag-1;
327 if(lag>=0) m += wc_List[lag].esize(type);
328 else for(
size_t n=0;
n<nLag;
n++) m += wc_List[
n].esize(type);
333 inline void setRMS();
337 for(
size_t i=0;
i<nLag;
i++) wc_List[
i].delink();
345 bool getwave(
size_t,
size_t,
char=
'W');
356 bool getMRAwave(
size_t ID,
size_t lag,
char atype=
'S',
int mode=0,
bool tof=
false);
362 void getSkyArea(
size_t id,
size_t lag,
double T);
364 void getSkyArea(
size_t id,
size_t lag,
double T,
double rms);
385 a = 1-
Gamma(ifoList.size(),a*a*ifoList.size());
386 this->e2or =
iGamma(ifoList.size()-1,
a);
400 inline string getmdcList(
size_t n) {
return mdcListSize()>n ? mdcList[
n] :
"\n"; }
402 inline string getmdcType(
size_t n) {
return mdcTypeSize()>n ? mdcType[
n] :
"\n"; }
404 inline std::vector<double>*
getmdcTime() {
return &mdcTime; }
406 inline double getmdcTime(
size_t n) {
return mdcTimeSize()>n ? mdcTime[
n] : 0.; }
408 inline size_t getmdc__ID(
size_t n) {
return mdc__IDSize()>n ? mdc__ID[
n] : 0; }
410 inline double getliveTime(
size_t n) {
return livTimeSize()>n ? livTime[
n] : 0.; }
431 for(
size_t n=0;
n<wdmListSize();
n++)
if(M==(wdmList[
n]->maxLayer()+1))
return wdmList[
n];
446 this->
delta = d==0. ? 0.00001 : d;
453 inline void set2or(
double p) { this->e2or = p; }
458 double threshold(
double,
double);
474 double getDelay(
const char*
c=
"");
483 void setDelay(
const char* =
"L1");
489 void updateTDamp(
int,
float**,
float**);
500 static inline double sumx(
double*);
501 static inline double dotx(
double*,
double*);
502 static inline double dotx(
float*,
float*);
503 static inline double dot4(
double*,
double*);
504 static inline double dotx(
double*,
double*,
double*);
505 static inline double dotx(
float*,
float*,
float*);
506 static inline double dot4(
double*,
double*,
double*);
507 static inline double dotx(
double*,
double**,
size_t);
508 static inline double dotx(
double**,
size_t,
double*);
509 static inline double dotx(
double**,
size_t,
double**,
size_t);
510 static inline double dotx(
double**,
size_t,
double**,
size_t,
double*);
511 static inline double dotx(
double*,
double*,
double);
512 static inline double dotx(
float*,
float*,
float);
513 static inline void addx(
double*,
double*,
double*);
514 static inline void addx(
double*,
double**,
size_t,
double*);
515 static inline void addx(
double**,
size_t,
double**,
size_t,
double*);
516 static inline double dotx(
double*,
double**,
size_t,
double*);
517 static inline double dot32(std::vector<float>*,
double*, std::vector<short>*);
518 static inline double dot32(
double*,
double*,
int*);
519 static inline double divx(
double*,
double*);
520 static inline double rotx(
double*,
double,
double*,
double,
double*);
521 static inline double rotx(
float*,
float,
float*,
float,
float*);
522 static inline double rot4(
double*,
double,
double*,
double,
double*);
523 static inline float rots(
float*,
float,
float*,
float,
float*);
524 static inline void mulx(
double**,
size_t,
double**,
size_t,
double*);
525 static inline void mulx(
double*,
double,
double*);
526 static inline void mulx(
float*,
float,
float*);
527 static inline void mulx(
double*,
double);
528 static inline void mulx(
float*,
float);
529 static inline void inix(
double**,
size_t,
double*);
530 static inline void inix(
double*,
double);
531 static inline void inix(
float*,
float);
532 static inline int netx(
double*,
double,
double*,
double,
double);
533 static inline int netx(
float*,
float,
float*,
float,
float);
534 static inline void pnt_(
float**,
float**,
short**,
int,
int);
535 static inline void cpp_(
float*&,
float**);
536 static inline void cpf_(
float*&
a,
double** p);
537 static inline void cpf_(
float*&
a,
double** p,
size_t);
539 static inline void dpfx(
float* fp,
float* fx);
540 static inline void pnpx(
float* fp,
float* fx,
float* am,
float* AM,
float* u,
float*
v);
541 static inline void dspx(
float* u,
float*
v,
float* am,
float* AM);
542 static inline void dspx(
float* fp,
float* fx,
float* am,
float* AM,
float* u,
float*
v);
544 inline int _sse_MRA_ps(
float*,
float*,
float,
int);
545 inline int _sse_mra_ps(
float*,
float*,
float,
int);
546 inline wavearray<float> _avx_norm_ps(
float**,
float**, std::vector<float*> &,
int);
548 inline void _avx_saveGW_ps(
float**,
float**,
int);
550 void test_sse(
int,
int);
648 void like(
char _LIKE) {this->_LIKE=_LIKE;}
649 void wdm(
bool _WDM) {this->_WDM =_WDM; }
662 size_t n = ifoList.size();
663 if(!ifoList.size() || wc_List.size()!=nLag)
return;
664 for(
size_t i=0;
i<
n;
i++) {
665 for(
size_t j=0;
j<nLag;
j++) {
666 if(!getwc(
j)->size())
continue;
667 if(!getifo(
i)->setrms(getwc(
j),
i)) {
668 cout<<
"network::setRMS() error\n";
692 NETX(d+= a[0]*b[0]; ,
705 NETX(d+= a[0]*b[0]; ,
727 NETX(c[0] = a[0]*b[0]; d+=c[0]; ,
728 c[1] = a[1]*b[1]; d+=c[1]; ,
729 c[2] = a[2]*b[2]; d+=c[2]; ,
730 c[3] = a[3]*b[3]; d+=c[3]; ,
731 c[4] = a[4]*b[4]; d+=c[4]; ,
732 c[5] = a[5]*b[5]; d+=c[5]; ,
733 c[6] = a[6]*b[6]; d+=c[6]; ,
734 c[7] = a[7]*b[7]; d+=c[7]; )
740 NETX(c[0] = a[0]*b[0]; d+=c[0]; ,
741 c[1] = a[1]*b[1]; d+=c[1]; ,
742 c[2] = a[2]*b[2]; d+=c[2]; ,
743 c[3] = a[3]*b[3]; d+=c[3]; ,
744 c[4] = a[4]*b[4]; d+=c[4]; ,
745 c[5] = a[5]*b[5]; d+=c[5]; ,
746 c[6] = a[6]*b[6]; d+=c[6]; ,
747 c[7] = a[7]*b[7]; d+=c[7]; )
753 c[0] = a[0]*b[0]; d+=c[0];
754 c[1] = a[1]*b[1]; d+=c[1];
755 c[2] = a[2]*b[2]; d+=c[2];
756 c[3] = a[3]*b[3]; d+=c[3];
762 NETX(d+= a[0]*b[0][j]; ,
775 NETX(d+= a[0][i]*b[0]; ,
788 NETX(d+= a[0][i]*b[0][j]; ,
789 d+= a[1][
i]*b[1][
j]; ,
790 d+= a[2][
i]*b[2][
j]; ,
791 d+= a[3][
i]*b[3][
j]; ,
792 d+= a[4][
i]*b[4][
j]; ,
793 d+= a[5][
i]*b[5][
j]; ,
794 d+= a[6][
i]*b[6][
j]; ,
795 d+= a[7][
i]*b[7][
j]; )
801 NETX(d+= a[0]/b[0]; ,
814 double q = (1.-
g)*um;
815 NETX(d+=
int(u[0]*u[0]>q) -
int((u[0]*u[0]/um+v[0]*v[0]/vm)>g); ,
816 d+=
int(u[1]*u[1]>q) -
int((u[1]*u[1]/um+v[1]*v[1]/vm)>g); ,
817 d+=
int(u[2]*u[2]>q) -
int((u[2]*u[2]/um+v[2]*v[2]/vm)>g); ,
818 d+=
int(u[3]*u[3]>q) -
int((u[3]*u[3]/um+v[3]*v[3]/vm)>g); ,
819 d+=
int(u[4]*u[4]>q) -
int((u[4]*u[4]/um+v[4]*v[4]/vm)>g); ,
820 d+=
int(u[5]*u[5]>q) -
int((u[5]*u[5]/um+v[5]*v[5]/vm)>g); ,
821 d+=
int(u[6]*u[6]>q) -
int((u[6]*u[6]/um+v[6]*v[6]/vm)>g); ,
822 d+=
int(u[7]*u[7]>q) -
int((u[7]*u[7]/um+v[7]*v[7]/vm)>g); )
829 NETX(d+=
int(u[0]*u[0]>q) -
int((u[0]*u[0]/um+v[0]*v[0]/vm)>g); ,
830 d+=
int(u[1]*u[1]>q) -
int((u[1]*u[1]/um+v[1]*v[1]/vm)>g); ,
831 d+=
int(u[2]*u[2]>q) -
int((u[2]*u[2]/um+v[2]*v[2]/vm)>g); ,
832 d+=
int(u[3]*u[3]>q) -
int((u[3]*u[3]/um+v[3]*v[3]/vm)>g); ,
833 d+=
int(u[4]*u[4]>q) -
int((u[4]*u[4]/um+v[4]*v[4]/vm)>g); ,
834 d+=
int(u[5]*u[5]>q) -
int((u[5]*u[5]/um+v[5]*v[5]/vm)>g); ,
835 d+=
int(u[6]*u[6]>q) -
int((u[6]*u[6]/um+v[6]*v[6]/vm)>g); ,
836 d+=
int(u[7]*u[7]>q) -
int((u[7]*u[7]/um+v[7]*v[7]/vm)>g); )
842 NETX(p[0] = a[0][i]*b[0][j];d+=p[0]; ,
843 p[1] = a[1][
i]*b[1][
j];d+=p[1]; ,
844 p[2] = a[2][
i]*b[2][
j];d+=p[2]; ,
845 p[3] = a[3][
i]*b[3][
j];d+=p[3]; ,
846 p[4] = a[4][
i]*b[4][
j];d+=p[4]; ,
847 p[5] = a[5][
i]*b[5][
j];d+=p[5]; ,
848 p[6] = a[6][
i]*b[6][
j];d+=p[6]; ,
849 p[7] = a[7][
i]*b[7][
j];d+=p[7]; )
855 NETX(d+= a[0]*b[0]; ,
868 NETX(d+= a[0]*b[0]; ,
881 NETX(p[0] = a[0]*b[0][j];d+=p[0]; ,
882 p[1] = a[1]*b[1][
j];d+=p[1]; ,
883 p[2] = a[2]*b[2][
j];d+=p[2]; ,
884 p[3] = a[3]*b[3][
j];d+=p[3]; ,
885 p[4] = a[4]*b[4][
j];d+=p[4]; ,
886 p[5] = a[5]*b[5][
j];d+=p[5]; ,
887 p[6] = a[6]*b[6][
j];d+=p[6]; ,
888 p[7] = a[7]*b[7][
j];d+=p[7]; )
893 NETX(p[0] += a[0]*b[0]; ,
905 NETX(p[0] += a[0]*b[0][j]; ,
906 p[1] += a[1]*b[1][
j]; ,
907 p[2] += a[2]*b[2][
j]; ,
908 p[3] += a[3]*b[3][
j]; ,
909 p[4] += a[4]*b[4][
j]; ,
910 p[5] += a[5]*b[5][
j]; ,
911 p[6] += a[6]*b[6][
j]; ,
912 p[7] += a[7]*b[7][
j]; )
917 NETX(p[0] += a[0][i]*b[0][j]; ,
918 p[1] += a[1][
i]*b[1][
j]; ,
919 p[2] += a[2][
i]*b[2][
j]; ,
920 p[3] += a[3][
i]*b[3][
j]; ,
921 p[4] += a[4][
i]*b[4][
j]; ,
922 p[5] += a[5][
i]*b[5][
j]; ,
923 p[6] += a[6][
i]*b[6][
j]; ,
924 p[7] += a[7][
i]*b[7][
j]; )
929 NETX(p[0] = a[0][i]*b[0][j]; ,
930 p[1] = a[1][
i]*b[1][
j]; ,
931 p[2] = a[2][
i]*b[2][
j]; ,
932 p[3] = a[3][
i]*b[3][
j]; ,
933 p[4] = a[4][
i]*b[4][
j]; ,
934 p[5] = a[5][
i]*b[5][
j]; ,
935 p[6] = a[6][
i]*b[6][
j]; ,
936 p[7] = a[7][
i]*b[7][
j]; )
941 NETX(p[0] = a[0]*b; ,
953 NETX(p[0] = a[0]*b; ,
989 NETX(p[0] = a[0][j]; ,
1026 NETX(e[0] = u[0]*c + v[0]*s;d+=e[0]*e[0]; ,
1027 e[1] = u[1]*c + v[1]*
s;d+=e[1]*e[1]; ,
1028 e[2] = u[2]*c + v[2]*
s;d+=e[2]*e[2]; ,
1029 e[3] = u[3]*c + v[3]*
s;d+=e[3]*e[3]; ,
1030 e[4] = u[4]*c + v[4]*
s;d+=e[4]*e[4]; ,
1031 e[5] = u[5]*c + v[5]*
s;d+=e[5]*e[5]; ,
1032 e[6] = u[6]*c + v[6]*
s;d+=e[6]*e[6]; ,
1033 e[7] = u[7]*c + v[7]*
s;d+=e[7]*e[7]; )
1039 NETX(e[0] = u[0]*c + v[0]*s;d+=e[0]*e[0]; ,
1040 e[1] = u[1]*c + v[1]*
s;d+=e[1]*e[1]; ,
1041 e[2] = u[2]*c + v[2]*
s;d+=e[2]*e[2]; ,
1042 e[3] = u[3]*c + v[3]*
s;d+=e[3]*e[3]; ,
1043 e[4] = u[4]*c + v[4]*
s;d+=e[4]*e[4]; ,
1044 e[5] = u[5]*c + v[5]*
s;d+=e[5]*e[5]; ,
1045 e[6] = u[6]*c + v[6]*
s;d+=e[6]*e[6]; ,
1046 e[7] = u[7]*c + v[7]*
s;d+=e[7]*e[7]; )
1052 e[0] = u[0]*c + v[0]*
s;d+=e[0]*e[0];
1053 e[1] = u[1]*c + v[1]*
s;d+=e[1]*e[1];
1054 e[2] = u[2]*c + v[2]*
s;d+=e[2]*e[2];
1055 e[3] = u[3]*c + v[3]*
s;d+=e[3]*e[3];
1061 NETX(e[0] -= u[0]*c + v[0]*s; d+=e[0]*e[0]; ,
1062 e[1] -= u[1]*c + v[1]*
s; d+=e[1]*e[1]; ,
1063 e[2] -= u[2]*c + v[2]*
s; d+=e[2]*e[2]; ,
1064 e[3] -= u[3]*c + v[3]*
s; d+=e[3]*e[3]; ,
1065 e[4] -= u[4]*c + v[4]*
s; d+=e[4]*e[4]; ,
1066 e[5] -= u[5]*c + v[5]*
s; d+=e[5]*e[5]; ,
1067 e[6] -= u[6]*c + v[6]*
s; d+=e[6]*e[6]; ,
1068 e[7] -= u[7]*c + v[7]*
s; d+=e[7]*e[7]; )
1075 NETX(q[0] = (p[0] + m[0][l]*n);,
1076 q[1] = (p[1] + m[1][
l]*
n);,
1077 q[2] = (p[2] + m[2][
l]*
n);,
1078 q[3] = (p[3] + m[3][
l]*
n);,
1079 q[4] = (p[4] + m[4][
l]*
n);,
1080 q[5] = (p[5] + m[5][
l]*
n);,
1081 q[6] = (p[6] + m[6][
l]*
n);,
1082 q[7] = (p[7] + m[7][
l]*
n);)
1088 for(
int i=0;
i<
XIFO;
i++) *(a++) = *p[
i];
1095 for(
int k=0;
k<
XIFO;
k++) *(a++) = p[
k][
i];
1102 for(
int i=0;
i<
XIFO;
i++) *(a++) = *p[
i]++;
1108 return (*F)[0] *p[(*J)[0]] + (*F)[1] *p[(*J)[1]]
1109 + (*F)[2] *p[(*J)[2]] + (*F)[3] *p[(*J)[3]]
1110 + (*F)[4] *p[(*J)[4]] + (*F)[5] *p[(*J)[5]]
1111 + (*F)[6] *p[(*J)[6]] + (*F)[7] *p[(*J)[7]]
1112 + (*F)[8] *p[(*J)[8]] + (*F)[9] *p[(*J)[9]]
1113 + (*F)[10]*p[(*J)[10]] + (*F)[11]*p[(*J)[11]]
1114 + (*F)[12]*p[(*J)[12]] + (*F)[13]*p[(*J)[13]]
1115 + (*F)[14]*p[(*J)[14]] + (*F)[15]*p[(*J)[15]]
1116 + (*F)[16]*p[(*J)[16]] + (*F)[17]*p[(*J)[17]]
1117 + (*F)[18]*p[(*J)[18]] + (*F)[19]*p[(*J)[19]]
1118 + (*F)[20]*p[(*J)[20]] + (*F)[21]*p[(*J)[21]]
1119 + (*F)[22]*p[(*J)[22]] + (*F)[23]*p[(*J)[23]]
1120 + (*F)[24]*p[(*J)[24]] + (*F)[25]*p[(*J)[25]]
1121 + (*F)[26]*p[(*J)[26]] + (*F)[27]*p[(*J)[27]]
1122 + (*F)[28]*p[(*J)[28]] + (*F)[29]*p[(*J)[29]]
1123 + (*F)[30]*p[(*J)[30]] + (*F)[31]*p[(*J)[31]];
1127 return F[0] *p[J[0]] + F[1] *p[J[1]] + F[2] *p[J[2]] + F[3] *p[J[3]]
1128 + F[4] *p[J[4]] + F[5] *p[J[5]] + F[6] *p[J[6]] + F[7] *p[J[7]]
1129 + F[8] *p[J[8]] + F[9] *p[J[9]] + F[10]*p[J[10]] + F[11]*p[J[11]]
1130 + F[12]*p[J[12]] + F[13]*p[J[13]] + F[14]*p[J[14]] + F[15]*p[J[15]]
1131 + F[16]*p[J[16]] + F[17]*p[J[17]] + F[18]*p[J[18]] + F[19]*p[J[19]]
1132 + F[20]*p[J[20]] + F[21]*p[J[21]] + F[22]*p[J[22]] + F[23]*p[J[23]]
1133 + F[24]*p[J[24]] + F[25]*p[J[25]] + F[26]*p[J[26]] + F[27]*p[J[27]]
1134 + F[28]*p[J[28]] + F[29]*p[J[29]] + F[30]*p[J[30]] + F[31]*p[J[31]];
1144 __m128* _t = (__m128*) t; __m128* _T = (__m128*) T;
1145 __m128* _fp = (__m128*) fp; __m128* _fx = (__m128*) fx;
1153 inline void network::pnpx(
float* fp,
float* fx,
float* am,
float* AM,
float* u,
float*
v) {
1158 __m128* _fp = (__m128*) fp; __m128* _fx = (__m128*) fx;
1159 __m128* _am = (__m128*) am; __m128* _AM = (__m128*) AM;
1160 __m128* _u = (__m128*) u; __m128* _v = (__m128*) v;
1175 __m128* _am = (__m128*) am; __m128* _AM = (__m128*) AM;
1176 __m128* _u = (__m128*) u; __m128* _v = (__m128*) v;
1177 __m128* _U = (__m128*) U; __m128* _V = (__m128*) V;
1185 inline void network::dspx(
float* fp,
float* fx,
float* am,
float* AM,
float* u,
float*
v) {
1190 pnpx(fp, fx, am, AM, u, v);
1208 int V = (
int)this->rNRG.size();
1209 float* ee = this->rNRG.data;
1210 float* pp = this->pNRG.data;
1216 for(j=0; j<V; ++j) if(ee[j]>Eo) pp[j]=0;
1218 __m128* _m00 = (__m128*) mam;
1219 __m128* _m90 = (__m128*) mAM;
1220 __m128* _amp = (__m128*) amp;
1221 __m128* _AMP = (__m128*) AMP;
1222 __m128* _a00 = (__m128*) a_00.data;
1223 __m128* _a90 = (__m128*) a_90.data;
1227 for(j=0; j<V; ++j) if(ee[j]>ee[
m]) m=j;
1228 if(ee[m]<=Eo)
break; mm = m*
f;
1231 int J = wdmMRA.size(m);
1232 float*
c = wdmMRA.getXTalk(m);
1234 if(E/EE < 0.01)
break;
1241 for(j=0; j<J; j++) {
1272 int V = (
int)this->rNRG.size();
1273 float* ee = this->rNRG.data;
1274 float* pp = this->pNRG.data;
1281 __m128* _m00 = (__m128*) mam;
1282 __m128* _m90 = (__m128*) mAM;
1283 __m128* _amp = (__m128*) amp;
1284 __m128* _AMP = (__m128*) AMP;
1285 __m128* _a00 = (__m128*) a_00.data;
1286 __m128* _a90 = (__m128*) a_90.data;
1292 for(j=0; j<V; ++j) if(ee[j]>ee[
m]) m=j;
1293 if(ee[m]<=Eo)
break; mm = m*
f;
1303 c = wdmMRA.getXTalk(m);
1304 for(j=0; j<J; j++) {
1310 if(ee[n]>0) cc += c[1];
1325 c = wdmMRA.getXTalk(m);
1326 for(j=0; j<J; j++) {
1328 if(E<E2 && n!=m) {c+=8;
continue;}
1348 std::vector<float*> &pAVX,
int I) {
1350 float*
g = norm.
data+1; norm=0.;
1359 float* mk = pAVX[1];
1360 float* rn = pAVX[22];
1362 float*
t = tmp.
data; tmp=0.;
1365 float am[4*8]
_ALIGNED; __m128* _am = (__m128*)am;
1366 float x[4*8]
_ALIGNED; __m128* _x = (__m128*)
x;
1367 float h[4*8]
_ALIGNED; __m128* _h = (__m128*)
h;
1369 float* an = pAVX[17];
1370 NETX(__m128* _a0 = (__m128*)(an+4*M*0);, __m128* _a1 = (__m128*)(an+4*M*1);,
1371 __m128* _a2 = (__m128*)(an+4*M*2);, __m128* _a3 = (__m128*)(an+4*M*3);,
1372 __m128* _a4 = (__m128*)(an+4*M*4);, __m128* _a5 = (__m128*)(an+4*M*5);,
1373 __m128* _a6 = (__m128*)(an+4*M*6);, __m128* _a7 = (__m128*)(an+4*M*7);)
1375 for(m=0; m<
M; m++) {
1377 NETX(_a0[m] = _mm_set_ps(q[0][m],q[0][m],p[0][m],p[0][m]); q[0][M+
m]=0; ,
1378 _a1[
m] = _mm_set_ps(q[1][m],q[1][m],p[1][m],p[1][m]); q[1][M+
m]=0; ,
1379 _a2[
m] = _mm_set_ps(q[2][m],q[2][m],p[2][m],p[2][m]); q[2][M+
m]=0; ,
1380 _a3[
m] = _mm_set_ps(q[3][m],q[3][m],p[3][m],p[3][m]); q[3][M+
m]=0; ,
1381 _a4[
m] = _mm_set_ps(q[4][m],q[4][m],p[4][m],p[4][m]); q[4][M+
m]=0; ,
1382 _a5[
m] = _mm_set_ps(q[5][m],q[5][m],p[5][m],p[5][m]); q[5][M+
m]=0; ,
1383 _a6[
m] = _mm_set_ps(q[6][m],q[6][m],p[6][m],p[6][m]); q[6][M+
m]=0; ,
1384 _a7[
m] = _mm_set_ps(q[7][m],q[7][m],p[7][m],p[7][m]); q[7][M+
m]=0; )
1387 for(m=0; m<
M; m++) {
1388 if(mk[m]<=0.)
continue;
1390 int J = wdmMRA.size(m)*2;
1392 float*
c = wdmMRA.getXTalk(m);
1393 __m128* _c = (__m128*)(c+4);
1395 NETX(u=p[0][m]; v=q[0][
m]; _am[0]=_mm_set_ps(v,u,v,u); _x[0]=_mm_setzero_ps(); ,
1396 u=p[1][
m]; v=q[1][
m]; _am[1]=_mm_set_ps(v,u,v,u); _x[1]=_mm_setzero_ps(); ,
1397 u=p[2][
m]; v=q[2][
m]; _am[2]=_mm_set_ps(v,u,v,u); _x[2]=_mm_setzero_ps(); ,
1398 u=p[3][
m]; v=q[3][
m]; _am[3]=_mm_set_ps(v,u,v,u); _x[3]=_mm_setzero_ps(); ,
1399 u=p[4][
m]; v=q[4][
m]; _am[4]=_mm_set_ps(v,u,v,u); _x[4]=_mm_setzero_ps(); ,
1400 u=p[5][
m]; v=q[5][
m]; _am[5]=_mm_set_ps(v,u,v,u); _x[5]=_mm_setzero_ps(); ,
1401 u=p[6][
m]; v=q[6][
m]; _am[6]=_mm_set_ps(v,u,v,u); _x[6]=_mm_setzero_ps(); ,
1402 u=p[7][
m]; v=q[7][
m]; _am[7]=_mm_set_ps(v,u,v,u); _x[7]=_mm_setzero_ps(); )
1404 for(j=0; j<J; j+=2) {
1406 NETX(_x[0]=_mm_add_ps(_x[0],_mm_mul_ps(_c[j],_a0[n]));,
1407 _x[1]=_mm_add_ps(_x[1],_mm_mul_ps(_c[j],_a1[n]));,
1408 _x[2]=_mm_add_ps(_x[2],_mm_mul_ps(_c[j],_a2[n]));,
1409 _x[3]=_mm_add_ps(_x[3],_mm_mul_ps(_c[j],_a3[n]));,
1410 _x[4]=_mm_add_ps(_x[4],_mm_mul_ps(_c[j],_a4[n]));,
1411 _x[5]=_mm_add_ps(_x[5],_mm_mul_ps(_c[j],_a5[n]));,
1412 _x[6]=_mm_add_ps(_x[6],_mm_mul_ps(_c[j],_a6[n]));,
1413 _x[7]=_mm_add_ps(_x[7],_mm_mul_ps(_c[j],_a7[n]));)
1416 NETX(_h[0]=_mm_mul_ps(_x[0],_am[0]);,
1417 _h[1]=_mm_mul_ps(_x[1],_am[1]);,
1418 _h[2]=_mm_mul_ps(_x[2],_am[2]);,
1419 _h[3]=_mm_mul_ps(_x[3],_am[3]);,
1420 _h[4]=_mm_mul_ps(_x[4],_am[4]);,
1421 _h[5]=_mm_mul_ps(_x[5],_am[5]);,
1422 _h[6]=_mm_mul_ps(_x[6],_am[6]);,
1423 _h[7]=_mm_mul_ps(_x[7],_am[7]);)
1425 NETX(t[0]=
h[ 0]+
h[ 1]+
h[ 2]+
h[ 3]; t[0]=t[0]>0?t[0]:0; g[0]+=t[0]; ,
1426 t[1]=
h[ 4]+
h[ 5]+
h[ 6]+
h[ 7]; t[1]=t[1]>0?t[1]:0; g[1]+=t[1]; ,
1427 t[2]=h[ 8]+h[ 9]+h[10]+h[11]; t[2]=t[2]>0?t[2]:0; g[2]+=t[2]; ,
1428 t[3]=h[12]+h[13]+h[14]+h[15]; t[3]=t[3]>0?t[3]:0; g[3]+=t[3]; ,
1429 t[4]=h[16]+h[17]+h[18]+h[19]; t[4]=t[4]>0?t[4]:0; g[4]+=t[4]; ,
1430 t[5]=h[20]+h[21]+h[22]+h[23]; t[5]=t[5]>0?t[5]:0; g[5]+=t[5]; ,
1431 t[6]=h[24]+h[25]+h[26]+h[27]; t[6]=t[6]>0?t[6]:0; g[6]+=t[6]; ,
1432 t[7]=h[28]+h[29]+h[30]+h[31]; t[7]=t[7]>0?t[7]:0; g[7]+=t[7]; )
1436 NETX(u=p[0][m]; v=q[0][
m]; e=(u*u+v*
v)/(t[0]+o); q[0][M+
m]=(e>=1)?0:e; ,
1437 u=p[1][
m]; v=q[1][
m]; e=(u*u+v*
v)/(t[1]+o); q[1][M+
m]=(e>=1)?0:e; ,
1438 u=p[2][
m]; v=q[2][
m]; e=(u*u+v*
v)/(t[2]+o); q[2][M+
m]=(e>=1)?0:e; ,
1439 u=p[3][
m]; v=q[3][
m]; e=(u*u+v*
v)/(t[3]+o); q[3][M+
m]=(e>=1)?0:e; ,
1440 u=p[4][
m]; v=q[4][
m]; e=(u*u+v*
v)/(t[4]+o); q[4][M+
m]=(e>=1)?0:e; ,
1441 u=p[5][
m]; v=q[5][
m]; e=(u*u+v*
v)/(t[5]+o); q[5][M+
m]=(e>=1)?0:e; ,
1442 u=p[6][
m]; v=q[6][
m]; e=(u*u+v*
v)/(t[6]+o); q[6][M+
m]=(e>=1)?0:e; ,
1443 u=p[7][
m]; v=q[7][
m]; e=(u*u+v*
v)/(t[7]+o); q[7][M+
m]=(e>=1)?0:e; )
1445 NETX(u=
x[ 0]+
x[ 2]; v=
x[ 1]+
x[ 3]; rn[
m]+=u*u+v*
v; ,
1446 u=x[ 4]+x[ 6]; v=x[ 5]+x[ 7]; rn[
m]+=u*u+v*
v; ,
1447 u=x[ 8]+x[10]; v=x[ 9]+x[11]; rn[
m]+=u*u+v*
v; ,
1448 u=x[12]+x[14]; v=x[13]+x[15]; rn[
m]+=u*u+v*
v; ,
1449 u=x[16]+x[18]; v=x[17]+x[19]; rn[
m]+=u*u+v*
v; ,
1450 u=x[20]+x[22]; v=x[21]+x[23]; rn[
m]+=u*u+v*
v; ,
1451 u=x[24]+x[26]; v=x[25]+x[27]; rn[
m]+=u*u+v*
v; ,
1452 u=x[28]+x[30]; v=x[29]+x[31]; rn[
m]+=u*u+v*
v; )
1456 for(n=1; n<=
XIFO; n++) {
1458 e = q[n-1][II+4]*q[n-1][II+4];
1460 q[n-1][II+5] = norm.
data[
n];
1476 float* nn = norm.
data;
1478 for(
int n=1; n<=
XIFO; n++) {
1479 nn[
n] = q[n-1][II+5];
1480 p[n-1][II+5] = nn[
n];
1481 e = p[n-1][II+4]*p[n-1][II+4];
1484 for(
int i=0; i<I; i++) {
1485 p[n-1][I+
i] = ec[
i]>0 ? q[n-1][I+
i] : 0.;
1497 for(
int i=0; i<I; i++) {
1498 for(
int n=0; n<
NIFO; n++) {
1499 a_00[i*NIFO+
n]=p[
n][
i];
1500 a_90[i*NIFO+
n]=q[
n][
i];
1506 #endif // NETWORK_HH
wavearray< double > t(hp.size())
monster wdmMRA
list of pixel pointers for MRA
std::vector< char * > ifoName
detector * getifo(size_t n)
param: detector index
static float _sse_abs_ps(__m128 *_a)
static double g(double e)
wavearray< double > skyHole
void like(char _LIKE)
buffer for standard response 90 ampl
WSeries< double > pixeLHood
std::vector< netcluster > wc_List
void set2or(double p)
param: threshold
std::vector< double > * getmdcTime()
static void cpp_(float *&, float **)
static float rots(float *, float, float *, float, float *)
std::vector< pixel > cluster
wavearray< double > a(hp.size())
size_t setSkyMask(double, char *, bool *, int)
static void mulx(double **, size_t, double **, size_t, double *)
static void inix(double **, size_t, double *)
static void pnpx(float *fp, float *fx, float *am, float *AM, float *u, float *v)
static void _sse_zero_ps(__m128 *_p)
std::vector< vectorD > NDM
wavearray< float > a_90
buffer for cluster sky 00 amplitude
double iGamma(double r, double p)
static void _sse_dsp4_ps(__m128 *u, __m128 *v, __m128 *_am, __m128 *_AM, __m128 *_u, __m128 *_v)
double getliveTime(size_t n)
double getmdcTime(size_t n)
cout<< endl;cout<< "ts size = "<< ts.size()<< " ts rate = "<< ts.rate()<< endl;tf.Forward(ts, wdm);int levels=tf.getLevel();cout<< "tf size = "<< tf.size()<< endl;double dF=tf.resolution();double dT=1./(2 *dF);cout<< "rate(hz) : "<< RATE<< "\t layers : "<< nLAYERS<< "\t dF(hz) : "<< dF<< "\t dT(ms) : "<< dT *1000.<< endl;int itime=TIME_PIXEL_INDEX;int ifreq=FREQ_PIXEL_INDEX;int index=(levels+1)*itime+ifreq;double time=itime *dT;double freq=(ifreq >0)?ifreq *dF:dF/4;cout<< endl;cout<< "PIXEL TIME = "<< time<< " sec "<< endl;cout<< "PIXEL FREQ = "<< freq<< " Hz "<< endl;cout<< endl;wavearray< double > x
static void _sse_cpf4_ps(__m128 *_aa, __m128 *_pp)
int _sse_mra_ps(float *, float *, float, int)
WDM< double > wdm(nLAYERS, nLAYERS, 6, 10)
std::vector< netpixel * > pList
string getmdcType(size_t n)
static float _sse_rotsub_ps(__m128 *_u, float c, __m128 *_v, float s, __m128 *_a)
void _avx_saveGW_ps(float **, float **, int)
std::vector< std::string > mdcType
string getmdcList(size_t n)
static void _sse_cpf_ps(float *a, __m128 *_p)
std::vector< size_t > mdc__ID
static double rotx(double *, double, double *, double, double *)
std::vector< double > mdcTime
size_t cluster(int kt, int kf)
param: time gap in pixels return: number of reconstructed clusters
void setRunID(size_t n)
param: run
std::vector< double > vectorD
std::vector< detector * > ifoList
static void _sse_rotadd_ps(__m128 *_u, float c, __m128 *_v, float s, __m128 *_a)
long subNetCut(network *net, int lag, float snc, TH2F *hist)
SUPERCLUSTER.
wavearray< double > skyENRG
size_t events(int type, int lag=-1)
virtual void Browse(TBrowser *)
void constraint(double d=1., double g=0.0001)
param: constraint parameter, p=0 - no constraint
WDM< double > * getwdm(size_t M)
param: number of wdm layers
static int netx(double *, double, double *, double, double)
wavearray< float > pNRG
buffers for cluster residual energy
std::vector< delayFilter > filter90
static double divx(double *, double *)
void Forward(size_t k)
param: number of steps
int _sse_MRA_ps(network *net, float *amp, float *AMP, float Eo, int K)
static void dspx(float *u, float *v, float *am, float *AM)
std::vector< double > livTime
static void addx(double *, double *, double *)
std::vector< std::string > mdcList
static void dpfx(float *fp, float *fx)
wavearray< float > a_00
wdm multi-resolution analysis
std::vector< WDM< double > * > wdmList
static void cpf_(float *&a, double **p)
static void pnt_(float **, float **, short **, int, int)
void setOffset(double t)
param: run
std::vector< waveSegment > segList
int _sse_mra_ps(network *NET, float *amp, float *AMP, float Eo, int K)
static double dot32(std::vector< float > *, double *, std::vector< short > *)
int getIndex(double theta, double phi)
param: theta [deg] param: phi [deg]
wavearray< double > skyProb
static double rot4(double *, double, double *, double, double *)
std::vector< delayFilter > filter
WSeries< double > pixeLNull
static void _sse_add_ps(__m128 *_a, __m128 *_b)
static double dot4(double *, double *)
double getNDM(size_t i, size_t j)
param: first detector param: second detector
wavearray< float > _avx_norm_ps(float **, float **, std::vector< float * > &, int)
size_t add(WDM< double > *wdm)
param: pointer to wdm return number of wdm tronsforms in the list
wavearray< double > skyMaskCC
netcluster * getwc(size_t n)
param: delay index
skymap nSensitivity
list of wdm tranformations
int _sse_MRA_ps(float *, float *, float, int)
static double dotx(double *, double *)
static void _sse_mul_ps(__m128 *_a, float b)
wavearray< float > rNRG
buffer for cluster sky 90 amplitudes
void setMRAcatalog(char *fn)
wavearray< short > skyMask
size_t getmdc__ID(size_t n)
static void _sse_dpf4_ps(__m128 *_Fp, __m128 *_Fx, __m128 *_fp, __m128 *_fx)
static double sumx(double *)
static void _sse_pnp4_ps(__m128 *_fp, __m128 *_fx, __m128 *_am, __m128 *_AM, __m128 *_u, __m128 *_v)