epic-lfhcal-tbana/NewStructure/HGCROC_Convert.cc

0001 #include "HGCROC_Convert.h"
0002
0003 #include "Analyses.h"
0004 #include "Event.h"
0005 #include "Tile.h"
0006 #include "HGCROC.h"
0007 #include "CommonHelperFunctions.h"
0008 #include "Setup.h"
0009
0010 #include "include/h2g_decode/run_decoder.h"
0011
0012
0013 int run_hgcroc_conversion(Analyses *analysis, waveform_fit_base *waveform_builder) {
0014     #ifdef DECODE_HGCROC // Built for HGCROC decoding
0015     std::cout << "Setting up event parameters for HGCROC data" << std::endl;
0016     // Make sure we actually have a waveform builder
0017     if (waveform_builder == nullptr) {
0018         std::cout << "No waveform builder specified" << std::endl;
0019         return 1;
0020     }
0021
0022     // Check mapping file
0023     if (analysis->MapInputName.IsNull()) {
0024         std::cout << "No mapping file specified" << std::endl;
0025         // return 1;
0026     }
0027     analysis->setup->Initialize(analysis->MapInputName.Data(), analysis->debug);
0028
0029     // Check run list file - Used to get run parameters automatically
0030     if (analysis->RunListInputName.IsNull()) {
0031         std::cout << "No run list file specified" << std::endl;
0032         // return 1;
0033     }
0034     std::map<int, RunInfo> ri = readRunInfosFromFile(analysis->RunListInputName.Data(), analysis->debug, 0);
0035
0036     // Set up the file to convert
0037     if (analysis->ASCIIinputName.IsNull()) {
0038         std::cout << "No input file specified" << std::endl;
0039         // return 1;
0040     }
0041
0042     // Set up the static event parameters
0043     // Clean up file headers'
0044     // THIS WILL HAVE TO CHANGE
0045     TObjArray* tok = analysis->ASCIIinputName.Tokenize("/");
0046     // get file name
0047     TString file = ((TObjString*)tok->At(tok->GetEntries() - 1))->String();
0048     delete tok;
0049     tok=file.Tokenize("_");
0050     TString header = ((TObjString*)tok->At(0))->String();
0051     delete tok;
0052
0053     // Get run number from file & find necessary run infos
0054     TString RunString=header(3,header.Sizeof());
0055     std::map<int,RunInfo>::iterator it=ri.find(RunString.Atoi());
0056     std::cout<<RunString.Atoi()<<"\t"<<it->second.species<<std::endl;
0057     analysis->event.SetRunNumber(RunString.Atoi());
0058     analysis->event.SetBeamName(it->second.species);
0059     analysis->event.SetBeamID(it->second.pdg);
0060     analysis->event.SetBeamEnergy(it->second.energy);
0061     analysis->event.SetVop(it->second.vop);
0062     analysis->event.SetVov(it->second.vop-it->second.vbr);
0063     analysis->event.SetBeamPosX(it->second.posX);
0064     analysis->event.SetBeamPosY(it->second.posY);
0065     analysis->calib.SetRunNumber(RunString.Atoi());
0066     analysis->calib.SetVop(it->second.vop);
0067     analysis->calib.SetVov(it->second.vop-it->second.vbr);
0068     analysis->event.SetROtype(ReadOut::Type::Hgcroc);
0069
0070     // Run the event builder
0071     std::list<aligned_event*> *events = run_event_builder((char*)analysis->ASCIIinputName.Data());
0072
0073     std::cout << "completed HGCROC event builder!\n" << std::endl;
0074
0075
0076     // convert from the aligned_events datatype to the Event datatype
0077     int event_number = 0;
0078     for (auto it = events->begin(); it != events->end(); it++) {
0079         if (true || event_number % 100 == 0) {
0080             std::cout << "\rFitting event " << event_number << std::flush;
0081         }
0082         aligned_event *ae = *it;
0083         analysis->event.SetEventID(event_number);
0084         event_number++;
0085         // Loop over each tile
0086         for (int i = 0; i < ae->get_num_fpga(); i++) {
0087             auto single_kcu = ae->get_event(i);
0088             for (int j = 0; j < ae->get_channels_per_fpga(); j++) {
0089                 int channel_number = i * ae->get_channels_per_fpga() + j;
0090                 int x, y, z;
0091                 if (decode_position(channel_number, x, y, z)) {
0092                     Hgcroc *tile = new Hgcroc();
0093                     auto cell_id = analysis->setup->GetCellID(y, x, z, 0);  // needs to be adapted once we have multiple modules
0094                     tile->SetCellID(cell_id);        // TODO: This is not the same cell ID as Fredi and Vincent set up
0095                     tile->SetROtype(ReadOut::Type::Hgcroc);
0096                     tile->SetLocalTriggerBit(0);            // What are these supposed to be?
0097                     tile->SetLocalTriggerPrimitive(0);
0098                     tile->SetE(0);                          // Need to process waveform to get this
0099                     tile->SetTOA(0);                        // Need to process waveform to get this
0100                     tile->SetTOT(0);                        // Need to process waveform to get this
0101
0102                     tile->SetNsample(single_kcu->get_n_samples());
0103                     for (int sample = 0; sample < single_kcu->get_n_samples(); sample++) {
0104                         tile->AppendWaveformADC(single_kcu->get_sample_adc(j, sample));
0105                         tile->AppendWaveformTOA(single_kcu->get_sample_toa(j, sample));
0106                         tile->AppendWaveformTOT(single_kcu->get_sample_tot(j, sample));
0107                     }
0108
0109                     // process tile waveform
0110                     waveform_builder->set_waveform(tile->GetADCWaveform());
0111                     waveform_builder->fit();
0112                     tile->SetE(waveform_builder->get_E());
0113                     tile->SetPedestal(waveform_builder->get_pedestal());
0114
0115                     analysis->event.AddTile(tile);
0116                 }
0117             }
0118         // Fill the event
0119         analysis->TdataOut->Fill();
0120         analysis->event.ClearTiles();
0121         }
0122     }
0123     std::cout << "\nFinished converting events\n" << std::endl;
0124     analysis->RootOutput->cd();
0125     // setup
0126     RootSetupWrapper rswtmp=RootSetupWrapper(analysis->setup);
0127     analysis->rsw=rswtmp;
0128     analysis->TsetupOut->Fill();
0129     // calib
0130     analysis->TcalibOut->Fill();
0131     analysis->TcalibOut->Write();
0132     // event data
0133     analysis->TdataOut->Fill();
0134     analysis->TsetupOut->Write();
0135     analysis->TdataOut->Write();
0136
0137     analysis->RootOutput->Close();
0138     return true;
0139
0140     #else
0141     std::cout << "This code is not built for HGCROC decoding" << std::endl;
0142     analysis->RootOutput->cd();
0143     // setup
0144     RootSetupWrapper rswtmp=RootSetupWrapper(analysis->setup);
0145     analysis->rsw=rswtmp;
0146     analysis->TsetupOut->Fill();
0147     // calib
0148     analysis->TcalibOut->Fill();
0149     analysis->TcalibOut->Write();
0150     // event data
0151     analysis->TdataOut->Fill();
0152     analysis->TsetupOut->Write();
0153     analysis->TdataOut->Write();
0154
0155     analysis->RootOutput->Close();
0156     return false;
0157
0158     #endif
0159 }
0160
0161 bool decode_position(int channel, int &x, int &y, int &z) {
0162     int channel_map[72] = {64, 63, 66, 65, 69, 70, 67, 68,  // this goes from 0 to 63 in lhfcal space to 0 to 71 in asic space
0163                            55, 56, 57, 58, 62, 61, 60, 59,  // So channel 0 on the detector is channel 64 on the asic
0164                            45, 46, 47, 48, 52, 51, 50, 49,  // What we want is the reverse of this, going from 0 to 71 in asic
0165                            37, 36, 39, 38, 42, 43, 40, 41,  // space to 0 to 63 in lhfcal space
0166                            34, 33, 32, 31, 27, 28, 29, 30,
0167                            25, 26, 23, 24, 20, 19, 22, 21,
0168                            16, 14, 15, 12,  9, 11, 10, 13,
0169                             7,  6,  5,  4,  0,  1,  2,  3,
0170                             -1, -1, -1, -1, -1, -1, -1, -1};
0171
0172     int fpga_factor[4] = {1, 3, 0, 2};
0173
0174     int asic_channel = channel % 72;
0175     // find the lhfcal channel
0176     int lhfcal_channel = -1;
0177     for (int i = 0; i < 72; i++) {
0178         if (channel_map[i] == asic_channel) {
0179             lhfcal_channel = i;
0180             break;
0181         }
0182     }
0183
0184     if (lhfcal_channel == -1) {
0185         // These are the empty channels
0186         x = -1;
0187         y = -1;
0188         z = -1;
0189         return false;
0190     }
0191
0192     int candy_bar_index = lhfcal_channel % 8;
0193     if (candy_bar_index < 4) {
0194         y = 1;
0195     } else {
0196         y = 0;
0197     }
0198
0199     if (candy_bar_index == 0 || candy_bar_index == 7) {
0200         x = 0;
0201     } else if (candy_bar_index == 1 || candy_bar_index == 6) {
0202         x = 1;
0203     } else if (candy_bar_index == 2 || candy_bar_index == 5) {
0204         x = 2;
0205     } else {
0206         x = 3;
0207     }
0208
0209     int fpga = channel / 144;
0210     int asic = (channel % 144) / 72;
0211
0212     z = fpga_factor[fpga] * 16 + asic * 8 + lhfcal_channel / 8;
0213
0214     return true;
0215 }