Back to home page

sPhenix code displayed by LXR

 
 

    


File indexing completed on 2025-08-05 08:19:43

0001 //
0002 // Make the MBD LUT's
0003 //
0004 #include <cmath>
0005 #include <TMath.h>
0006 #include <mbd/MbdCalib.h>
0007 #include <mbd/MbdGeomV1.h>
0008 #include <mbd/MbdDefs.h>
0009 #include <cstdint>
0010 
0011 
0012 #if defined(__CLING__)
0013 R__LOAD_LIBRARY(libmbd.so)
0014 #endif
0015 
0016 MbdGeom *mbdgeom{nullptr};
0017 MbdCalib *mcal{nullptr};
0018 
0019 
0020 float tdclut[128][1024] = {};
0021 uint16_t adclut[128][1024] = {};
0022 uint16_t slewlut[128][4096] = {};
0023 double step = 0;                   // ns/LL1 tdc count
0024 double slewcorr_step[128] = {};
0025 
0026 void make_generic_adclut();
0027 void make_generic_slewlut();
0028 void make_adclut(const char *scorr_fname = "mbd_slewcorr.calib");
0029 void make_slewlut(const char *scorr_fname = "mbd_slewcorr.calib");
0030 
0031 void make_generic_adclut()
0032 {
0033   // write adc lut
0034   ofstream adclutfile("mbdadc.lut");
0035   for (int ipmt=0; ipmt<MbdDefs::MBD_N_PMT; ipmt++)
0036   {
0037     adclutfile << "adc " << ipmt << endl;
0038     for (int iaddr=0; iaddr<1024; iaddr++)
0039     {
0040       uint16_t qsum7 = iaddr>>3;
0041       uint16_t slew3 = iaddr>>7;
0042       adclut[ipmt][iaddr] = (slew3<<7) + qsum7;   // assuming slew3 are high order bits
0043 
0044       adclutfile << slew3*1000 + qsum7 << "\t";
0045 
0046       if ( iaddr%8 == 7 )
0047       {
0048         adclutfile << endl;
0049       }
0050     }
0051   }
0052   adclutfile.close();
0053 }
0054 
0055 void make_generic_slewlut()
0056 {
0057   // write slew lut
0058   ofstream slewlutfile("mbdslew.lut");
0059   for (int ipmt=0; ipmt<MbdDefs::MBD_N_PMT; ipmt++)
0060   {
0061     slewlutfile << "slew " << ipmt << endl;
0062     for (int islew=0; islew<8; islew++)   // 3-bit
0063     {
0064       uint16_t saddr = islew;
0065       for (int itdc=0; itdc<512; itdc++)  // 9-bit
0066       {
0067         uint16_t taddr = itdc;
0068         if ( taddr!=0 )
0069         {
0070           slewlut[ipmt][itdc] = 1000 + taddr;
0071           // 1000 represents the active bit in the lut file
0072           // when loading into the actual MBD ADC LUT, we should
0073           // we should set bit 10 = 1.
0074         }
0075         else
0076         {
0077           slewlut[ipmt][itdc] = taddr;
0078         }
0079 
0080         slewlutfile << slewlut[ipmt][itdc] << "\t";
0081 
0082         if ( itdc%8 == 7 )
0083         {
0084           slewlutfile << endl;
0085         }
0086       }
0087     }
0088   }
0089   slewlutfile.close();
0090 }
0091 
0092 void make_adclut(const char *scorr_fname = "mbd_slewcorr.calib")
0093 {
0094 
0095 
0096   // write adc lut
0097   ofstream adclutfile("mbdadc.lut");
0098   for (int ipmt=0; ipmt<MbdDefs::MBD_N_PMT; ipmt++)
0099   {
0100     int feech = mbdgeom->get_feech(ipmt,0);
0101 
0102     // get the max slew correction, which we set at ADC=40
0103     double max_slew = mcal->get_scorr(feech,40);
0104     slewcorr_step[ipmt] = max_slew/8;
0105     int nsteps = 0;
0106  
0107     adclutfile << "adc " << ipmt << endl;
0108     for (int iaddr=0; iaddr<1024; iaddr++)
0109     {
0110       uint16_t qsum7 = iaddr>>3;
0111 
0112       uint16_t actual_adc = (iaddr<<4) + 8; // 8 = 2^4/2, so actual_adc is at midpoint
0113       double scorr = mcal->get_scorr(feech,actual_adc);
0114       int scorr_step = static_cast<int>( scorr/slewcorr_step[ipmt] );
0115       if ( scorr_step < 0 )
0116       {
0117         scorr_step = 0;
0118       }
0119       else if ( scorr_step > 7 )
0120       {
0121         scorr_step = 7;
0122       }
0123 
0124       uint16_t slew3 = 7 - scorr_step;
0125 
0126       if ( ipmt==0 )
0127       {
0128         cout << "adcslewbin " << iaddr << "\t" << scorr << "\t" << scorr_step << "\t" << slew3 << endl;
0129       }
0130 
0131       adclut[ipmt][iaddr] = (slew3<<7) + qsum7;   // assuming slew3 are high order bits
0132 
0133       adclutfile << slew3*1000 + qsum7 << "\t";
0134 
0135       if ( iaddr%8 == 7 )
0136       {
0137         adclutfile << endl;
0138       }
0139     }
0140   }
0141   adclutfile.close();
0142 }
0143 
0144 void make_slewlut(const char *scorr_fname = "mbd_slewcorr.calib")
0145 {
0146 
0147   // write slew lut
0148   ofstream slewlutfile("mbdslew.lut");
0149   for (int ipmt=0; ipmt<MbdDefs::MBD_N_PMT; ipmt++)
0150   {
0151     slewlutfile << "slew " << ipmt << endl;
0152     for (int islew=0; islew<8; islew++)   // 3-bit
0153     {
0154       uint16_t saddr = islew;
0155       for (int itdc=0; itdc<512; itdc++)  // 9-bit
0156       {
0157         int taddr = itdc;
0158         if ( taddr!=0 )
0159         {
0160           // apply slew correction
0161           taddr -= static_cast<int>( (7-islew)*(slewcorr_step[ipmt]/step) );
0162           if ( ipmt==0 )
0163           {
0164             cout << "slewxxx " << islew << "\t" << itdc << "\t" << (7-islew) << "\t" << slewcorr_step[ipmt] << "\t" << step << "\t" << (slewcorr_step[ipmt]/step) << "\t" << taddr << endl;
0165           }
0166 
0167           if ( taddr<0 )
0168           {
0169             taddr = 1;
0170           }
0171           slewlut[ipmt][itdc] = 1000 + taddr;
0172           // 1000 represents the active bit in the lut file
0173           // when loading into the actual MBD ADC LUT, we should
0174           // we should set bit 10 = 1.
0175         }
0176         else  // 0 is no-hit for TDC
0177         {
0178           slewlut[ipmt][itdc] = taddr;
0179         }
0180 
0181         slewlutfile << slewlut[ipmt][itdc] << "\t";
0182 
0183         if ( itdc%8 == 7 )
0184         {
0185           slewlutfile << endl;
0186         }
0187       }
0188     }
0189   }
0190   slewlutfile.close();
0191 }
0192 
0193 void make_mbd_l1lut(const char *tcorr_fname = "mbd_timecorr.calib",
0194     const char *t0_fname = "mbd_tt_t0.calib",
0195     const char *scorr_fname = "mbd_slewcorr.calib",
0196     const char *tpranges_fname = "mbd_tpscan.txt")
0197 {
0198   mbdgeom = new MbdGeomV1();
0199 
0200   mcal = new MbdCalib();
0201   mcal->Download_TimeCorr( tcorr_fname );
0202   mcal->Download_TTT0( t0_fname );
0203   mcal->Download_SlewCorr( scorr_fname );
0204 
0205   ifstream tprangesfile(tpranges_fname);
0206   int temp_feech;
0207   double mintdc[128] = {0};
0208   double maxtdc[128] = {0};
0209   double mintime[128] = {0};
0210   double maxtime[128] = {0};
0211   double dtdc[128] = {0};
0212   double dtime[128] = {0};
0213 
0214   // Find the overlapping good time ranges for all the channels
0215   TGraph *g_mintdc = new TGraph();
0216   TGraph *g_maxtdc = new TGraph();
0217   g_mintdc->SetName("g_mintdc");
0218   g_maxtdc->SetName("g_maxtdc");
0219   TGraph *g_mintimes = new TGraph();
0220   TGraph *g_maxtimes = new TGraph();
0221   g_mintimes->SetName("g_mintdc");
0222   g_maxtimes->SetName("g_maxtdc");
0223 
0224   while ( tprangesfile >> temp_feech )
0225   {
0226     int pmtch = (double)mbdgeom->get_pmt(temp_feech);
0227 
0228     tprangesfile >> mintdc[pmtch] >> maxtdc[pmtch] >> dtdc[pmtch] >> mintime[pmtch] >> maxtime[pmtch] >> dtime[pmtch];
0229 
0230     //cout << "maxtime " << pmtch << "\t" << maxtime[pmtch] << "\t" << mcal->get_tcorr( temp_feech, 16 ) << endl;;
0231     maxtime[pmtch] = mcal->get_tcorr( temp_feech, 16 ) - mcal->get_tt0( pmtch );
0232 
0233     cout << "mintime " << pmtch << "\t" << mintime[pmtch] << "\t" << mcal->get_tcorr( temp_feech, static_cast<int>(std::round(maxtdc[pmtch])) ) << endl;;
0234     cout << mcal->get_tt0(pmtch) << endl;
0235     mintime[pmtch] = mcal->get_tcorr( temp_feech, static_cast<int>(round(maxtdc[pmtch])) ) - mcal->get_tt0( pmtch );
0236 
0237     // correct for bad channels
0238     if ( fabs(mcal->get_tt0( pmtch ))>100. )
0239     {
0240       maxtime[pmtch] = mcal->get_tcorr( temp_feech, 16 ) - 3.0;
0241       mintime[pmtch] = mcal->get_tcorr( temp_feech, static_cast<int>(round(maxtdc[pmtch])) ) - 3.0;
0242     }
0243 
0244     int n = g_mintdc->GetN();
0245     g_mintdc->SetPoint( n, pmtch, mintdc[pmtch] );
0246 
0247     n = g_maxtdc->GetN();
0248     g_maxtdc->SetPoint( n, pmtch, maxtdc[pmtch] );
0249 
0250     n = g_mintimes->GetN();
0251     g_mintimes->SetPoint( n, pmtch, mintime[pmtch] );
0252 
0253     n = g_maxtimes->GetN();
0254     g_maxtimes->SetPoint( n, pmtch, maxtime[pmtch] );
0255   }
0256 
0257   int n = g_mintimes->GetN();
0258   double lut_mintime = TMath::MaxElement(n,g_mintimes->GetY());
0259 
0260   n = g_maxtimes->GetN();
0261   double lut_maxtime = TMath::MinElement(n,g_maxtimes->GetY());
0262 
0263   cout << "max-min time: " << lut_mintime << "\t" << lut_maxtime << "\t" << lut_maxtime - lut_mintime << endl;
0264 
0265   g_mintimes->SetMarkerStyle(20);
0266   g_mintimes->SetMarkerColor(3);
0267   g_maxtimes->SetMarkerStyle(20);
0268   g_maxtimes->SetMarkerColor(4);
0269   g_mintimes->SetMinimum(lut_mintime*1.5);
0270   g_mintimes->SetMaximum(lut_maxtime*1.2);
0271   g_mintimes->Draw("ap");
0272   g_maxtimes->Draw("p");
0273 
0274   double max_range = lut_maxtime - lut_mintime;
0275   step = max_range/511.;           // ns/LL1 tdc count
0276   double zstep = (2.0/29.9792458)*(1.0/step);    
0277   cout << "L1 trig time unit is " << step << " ns/count" << endl;
0278   cout << "L1 trig time unit is " << zstep << " count/cm" << endl;
0279 
0280   for (int ifeech=0; ifeech<MbdDefs::MBD_N_FEECH; ifeech++)
0281   {
0282     if ( mbdgeom->get_type(ifeech) == 1 ) continue; // skip q-ch's
0283     int pmtch = mbdgeom->get_pmt(ifeech);
0284 
0285     for (int iaddr=0; iaddr<1024; iaddr++ )
0286     {
0287       // could average over the 16 values
0288       int tdc = (iaddr<<4) + 8;
0289       float true_time = mcal->get_tcorr(ifeech,tdc) - mcal->get_tt0(pmtch) - lut_mintime;
0290 
0291       if ( mcal->get_tt0(pmtch)<-100. )
0292       {
0293         // set bad channels to 0;
0294         //tdclut[pmtch][iaddr] = 0;
0295         //continue;
0296 
0297         true_time = mcal->get_tcorr(ifeech,tdc) - 3.0 - lut_mintime;
0298       }
0299 
0300       if ( true_time<0. || true_time>max_range || isnan(true_time) || true_time>22.5 )
0301       {
0302         tdclut[pmtch][iaddr] = 0;
0303       }
0304       else
0305       {
0306         if ( isnan(true_time) )
0307         {
0308           cout << "isnan " << ifeech << " " << mcal->get_tcorr(ifeech,tdc) 
0309             << "\t" << mcal->get_tt0(pmtch) << endl;;
0310         }
0311 
0312         //tdclut[pmtch][iaddr] = static_cast<int>( round(true_time/step) );
0313         tdclut[pmtch][iaddr] = round(true_time/step);
0314       }
0315 
0316       if ( ifeech==0 && (iaddr<24 || iaddr>(1023-96) ) )
0317       {
0318         if (iaddr%8==0) cout << iaddr << "\t";
0319         cout << tdclut[pmtch][iaddr] << "\t";   // be sure to set upper bit to 1 for non-zero values when writing to actual LUT in ADC
0320         if (iaddr%8==7) cout << endl;
0321       }
0322     }
0323   }
0324 
0325   // write tdc lut
0326   ofstream tdclutfile("mbdtdc.lut");
0327   for (int ifeech=0; ifeech<MbdDefs::MBD_N_FEECH; ifeech++)
0328   {
0329     if ( mbdgeom->get_type(ifeech) == 1 ) continue; // skip q-ch's
0330     int pmtch = mbdgeom->get_pmt(ifeech);
0331 
0332     tdclutfile << "tdc " << pmtch << endl;
0333     for (int iaddr=0; iaddr<1024; iaddr++)
0334     {
0335 
0336       tdclutfile << tdclut[pmtch][iaddr] << "\t";
0337 
0338       if ( iaddr%8 == 7 )
0339       {
0340         tdclutfile << endl;
0341       }
0342     }
0343   }
0344   tdclutfile.close();
0345 
0346 
0347   make_adclut();
0348   make_slewlut();
0349 }
0350