/[PAMELA software]/DarthVader/TrackerLevel2/src/TrkLevel1.cpp
ViewVC logotype

Diff of /DarthVader/TrackerLevel2/src/TrkLevel1.cpp

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

revision 1.11 by pam-fi, Thu Jan 11 10:20:58 2007 UTC revision 1.20 by pam-fi, Mon Aug 20 15:47:57 2007 UTC
# Line 239  Bool_t TrkCluster::IsBad(Int_t nbad){ Line 239  Bool_t TrkCluster::IsBad(Int_t nbad){
239      il = indmax;      il = indmax;
240      ir = indmax;      ir = indmax;
241      for(Int_t i=1; i<nbad; i++){      for(Int_t i=1; i<nbad; i++){
242          if (ir == CLlength && il == 0)break;          if (ir == CLlength-1 && il == 0)break;
243          else if (ir == CLlength && il != 0)il--;          else if (ir == CLlength-1 && il != 0)il--;
244          else if (ir != CLlength && il == 0)ir++;          else if (ir != CLlength-1 && il == 0)ir++;
245          else{          else{
246              if(clsignal[il-1] > clsignal[ir+1])il--;              if(clsignal[il-1] > clsignal[ir+1])il--;
247              else ir++;              else ir++;
# Line 264  Bool_t TrkCluster::IsSaturated(Int_t nba Line 264  Bool_t TrkCluster::IsSaturated(Int_t nba
264      il = indmax;      il = indmax;
265      ir = indmax;      ir = indmax;
266      for(Int_t i=1; i<nbad; i++){      for(Int_t i=1; i<nbad; i++){
267          if (ir == CLlength && il == 0)break;          if (ir == CLlength-1 && il == 0)break;
268          else if (ir == CLlength && il != 0)il--;          else if (ir == CLlength-1 && il != 0)il--;
269          else if (ir != CLlength && il == 0)ir++;          else if (ir != CLlength-1 && il == 0)ir++;
270          else{          else{
271              if(clsignal[il-1] > clsignal[ir+1])il--;              if(clsignal[il-1] > clsignal[ir+1])il--;
272              else ir++;              else ir++;
# Line 291  void TrkCluster::Dump(){ Line 291  void TrkCluster::Dump(){
291      cout << "Position of maximun "<< maxs <<endl;      cout << "Position of maximun "<< maxs <<endl;
292      cout << "Multiplicity        "<< GetMultiplicity() <<endl;      cout << "Multiplicity        "<< GetMultiplicity() <<endl;
293      cout << "Tot signal          "<< GetSignal() << " (ADC channels)"<<endl ;      cout << "Tot signal          "<< GetSignal() << " (ADC channels)"<<endl ;
294      cout << "Signal/Noise        "<< GetSignalToNoise();      cout << "Signal/Noise        "<< GetSignalToNoise()<<endl;
295      cout <<endl<< "Strip signals       ";      cout << "COG                 "<< GetCOG(0)<<endl;;
296        cout << "Strip signals       ";
297      for(Int_t i =0; i<CLlength; i++)cout << " " <<clsignal[i];      for(Int_t i =0; i<CLlength; i++)cout << " " <<clsignal[i];
298      cout <<endl<< "Strip sigmas        ";      cout <<endl<< "Strip sigmas        ";
299      for(Int_t i =0; i<CLlength; i++)cout << " " <<clsigma[i];      for(Int_t i =0; i<CLlength; i++)cout << " " <<clsigma[i];
# Line 313  void TrkCluster::Dump(){ Line 314  void TrkCluster::Dump(){
314  /**  /**
315   * Method to fill a level1 struct with only one cluster (done to use F77 p.f.a. routines on a cluster basis).   * Method to fill a level1 struct with only one cluster (done to use F77 p.f.a. routines on a cluster basis).
316   */   */
317  cTrkLevel1* TrkCluster::GetLevel1Struct(){  void TrkCluster::GetLevel1Struct(cTrkLevel1* l1){
318                                    
319  //    cTrkLevel1* l1 = new cTrkLevel1;  //    cTrkLevel1* l1 = new cTrkLevel1;
320    
321      cTrkLevel1* l1 = &level1event_ ;  //    cTrkLevel1* l1 = &level1event_ ;
322                    
323      l1->nclstr1 = 1;      l1->nclstr1 = 1;
324      l1->view[0] = view;      l1->view[0] = view;
# Line 335  cTrkLevel1* TrkCluster::GetLevel1Struct( Line 336  cTrkLevel1* TrkCluster::GetLevel1Struct(
336          l1->clbad[i] = clbad[i];          l1->clbad[i] = clbad[i];
337      };      };
338            
339      return l1;  //    return l1;
340  };  };
341  //--------------------------------------  //--------------------------------------
342  //  //
# Line 397  Float_t TrkCluster::GetETA(Int_t neta, f Line 398  Float_t TrkCluster::GetETA(Int_t neta, f
398  //    cout << "GetETA(neta,angle) "<< neta << " "<< angle;  //    cout << "GetETA(neta,angle) "<< neta << " "<< angle;
399  //      LoadPfaParam();  //      LoadPfaParam();
400    
401        TrkParams::Load(4);
402        if( !TrkParams::IsLoaded(4) ){
403            cout << "int Trajectory::DoTrack2(float* al) --- ERROR --- p.f.a. parameters  not loaded"<<endl;
404            return 0;
405        }
406    
407      float ax = angle;      float ax = angle;
408      int ic = 1;      int ic = 1;
409      GetLevel1Struct();      GetLevel1Struct();
# Line 425  TrkLevel1::TrkLevel1(){ Line 432  TrkLevel1::TrkLevel1(){
432              cnn[j][i]=0;              cnn[j][i]=0;
433          };          };
434      };      };
435        TrkParams::SetTrackingMode();
436        TrkParams::SetPrecisionFactor();
437        TrkParams::SetStepMin();
438        TrkParams::SetPFA();
439    }
440    //--------------------------------------
441    //
442    //
443    //--------------------------------------
444    void TrkLevel1::Set(){
445        if(!Cluster)Cluster = new TClonesArray("TrkCluster");
446  }  }
447  //--------------------------------------  //--------------------------------------
448  //  //
# Line 450  void TrkLevel1::Dump(){ Line 468  void TrkLevel1::Dump(){
468      for(int i=0; i<this->nclstr(); i++)     ((TrkCluster *)t[i])->Dump();      for(int i=0; i<this->nclstr(); i++)     ((TrkCluster *)t[i])->Dump();
469            
470  }  }
471    /**
472     * \brief Dump processing status
473     */
474    void TrkLevel1::StatusDump(int view){
475        cout << "DSP n. "<<view+1<<" (level1-)status: "<<hex<<showbase<<good[view]<<dec<<endl;    
476    };
477    /**
478     * \brief Check event status
479     *
480     * Check the event status, according to a flag-mask given as input.
481     * Return true if the view passes the check.
482     *
483     * @param view View number (0-11)
484     * @param flagmask Mask of flags to check (eg. flagmask=0x111 no missing packet,
485     *  no crc error, no software alarm)
486     *
487     * @see TrkLevel2 class definition to know how the status flag is defined
488     *
489     */
490    Bool_t TrkLevel1::StatusCheck(int view, int flagmask){
491    
492        if( view<0 || view >= 12)return false;
493        return !(good[view]&flagmask);
494    
495    };
496    
497    
498  //--------------------------------------  //--------------------------------------
499  //  //
500  //  //
# Line 459  void TrkLevel1::Dump(){ Line 504  void TrkLevel1::Dump(){
504   */   */
505  void TrkLevel1::SetFromLevel1Struct(cTrkLevel1 *l1, Bool_t full){  void TrkLevel1::SetFromLevel1Struct(cTrkLevel1 *l1, Bool_t full){
506    
507    //    cout << "void TrkLevel1::SetFromLevel1Struct(cTrkLevel1 *l1, Bool_t full)"<<endl;
508        
509        Clear();
510      //  ---------------      //  ---------------
511      //  *** CLUSTER ***      //  *** CLUSTER ***
512      //  ---------------      //  ---------------
# Line 484  void TrkLevel1::SetFromLevel1Struct(cTrk Line 532  void TrkLevel1::SetFromLevel1Struct(cTrk
532              t_cl->clsigma  = new Float_t[t_cl->CLlength];              t_cl->clsigma  = new Float_t[t_cl->CLlength];
533              t_cl->cladc    = new Int_t[t_cl->CLlength];              t_cl->cladc    = new Int_t[t_cl->CLlength];
534              t_cl->clbad    = new Bool_t[t_cl->CLlength];              t_cl->clbad    = new Bool_t[t_cl->CLlength];
535    
536              Int_t index = 0;              Int_t index = 0;
537              for(Int_t is = from; is < to; is++ ){              for(Int_t is = from; is < to; is++ ){
538                  t_cl->clsignal[index] = (Float_t) l1->clsignal[is];                  t_cl->clsignal[index] = (Float_t) l1->clsignal[is];
# Line 515  void TrkLevel1::SetFromLevel1Struct(cTrk Line 564  void TrkLevel1::SetFromLevel1Struct(cTrk
564   * Fills a struct cTrkLevel1 with values from a TrkLevel1 object (to put data into a F77 common).   * Fills a struct cTrkLevel1 with values from a TrkLevel1 object (to put data into a F77 common).
565   */   */
566    
567  cTrkLevel1* TrkLevel1::GetLevel1Struct() {  void TrkLevel1::GetLevel1Struct(cTrkLevel1* l1) {
568            
569      cTrkLevel1 *l1=0;  //    cTrkLevel1* l1 = &level1event_ ;
570      //      
571      for(Int_t i=0; i<12 ; i++){      for(Int_t i=0; i<12 ; i++){
572          l1->good[i] = good[i];          l1->good[i] = good[i];
573          for(Int_t j=0; j<24 ; j++){          for(Int_t j=0; j<24 ; j++){
574              l1->cnev[j][i]    = cn[j][i];              l1->cnev[j][i]    = cn[j][i]  ;
575  //          l1->cnrmsev[j][i] = cnrms[j][i];              l1->cnnev[j][i]   = cnn[j][i] ;
576              l1->cnnev[j][i]   = cnn[j][i];              l1->cnrmsev[j][i] = 0. ;
577          };          };
578            l1->fshower[i] = 0;
579      };      };
580        
581  //  *** CLUSTERS ***      l1->nclstr1=0;
582        l1->totCLlength=0;
583        Int_t index=0;
584      if(Cluster){      if(Cluster){
585          l1->nclstr1 =  Cluster->GetEntries();          Int_t i=0;
586          for(Int_t i=0;i<l1->nclstr1;i++){          for(Int_t ii=0;ii<Cluster->GetEntries();ii++){
587                        TrkCluster *clu = GetCluster(ii);
588              l1->view[i]     = ((TrkCluster *)Cluster->At(i))->view;              // ----------------------------------------
589              l1->maxs[i]     = ((TrkCluster *)Cluster->At(i))->maxs;              // attenzione!!
590              // COMPLETARE //              // se il cluster non e` salvato (view = 0)
591              // COMPLETARE //              // DEVE essere escluso dal common F77
592              // COMPLETARE //              // ----------------------------------------
593              // COMPLETARE //              if(clu->view != 0 ){
594              // COMPLETARE //                  l1->view[i]     = clu->view;
595              // COMPLETARE //                  l1->ladder[i]   = clu->GetLadder();
596                            l1->maxs[i]     = clu->maxs;
597                    l1->mult[i]     = clu->GetMultiplicity();
598                    l1->dedx[i]     = clu->GetSignal();
599                    l1->indstart[i] = index+1;
600                    l1->indmax[i]   = l1->indstart[i] + clu->indmax;
601                    l1->totCLlength += clu->CLlength;
602                    for(Int_t iw=0; iw < clu->CLlength; iw++){
603                        l1->clsignal[index] = clu->clsignal[iw];
604                        l1->clsigma[index]  = clu->clsigma[iw];
605                        l1->cladc[index]    = clu->cladc[iw];
606                        l1->clbad[index]    = clu->clbad[iw];
607                        index++;
608                    }
609                    i++;
610                }
611          }          }
612          // COMPLETARE //          l1->nclstr1 =  i;      
         // COMPLETARE //  
         // COMPLETARE //  
         // COMPLETARE //  
         // COMPLETARE //  
         // COMPLETARE //  
613      }      }
614      return l1;  
615    //    return l1;
616  }  }
617  //--------------------------------------  //--------------------------------------
618  //  //
# Line 601  TrkCluster *TrkLevel1::GetCluster(int is Line 663  TrkCluster *TrkLevel1::GetCluster(int is
663  //  //
664  //  //
665  //--------------------------------------  //--------------------------------------
666    // /**
667    //  * Load Position-Finding-Algorythm parameters (call the F77 routine).
668    //  *
669    //  */
670    // int TrkLevel1::LoadPfaParam(TString path){
671            
672    //     if( path.IsNull() ){
673    //      path = gSystem->Getenv("PAM_CALIB");
674    //      if(path.IsNull()){
675    //          cout << " TrkLevel1::LoadPfaParam() ==> No PAMELA environment variables defined "<<endl;
676    //          return 0;
677    //      }
678    //      path.Append("/trk-param/eta_param-0/");
679    //     }
680    
681    //     strcpy(path_.path,path.Data());
682    //     path_.pathlen = path.Length();
683    //     path_.error   = 0;
684    //     cout <<"Loading p.f.a. parameters: "<<path<<endl;
685    //     return readetaparam_();
686    // }
687    
688    // /**
689    //  * Load magnetic field parameters (call the F77 routine).
690    //  *
691    //  */
692    // int TrkLevel1::LoadFieldParam(TString path){
693            
694    // //    if( strcmp(path_.path,path.Data()) ){
695    //     if( path.IsNull() ){
696    //      path = gSystem->Getenv("PAM_CALIB");
697    //      if(path.IsNull()){
698    //          cout << " TrkLevel1::LoadFieldParam() ==> No PAMELA environment variables defined "<<endl;
699    //          return 0;
700    //      }
701    //      path.Append("/trk-param/field_param-0/");
702    //     }
703    //     cout <<"Loading magnetic field "<<path<<endl;
704    //     strcpy(path_.path,path.Data());
705    //     path_.pathlen = path.Length();
706    //     path_.error   = 0;
707    //     return readb_();
708    // //    }      
709    // //    return 0;
710    // }
711    // /**
712    //  * Load magnetic field parameters (call the F77 routine).
713    //  *
714    //  */
715    // int TrkLevel1::LoadChargeParam(TString path){
716            
717    // //    if( strcmp(path_.path,path.Data()) ){
718    //     if( path.IsNull() ){
719    //      path = gSystem->Getenv("PAM_CALIB");
720    //      if(path.IsNull()){
721    //          cout << " TrkLevel1::LoadChargeParam() ==> No PAMELA environment variables defined "<<endl;
722    //          return 0;
723    //      }
724    //      path.Append("/trk-param/charge_param-1/");
725    //     }
726    //     cout <<"Loading charge-correlation parameters: "<<path<<endl;
727    //     strcpy(path_.path,path.Data());
728    //     path_.pathlen = path.Length();
729    //     path_.error   = 0;
730    //     return readchargeparam_();
731    // //    }      
732    // //    return 0;
733    // }
734    // /**
735    //  * Load magnetic field parameters (call the F77 routine).
736    //  *
737    //  */
738    // int TrkLevel1::LoadAlignmentParam(TString path){
739            
740    // //    if( strcmp(path_.path,path.Data()) ){
741    //     if( path.IsNull() ){
742    //      path = gSystem->Getenv("PAM_CALIB");
743    //      if(path.IsNull()){
744    //          cout << " TrkLevel1::LoadAlignmentParam() ==> No PAMELA environment variables defined "<<endl;
745    //          return 0;
746    //      }
747    //      path.Append("/trk-param/align_param-0/");
748    //     }
749    //     cout <<"Loading alignment parameters: "<<path<<endl;
750    //     strcpy(path_.path,path.Data());
751    //     path_.pathlen = path.Length();
752    //     path_.error   = 0;
753    //     return readalignparam_();
754    // //    }      
755    // //    return 0;
756    // }
757    // /**
758    //  * Load magnetic field parameters (call the F77 routine).
759    //  *
760    //  */
761    // int TrkLevel1::LoadMipParam(TString path){
762            
763    // //    if( strcmp(path_.path,path.Data()) ){
764    //     if( path.IsNull() ){
765    //      path = gSystem->Getenv("PAM_CALIB");
766    //      if(path.IsNull()){
767    //          cout << " TrkLevel1::LoadMipParam() ==> No PAMELA environment variables defined "<<endl;
768    //          return 0;
769    //      }
770    //      path.Append("/trk-param/mip_param-0/");
771    //     }
772    //     cout <<"Loading ADC-to-MIP conversion parameters: "<<path<<endl;
773    //     strcpy(path_.path,path.Data());
774    //     path_.pathlen = path.Length();
775    //     path_.error   = 0;
776    //     return readmipparam_();
777    // //    }      
778    // //    return 0;
779    // }
780    // /**
781    //  * Load magnetic field parameters (call the F77 routine).
782    //  *
783    //  */
784    // int TrkLevel1::LoadVKMaskParam(TString path){
785            
786    // //    if( strcmp(path_.path,path.Data()) ){
787    //     if( path.IsNull() ){
788    //      path = gSystem->Getenv("PAM_CALIB");
789    //      if(path.IsNull()){
790    //          cout << " TrkLevel1::LoadVKMaskParam() ==> No PAMELA environment variables defined "<<endl;
791    //          return 0;
792    //      }
793    //      path.Append("/trk-param/mask_param-1/");
794    //     }
795    //     cout <<"Loading VK-mask parameters: "<<path<<endl;
796    //     strcpy(path_.path,path.Data());
797    //     path_.pathlen = path.Length();
798    //     path_.error   = 0;
799    //     return readvkmask_();
800    // //    }      
801    // //    return 0;
802    // }
803    
804    // /**
805    //  * Load all (default) parameters. Environment variable must be defined.
806    //  *
807    //  */
808    // int TrkLevel1::LoadParams(){
809    
810    //     int result=0;
811        
812    //     result = result * LoadFieldParam();
813    //     result = result * LoadPfaParam();
814    //     result = result * LoadChargeParam();
815    //     result = result * LoadAlignmentParam();
816    //     result = result * LoadMipParam();
817    //     result = result * LoadVKMaskParam();
818    
819    //     return result;
820    // }
821    
822    
823    
824    int TrkLevel1::GetPfaNbinsAngle(){
825        TrkParams::Load(4);
826        if( !TrkParams::IsLoaded(4) ){
827            cout << "int TrkLevel1::GetPfaNbinsAngle() --- ERROR --- p.f.a. parameters  not loaded"<<endl;
828            return 0;
829        }
830        return pfa_.nangbin;
831    };
832    
833    int TrkLevel1::GetPfaNbinsETA(){
834        TrkParams::Load(4);
835        if( !TrkParams::IsLoaded(4) ){
836            cout << "int TrkLevel1::GetPfaNbinsETA() --- ERROR --- p.f.a. parameters  not loaded"<<endl;
837            return 0;
838        }
839        return pfa_.netaval;
840    };
841    
842  /**  /**
843   * Load Position-Finding-Algorythm parameters (call the F77 routine).   *
844   *   *
845   */   */
846  int TrkLevel1::LoadPfaParam(TString path){  float* TrkLevel1::GetPfaCoord(TString pfa, int nview, int nladder, int nang){
847            
848      if( strcmp(path_.path,path.Data()) ){      TrkParams::Load(4);
849          cout <<"Loading p.f.a. param::eters\n";      if( !TrkParams::IsLoaded(4) ){
850          strcpy(path_.path,path.Data());          cout << "float* TrkLevel1::GetPfaCoord(TString pfa, int nview, int nladder, int nang) --- ERROR --- p.f.a. parameters  not loaded"<<endl;
851          path_.pathlen = path.Length();          return 0;
852          path_.error   = 0;      }
853          return readetaparam_();    
854      }        int nbins = GetPfaNbinsETA();
855      return 0;      if(!nbins)return 0;
856    
857        float *fcorr = new float [nbins];
858    
859        if(!pfa.CompareTo("ETA2",TString::kIgnoreCase)){
860            for(int ib=0; ib<nbins; ib++){
861                fcorr[ib] = pfa_.feta2[nang][nladder][nview][ib];
862                cout << pfa_.eta2[nang][ib] << " - " <<  pfa_.feta2[nang][nladder][nview][ib]<<endl;;
863            }
864        }else if (!pfa.CompareTo("ETA3",TString::kIgnoreCase)){
865            for(int ib=0; ib<nbins; ib++)fcorr[ib] = pfa_.feta3[nang][nladder][nview][ib];
866        }else if (!pfa.CompareTo("ETA4",TString::kIgnoreCase)){
867            for(int ib=0; ib<nbins; ib++)fcorr[ib] = pfa_.feta4[nang][nladder][nview][ib];
868        }else{
869            cout << pfa<<" pfa parameters not implemented "<<endl;
870            return 0;
871        }    
872    
873        return fcorr;
874    
875    };
876    
877    float* TrkLevel1::GetPfaAbs(TString pfa, int nang){
878      
879        TrkParams::Load(4);
880        if( !TrkParams::IsLoaded(4) ){
881            cout << "float* TrkLevel1::GetPfaAbs(TString pfa, int nang) --- ERROR --- p.f.a. parameters  not loaded"<<endl;
882            return 0;
883        }
884    
885        int nbins = GetPfaNbinsETA();
886        if(!nbins)return 0;
887    
888        float *fcorr = new float [nbins];
889    
890        if(!pfa.CompareTo("ETA2",TString::kIgnoreCase)){
891            for(int ib=0; ib<nbins; ib++)fcorr[ib] = pfa_.eta2[nang][ib];
892        }else if (!pfa.CompareTo("ETA3",TString::kIgnoreCase)){
893            for(int ib=0; ib<nbins; ib++)fcorr[ib] = pfa_.eta3[nang][ib];
894        }else if (!pfa.CompareTo("ETA4",TString::kIgnoreCase)){
895            for(int ib=0; ib<nbins; ib++)fcorr[ib] = pfa_.eta4[nang][ib];
896        }else{
897            cout << pfa<<" pfa parameters not implemented "<<endl;
898            return 0;
899        }    
900    
901        return fcorr;
902    
903    };
904    
905    /**
906     * Method to call the F77 routine that performs level1->level2 processing.
907     * The level2 output is stored in a common block, which can be retrieved
908     * by mean of the method TrkLevel2::SetFromLevel2Struct().
909     * NB If the TrkLevel1 object is readout from a tree, and the
910     * TrkLevel1::ProcessEvent(int pfa) is used to reprocess the event, attention
911     * should be payed to the fact that single clusters (clusters not associated
912     * with any track) might not be stored. Full reprocessing should be done starting
913     * from level0 data.
914     */
915    //int TrkLevel1::ProcessEvent(int pfa){
916    int TrkLevel1::ProcessEvent(){
917    
918    //    cout << "int TrkLevel1::ProcessEvent()" << endl;
919        TrkParams::Load( );
920        if( !TrkParams::IsLoaded() )return 0;
921    
922        GetLevel1Struct();
923    
924    //    analysisflight_(&pfa);
925    //    TrkParams::SetPFA(pfa);
926        analysisflight_();
927    
928        return 1;
929    
930  }  }
931    
932    
Legend:
Removed from v.1.11  
changed lines
  Added in v.1.20
  ViewVC Help
Powered by ViewVC 1.1.23