/[PAMELA software]/DarthVader/CalorimeterLevel2/src/CaloLevel0.cpp
ViewVC logotype

Diff of /DarthVader/CalorimeterLevel2/src/CaloLevel0.cpp

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1.15 by mocchiut, Mon Dec 3 20:55:58 2007 UTC revision 1.20 by mocchiut, Fri Apr 18 18:55:47 2008 UTC
# Line 728  Int_t CaloLevel0::CalcCrossTalkCorr(GL_T Line 728  Int_t CaloLevel0::CalcCrossTalkCorr(GL_T
728    return(0);    return(0);
729  }  }
730    
731    void CaloLevel0::FindBaseCompress(Int_t l, Int_t m, Int_t pre){
732      Int_t n = 0;
733      Float_t q = 0;
734      this->FindBaseCompress(l,m,pre,n,q);
735    }
736    
737    void CaloLevel0::FindBaseCompress(Int_t l, Int_t m, Int_t pre, Int_t &nst, Float_t &qp){
738      for (Int_t e = pre*16; e < (pre+1)*16 ; e++){
739        dexy[l][m][e] = dexyc[l][m][e];
740      };  
741      this->FindBaseRaw(l,m,pre,nst,qp);
742    }
743    
744  void CaloLevel0::FindBaseRaw(Int_t l, Int_t m, Int_t pre){  void CaloLevel0::FindBaseRaw(Int_t l, Int_t m, Int_t pre){
745      Float_t minstrip = 100000.;    Int_t n = 0;
746      Float_t rms = 0.;    Float_t q = 0;
747      base[l][m][pre] = 0.;    this->FindBaseRaw(l,m,pre,n,q);
748    }
749    
750    void CaloLevel0::FindBaseRaw(Int_t l, Int_t m, Int_t pre, Int_t &nst, Float_t &qp){
751      //
752      Float_t minstrip = 100000.;
753      Float_t rms = 0.;
754      base[l][m][pre] = 0.;
755      qp = 0.;
756      //
757      for (Int_t e = pre*16; e < (pre+1)*16 ; e++){
758        if ( calgood[l][m][e] == 0. && obadmask[l][m][e] == 0 &&  dexy[l][m][e]-calped[l][m][e] < minstrip &&  dexy[l][m][e] > 0.) {
759          minstrip = dexy[l][m][e]-calped[l][m][e];
760          rms = calthr[l][m][pre];
761        };
762        qp += (dexy[l][m][e]-calped[l][m][e]-sbase[l][m][e]);
763      };
764      //
765      if ( debug && l==1 ){
766        printf("\n BASELINE CALCULATION for view %i pl %i pre %i: \n => minstrip %f rms %f \n => qp = %f \n",l,m,pre,minstrip,rms,qp);
767      };
768      if ( minstrip != 100000. ) {
769        Float_t strip6s = 0.;
770      for (Int_t e = pre*16; e < (pre+1)*16 ; e++){      for (Int_t e = pre*16; e < (pre+1)*16 ; e++){
771          if ( calgood[l][m][e] == 0. && obadmask[l][m][e] == 0 &&  dexy[l][m][e]-calped[l][m][e] < minstrip &&  dexy[l][m][e] > 0.) {        if ( (dexy[l][m][e]-calped[l][m][e]) >= minstrip && (dexy[l][m][e]-calped[l][m][e]) <= (minstrip+rms) ) {
772              minstrip = dexy[l][m][e]-calped[l][m][e];          strip6s += 1.;
773              rms = calthr[l][m][pre];          base[l][m][pre] += (dexy[l][m][e] - calped[l][m][e]);
774          };        };
775          //
776          //  compression
777          //
778          if ( abs((int)(dexy[l][m][e]-calped[l][m][e])) <= (minstrip+rms) ) {
779            dexyc[l][m][e] = 0.;
780          } else {
781            dexyc[l][m][e] = dexy[l][m][e];
782          };
783      };      };
784      if ( minstrip != 100000. ) {      //
785          Float_t strip6s = 0.;      nst = (Int_t)strip6s;
786          for (Int_t e = pre*16; e < (pre+1)*16 ; e++){      //
787              if ( (dexy[l][m][e]-calped[l][m][e]) >= minstrip && (dexy[l][m][e]-calped[l][m][e]) <= (minstrip+rms) ) {      if ( debug && l==1 ){
788                  strip6s += 1.;        printf(" strip6s %f \n",strip6s);
789                  base[l][m][pre] += (dexy[l][m][e] - calped[l][m][e]);      };
790              };      //  if ( strip6s >= 9. ){      
791              //      if ( strip6s >= 2. ){          
792              //  compression        Double_t arro = base[l][m][pre]/strip6s;
793              //        Float_t deci = 1000.*((float)arro - float(int(arro)));                  
794              if ( abs((int)(dexy[l][m][e]-calped[l][m][e])) <= (minstrip+rms) ) {        if ( deci < 500. ) {
795                  dexyc[l][m][e] = 0.;          arro = double(int(arro));
796              } else {        } else {
797                  dexyc[l][m][e] = dexy[l][m][e];          arro = 1. + double(int(arro));
798              };        };
799          };        base[l][m][pre] = arro;
800          if ( strip6s >= 9. ){              //
801              Double_t arro = base[l][m][pre]/strip6s;        // if too few strips were used to determine the baseline check if it is comparable with the previous event, if not mark it as bad
802              Float_t deci = 1000.*((float)arro - float(int(arro)));                            //
803              if ( deci < 500. ) {        if ( debug && l==1 ) printf(" Calculated baseline: base %f sbase*1.02 %f \n",base[l][m][pre],1.02*sbase[l][m][pre]);
804                  arro = double(int(arro));        //
805              } else {        if ( strip6s < 4 && base[l][m][pre] > 1.02*sbase[l][m][pre] && sbase[l][m][pre] > 0. ){
806                  arro = 1. + double(int(arro));          if ( debug ) printf(" Suspicious calculated baseline: base %f sbase*1.02 %f strip6s %i \n",base[l][m][pre],1.02*sbase[l][m][pre],(Int_t)strip6s);
807              };          base[l][m][pre] = 31000.;      
808              base[l][m][pre] = arro;        };
         } else {  
             base[l][m][pre] = 31000.;  
             for (Int_t e = pre*16; e < (pre+1)*16 ; e++){  
                 dexyc[l][m][e] = dexy[l][m][e];  
             };  
         };  
809      } else {      } else {
810          base[l][m][pre] = 31000.;        base[l][m][pre] = 31000.;
811          for (Int_t e = pre*16; e < (pre+1)*16 ; e++){
812            dexyc[l][m][e] = dexy[l][m][e];
813          };
814      };      };
815      } else {
816        base[l][m][pre] = 31000.;
817      };
818  }  }
819    
820  Int_t CaloLevel0::Calibrate(Int_t ei){  Int_t CaloLevel0::Calibrate(Int_t ei){
# Line 829  Int_t CaloLevel0::Calibrate(Int_t ei){ Line 872  Int_t CaloLevel0::Calibrate(Int_t ei){
872    Float_t ener;    Float_t ener;
873    Int_t doneb = 0;    Int_t doneb = 0;
874    Int_t donec = 0;    Int_t donec = 0;
875    Int_t ck = 0;    Int_t ck[6] = {0,0,0,0,0,0};
876    Int_t ipre = 0;    Int_t ipre = 0;
877    Int_t ip[3] = {0};    //  Int_t ip[3] = {0};
878      Int_t ip[3] = {0,0,0};
879    Float_t base0, base1, base2;    Float_t base0, base1, base2;
880    base0 = 0.;    base0 = 0.;
881    base1 = 0.;    base1 = 0.;
# Line 907  Int_t CaloLevel0::Calibrate(Int_t ei){ Line 951  Int_t CaloLevel0::Calibrate(Int_t ei){
951          //          //
952          pre = -1;          pre = -1;
953          cbase0 = 0.;          cbase0 = 0.;
954            Int_t nstt[2];
955            Float_t rqp[2];
956          for (Int_t i = 0; i < 3; i++){          for (Int_t i = 0; i < 3; i++){
957              nstt[0] = 0;
958              nstt[1] = 0;
959              rqp[0] = 0.;
960              rqp[1] = 0.;
961            for (Int_t j = 0; j < 2; j++){            for (Int_t j = 0; j < 2; j++){
962              pre = j + i*2;              pre = j + i*2;
963              //              //
964              // baseline check and calculation              // baseline check and calculation
965              //              //
966              if ( !isRAW ) {              if ( !isRAW ){
967                base[l][m][pre] = de->base[l][m][pre] ;                  //
968                  // if it is a compress event with fully transmitted pre try to calculate the baseline
969                  //
970                  if ( de->base[l][m][pre] != 0. && de->base[l][m][pre]<31000. ) {
971                    base[l][m][pre] = de->base[l][m][pre] ;  
972                  } else {
973                    FindBaseCompress(l,m,pre,nstt[j],rqp[j]);
974                  };
975                cbase0 += base[l][m][pre];                cbase0 += base[l][m][pre];
976              } else {              } else {
977                //                //
978                // if it is a raw event and we haven't checked                // if it is a raw event calculate the baseline.
               // yet, calculate the baseline.  
979                //                //
980                FindBaseRaw(l,m,pre);                FindBaseRaw(l,m,pre,nstt[j],rqp[j]);
981                cbase0 += base[l][m][pre];                cbase0 += base[l][m][pre];
982              };              };
983            };            };
984              //
985              // if we are able to calculate the baseline with more than 3 strips on one pre and not in the other one choose the pre with more calculated strips
986              //
987              if ( nstt[0] < 4 && nstt[1] >= 4 ) base[l][m][pre-1] = 31000.;
988              if ( nstt[0] >= 4 && nstt[1] < 4 ) base[l][m][pre] = 31000.;
989    //        //
990    //        // if we are NOT able to calculate the baseline with more than 3 strips on both pres take the baseline (if any) of the one which has less energy
991    //        //
992    //        if ( nstt[0] < 4 && nstt[1] < 4 ){
993    //          if ( rqp[0] >= rqp[1] ) base[l][m][pre-1] = 31000.;
994    //          if ( rqp[0] < rqp[1] ) base[l][m][pre] = 31000.;
995    //        };
996          };          };
997          //          //
998          // run over strips          // run over strips
# Line 935  Int_t CaloLevel0::Calibrate(Int_t ei){ Line 1003  Int_t CaloLevel0::Calibrate(Int_t ei){
1003            ip[i] = 0;            ip[i] = 0;
1004            for (Int_t n = i*32 ; n < (i+1)*32 ; n++){                            for (Int_t n = i*32 ; n < (i+1)*32 ; n++){                
1005              if (n%16 == 0) {              if (n%16 == 0) {
               ck = 0;  
1006                done = 0;                done = 0;
1007                doneb = 0;                doneb = 0;
1008                donec = 0;                donec = 0;
1009                pre++;                pre++;
1010                  ck[pre] = 0;
1011                qpre[pre] = 0.;                qpre[pre] = 0.;
1012              };              };
1013              //              //
# Line 949  Int_t CaloLevel0::Calibrate(Int_t ei){ Line 1017  Int_t CaloLevel0::Calibrate(Int_t ei){
1017              //              //
1018              if ( !done ){              if ( !done ){
1019                if ( (base[l][m][pre] == 31000. || base[l][m][pre] == 0.) ){                if ( (base[l][m][pre] == 31000. || base[l][m][pre] == 0.) ){
1020                  ck = 1;                  ck[pre] = 1;
1021                  if (pre%2 == 0) {                  if (pre%2 == 0) {
1022                    ip[i] = pre + 1;                    ip[i] = pre + 1;
1023                  } else {                  } else {
# Line 957  Int_t CaloLevel0::Calibrate(Int_t ei){ Line 1025  Int_t CaloLevel0::Calibrate(Int_t ei){
1025                  };                  };
1026                  if ( (base[l][m][ip[i]] == 31000. || base[l][m][ip[i]] == 0. || !crosst ) ){                  if ( (base[l][m][ip[i]] == 31000. || base[l][m][ip[i]] == 0. || !crosst ) ){
1027                    //                    //
1028                    ck = 2;                    ck[pre] = 2;
1029                    if ( sbase[l][m][pre] == 31000. || sbase[l][m][pre] == 0. ) {                    if ( sbase[l][m][pre] == 31000. || sbase[l][m][pre] == 0. ) {
1030                      ck = 3;                      ck[pre] = 3;
1031                    };                    };
1032                  };                  };
1033                  done = 1;                  done = 1;
# Line 969  Int_t CaloLevel0::Calibrate(Int_t ei){ Line 1037  Int_t CaloLevel0::Calibrate(Int_t ei){
1037              // CALIBRATION ALGORITHM              // CALIBRATION ALGORITHM
1038              //              //
1039              if ( !doneb ){              if ( !doneb ){
1040                if ( debug ) printf(" ck is %i \n",ck);                if ( debug ) printf(" ck[pre] is %i \n",ck[pre]);
1041                switch (ck) {                switch (ck[pre]) {
1042                case 0:                case 0:
1043                  base0 = base[l][m][pre];                  base0 = base[l][m][pre];
1044                  base2 = calbase[l][m][pre];                  base2 = calbase[l][m][pre];
1045                  if ( debug ) printf(" base0 = base l m pre = %f base2 = calbase l m pre = %f \n",base[l][m][pre],calbase[l][m][pre]);                  if ( debug ) printf(" base0 = base l%i m%i pre%i = %f base2 = calbase l m pre = %f \n",l,m,pre,base[l][m][pre],calbase[l][m][pre]);
1046                  break;                  break;
1047                case 1:                case 1:
1048                  base0 = base[l][m][ip[i]];                  base0 = base[l][m][ip[i]];
1049                  base2 = calbase[l][m][ip[i]];                  base2 = calbase[l][m][ip[i]];
1050                  if ( debug ) printf(" base0 = base l m ip(i) = %f base2 = calbase l m ip(i) = %f \n",base[l][m][ip[i]],calbase[l][m][ip[i]]);                  if ( debug ) printf(" base0 = base l%i m%i ip(i)%i = %f base2 = calbase l m ip(i) = %f \n",l,m,ip[i],base[l][m][ip[i]],calbase[l][m][ip[i]]);
1051                  break;                  break;
1052                case 2:                case 2:
1053                  base0 = sbase[l][m][pre];                  base0 = sbase[l][m][pre];
1054                  base2 = calbase[l][m][pre];                      base2 = calbase[l][m][pre];    
1055                  if ( debug ) printf(" base0 = sbase l m pre = %f base2 = calbase l m pre = %f \n",sbase[l][m][pre],calbase[l][m][pre]);                  if ( debug ) printf(" base0 = sbase l%i m%i pre%i = %f base2 = calbase l m pre = %f \n",l,m,pre,sbase[l][m][pre],calbase[l][m][pre]);
1056                  break;                  break;
1057                case 3:                case 3:
1058                  base0 = calbase[l][m][pre];                  base0 = calbase[l][m][pre];
1059                  base2 = calbase[l][m][pre];                  base2 = calbase[l][m][pre];
1060                  if ( debug ) printf(" base0 = calbase l m pre = %f base2 = calbase l m pre = %f \n",calbase[l][m][pre],calbase[l][m][pre]);                  if ( debug ) printf(" base0 = calbase l%i m%i pre%i = %f base2 = calbase l m pre = %f \n",l,m,pre,calbase[l][m][pre],calbase[l][m][pre]);
1061                  break;                  break;
1062                };                };
1063                base1 = calbase[l][m][pre];                base1 = calbase[l][m][pre];
# Line 1017  Int_t CaloLevel0::Calibrate(Int_t ei){ Line 1085  Int_t CaloLevel0::Calibrate(Int_t ei){
1085              };              };
1086            };            };
1087            if ( crosst ){            if ( crosst ){
1088              if (ck == 1){              if (ck[pre] == 1 || ck[pre-1] == 1){
1089                if (ip[i]%2 == 0) {                if (ck[pre] == 1){
1090                  ipre = ip[i] + 1;                  ipre = pre;
1091                    ip[i] = pre - 1;
1092                } else {                } else {
1093                  ipre = ip[i] - 1;                  ipre = pre - 1;
1094                    ip[i] = pre;
1095                };                };
1096                  //              if (ip[i]%2 == 0) {
1097                  //                ipre = ip[i] + 1;
1098                  //              } else {
1099                  //                ipre = ip[i] - 1;
1100                  //              };
1101                for (Int_t j = ipre*16 ; j < (ipre+1)*16 ; j++){                for (Int_t j = ipre*16 ; j < (ipre+1)*16 ; j++){
1102                  if ( !ctground ){                  if ( !ctground ){
1103                    clevel1->estrip[j][m][l] += (qpre[ipre] - qpre[ip[i]]) * ctprecor[l][m][ipre];                    clevel1->estrip[j][m][l] += qpre[ipre] * ctprecor[l][m][ipre] - qpre[ip[i]] * ctprecor[l][m][ip[i]];
1104                  } else {                  } else {
1105                    clevel1->estrip[j][m][l] += (qpre[ipre] - qpre[ip[i]]) * 0.00478;                    clevel1->estrip[j][m][l] += (qpre[ipre] - qpre[ip[i]]) * 0.00478;
1106                  };                  };
1107                };                };
1108              };              };
1109              if (ck == 2){              if (ck[pre] == 2 && ck[pre-1] == 2){
1110                for (Int_t j = i*32 ; j < (i+1)*32 ; j++){                for (Int_t j = i*32 ; j < (i+1)*32 ; j++){
1111                  ipre = j/16 + 1;                  //              ipre = j/16 + 1;
1112                    ipre = j/16 ;
1113                  if ( !ctground ){                  if ( !ctground ){
1114                    clevel1->estrip[j][m][l] +=  qpre[ipre] * ctprecor[l][m][ipre];                    clevel1->estrip[j][m][l] +=  qpre[ipre] * ctprecor[l][m][ipre];
1115                  } else {                  } else {
Legend:
Removed from v.1.15  
changed lines
  Added in v.1.20
  ViewVC Help
Powered by ViewVC 1.1.23