607 HepMatrix Correlation(dignum,dignum,0), Inverse(dignum,dignum,0),
Zero(dignum,dignum,0);
608 HepMatrix PulseShape(dignum,1,0), DPulseShape(dignum,1,0);
609 HepMatrix
a(dignum,1,0),
b(dignum,1,0), WeightZero(dignum,1,0);
614 log<<MSG::INFO<<
" Calculating OF2 weights ";
616 log<<MSG::INFO<<
" Calculating OF1 weights ";
617 if (deltaCorrelation)
618 log<<MSG::INFO<<
"with Delta correlation matrix "<<
endmsg;
620 log<<MSG::INFO<<
"with correlation matrix obtained from data "<<
endmsg;
626 if (deltaCorrelation){
628 for (
int i=0;
i<dignum;
i++)
629 for (
int j=0;j<dignum;j++)
630 Correlation[
i][j]=
m_R[0][0][0][0][
i][j];
632 Inverse=Correlation.inverse(ierr);
634 for (
int pha=-100;pha<101;pha++){
637 for (
int i=0;
i<dignum;
i++){
639 PulseShape[
i][0]=LshapeForm[
i*25+100+pha];
640 DPulseShape[
i][0]=LdshapeForm[
i*25+100+pha];
642 PulseShape[
i][0]=HshapeForm[
i*25+100+pha];
643 DPulseShape[
i][0]=HdshapeForm[
i*25+100+pha];
645 g[pha+100][
i]=PulseShape[
i][0];
651 HepMatrix SystemMatrix(dignum+2,dignum+2,0);
652 HepVector
Result(dignum + 2,0);
653 HepVector IndependTermsAmp(dignum+2,0);
654 HepVector IndependTermsTime(dignum+2,0);
655 for (
int i=0;
i<dignum;
i++){
657 for (
int j=0; j<dignum; j++){
658 SystemMatrix[
i][j]=Correlation[
i][j];
660 SystemMatrix[dignum][
i] = PulseShape[
i][0];
661 SystemMatrix[
i][dignum] = PulseShape[
i][0];
663 SystemMatrix[dignum + 1][
i] = DPulseShape[
i][0];
664 SystemMatrix[
i][dignum + 1] = DPulseShape[
i][0];
668 IndependTermsAmp[dignum] = 1.;
669 IndependTermsTime[dignum+1] = -1.;
671 Result = solve(SystemMatrix,IndependTermsAmp);
673 for (
int ismp=0; ismp<dignum; ismp++) {
677 Result = solve(SystemMatrix,IndependTermsTime);
679 for (
int ismp=0; ismp<dignum; ismp++){
684 HepMatrix SystemMatrix(dignum+3,dignum+3,0);
685 HepVector
Result(dignum+3,0);
686 HepVector IndependTermsAmp(dignum+3,0);
687 HepVector IndependTermsTime(dignum+3,0);
688 HepVector IndependTermsPed(dignum+3,0);
690 for (
int i=0;
i<dignum;
i++){
691 for (
int j=0; j<dignum; j++){
692 SystemMatrix[
i][j]=Correlation[
i][j];
694 SystemMatrix[dignum][
i] = PulseShape[
i][0];
695 SystemMatrix[
i][dignum] = PulseShape[
i][0];
697 SystemMatrix[dignum + 1][
i] = DPulseShape[
i][0];
698 SystemMatrix[
i][dignum + 1] = DPulseShape[
i][0];
701 SystemMatrix[dignum + 2][
i] = 1.;
702 SystemMatrix[
i][dignum + 2] = 1.;
705 IndependTermsAmp[dignum] = 1.;
706 IndependTermsTime[dignum+1] = -1.;
707 IndependTermsPed[dignum+2] = 1.;
709 Result = solve(SystemMatrix,IndependTermsAmp);
711 for (
int ismp=0; ismp<dignum; ismp++) {
715 Result = solve(SystemMatrix,IndependTermsTime);
717 for (
int ismp=0; ismp<dignum; ismp++){
721 Result = solve(SystemMatrix,IndependTermsPed);
723 for (
int ismp=0; ismp<dignum; ismp++){
729 for (
int ismp=0; ismp<dignum; ismp++){
730 w_a[pha+100][ismp] = 0.;
731 w_b[pha+100][ismp] = 0.;
732 w_c[pha+100][ismp] = 0.;
737 printf(
"ros=%d drawer=%d channel=%d gain=%d N=%d \n",
741 for (
int i=0;
i<dignum;
i++)
742 for (
int j=0;j<dignum;j++)
745 Inverse=Correlation.inverse(ierr);
748 log<<MSG::WARNING<<
" Correlation matrix cannot be inverted for ros="<<
ros<<
" drawer="<<
drawer
750 log<<MSG::WARNING<<
" Weights set to zero for this channel..."<<
endmsg;
751 for (
int i=0;
i<dignum;
i++){
752 for (
int j=0; j<dignum; j++){
759 for (
int pha=-100;pha<101;pha++){
762 for (
int i=0;
i<dignum;
i++){
764 PulseShape[
i][0]=LshapeForm[
i*25+100+pha];
765 DPulseShape[
i][0]=LdshapeForm[
i*25+100+pha];
767 PulseShape[
i][0]=HshapeForm[
i*25+100+pha];
768 DPulseShape[
i][0]=HdshapeForm[
i*25+100+pha];
770 g[pha+100][
i]=PulseShape[
i][0];
775 HepMatrix SystemMatrix(dignum+2,dignum+2,0);
776 HepVector
Result(dignum + 2,0);
777 HepVector IndependTermsAmp(dignum+2,0);
778 HepVector IndependTermsTime(dignum+2,0);
780 for (
int i=0;
i<dignum;
i++){
781 for (
int j=0; j<dignum; j++){
782 SystemMatrix[
i][j]=Correlation[
i][j];
784 SystemMatrix[dignum][
i] = PulseShape[
i][0];
785 SystemMatrix[
i][dignum] = PulseShape[
i][0];
787 SystemMatrix[dignum + 1][
i] = DPulseShape[
i][0];
788 SystemMatrix[
i][dignum + 1] = DPulseShape[
i][0];
791 IndependTermsAmp[dignum] = 1.;
792 IndependTermsTime[dignum+1] = -1.;
794 Result = solve(SystemMatrix,IndependTermsAmp);
796 for (
int ismp=0; ismp<dignum; ismp++) {
800 Result = solve(SystemMatrix,IndependTermsTime);
802 for (
int ismp=0; ismp<dignum; ismp++){
807 HepMatrix SystemMatrix(dignum+3,dignum+3,0);
808 HepVector
Result(dignum + 3,0);
809 HepVector IndependTermsAmp(dignum+3,0);
810 HepVector IndependTermsTime(dignum+3,0);
811 HepVector IndependTermsPed(dignum+3,0);
813 for (
int i=0;
i<dignum;
i++){
814 for (
int j=0; j<dignum; j++){
815 SystemMatrix[
i][j]=Correlation[
i][j];
817 SystemMatrix[dignum][
i] = PulseShape[
i][0];
818 SystemMatrix[
i][dignum] = PulseShape[
i][0];
820 SystemMatrix[dignum + 1][
i] = DPulseShape[
i][0];
821 SystemMatrix[
i][dignum + 1] = DPulseShape[
i][0];
823 SystemMatrix[dignum + 2][
i] = 1.;
824 SystemMatrix[
i][dignum + 2] = 1.;
827 IndependTermsAmp[dignum] = 1.;
828 IndependTermsTime[dignum+1] = -1.;
829 IndependTermsPed[dignum+2] = 1.;
831 Result = solve(SystemMatrix,IndependTermsAmp);
833 for (
int ismp=0; ismp<dignum; ismp++) {
837 Result = solve(SystemMatrix,IndependTermsTime);
839 for (
int ismp=0; ismp<dignum; ismp++){
843 Result = solve(SystemMatrix,IndependTermsPed);
845 for (
int ismp=0; ismp<dignum; ismp++){
850 printf(
"PHASE: %d deter=%f\n",pha,SystemMatrix.determinant());
851 for (
int i=0;
i<dignum+2;
i++){
852 for (
int j=0; j<dignum+2; j++){
853 printf(
"%4.4f ",SystemMatrix[
i][j]);
857 for (
int i=0;
i<dignum;
i++){
858 for (
int j=0; j<dignum; j++){
863 for (
int ismp=0; ismp<dignum; ismp++)
864 printf(
"samp: %d ->w_a: %f w_b:%f\n",
865 ismp,
w_a[pha+100][ismp],
w_b[pha+100][ismp]);
870 for (
int ismp=0; ismp<dignum; ismp++){
871 w_a[pha+100][ismp] = 0.;
872 w_b[pha+100][ismp] = 0.;
873 w_c[pha+100][ismp] = 0.;