THcScintillatorPlane.cxx

	  timec_pos = tdc_pos*fScinTdcToTime
	    - fHodoPosInvAdcOffset[index]
	    - fHodoPosInvAdcAdc[index]/TMath::Sqrt(TMath::Max(20.0*.020,adcint_pos));
	  timec_neg = tdc_neg*fScinTdcToTime
	    - fHodoNegInvAdcOffset[index]
	    - fHodoNegInvAdcAdc[index]/TMath::Sqrt(TMath::Max(20.0*.020,adcint_neg));
	} else {		// FADC style
	  timec_pos =  tdc_pos*fScinTdcToTime -tw_corr_pos + fHodo_LCoeff[index];
	  timec_neg =  tdc_neg*fScinTdcToTime -tw_corr_neg- 2*fHodoCableFit[index] + fHodo_LCoeff[index];
	}

 	Double_t TWCorrDiff = fGoodNegTdcTimeWalkCorr.at(padnum-1) - 2*fHodoCableFit[index] - fGoodPosTdcTimeWalkCorr.at(padnum-1); 
     	
        Double_t fHitDistCorr = 0.5*TWCorrDiff*fHodoVelFit[index];  

	/*Debug*/
	//cout << "*****************" << endl;
	//cout << "fPlNum: " << fPlaneNum << endl;
	//cout << "*****************" << endl;
	//cout << "Paddle index: " << index << endl;
	//cout << "pos_sigma: " << fHodoSigma[index] << endl;
	//cout << "fZPos: " << fZpos << endl;
	//cout << "fDzPos: " << fDzpos << endl;
	//cout << "Zcorr = fZpos+(index%2)*fDzpos = " << fZpos+(index%2)*fDzpos << endl;

	//cout << Form("****fHodo_LCoeff[%d]", index) << fHodo_LCoeff[index] << endl;
	//cout << Form("****fHodoCableFit[%d]", index) << fHodoCableFit[index] << endl;
	//cout << Form("****fHodoVelFit[%d]", index) << fHodoVelFit[index] << endl;
	//cout << Form("****c1_Pos/Neg[%d]", index) <<  fHodoPos_c1[index] << " / " << fHodoNeg_c1[index] << endl;
	//cout << Form("****c2_Pos/Neg[%d]", index) <<  fHodoPos_c2[index] << " / " << fHodoNeg_c2[index] << endl;
	//cout << "TW Corr val. +/-: " << tw_corr_pos << " / " << tw_corr_neg << endl; 
	//cout << "TW UnCorr +/-: " << fGoodPosTdcTimeUnCorr.at(padnum-1) << " / " << fGoodNegTdcTimeUnCorr.at(padnum-1) << endl;
	//cout << "TW Corr +/-: " << fGoodPosTdcTimeWalkCorr.at(padnum-1) << " / " << fGoodNegTdcTimeWalkCorr.at(padnum-1) << endl;
	

	fGoodDiffDistTrack.at(index) =  fHitDistCorr;
	// Find hit position using ADCs
	// If postime larger, then hit was nearer negative side.
	
	Double_t vellight=fHodoVelLight[index]; //read from hodo_cuts.param, where it is set fixed to 15.0 
	//Double_t vellight=fHodoVelFit[index];   //use scin prop vel. values from hodo_calibVp_run#.param file
	
	Double_t dist_from_center=0.5*(timec_neg-timec_pos)*vellight;
	Double_t scint_center=0.5*(fPosLeft+fPosRight);
	Double_t hit_position=scint_center+dist_from_center;
	hit_position=TMath::Min(hit_position,fPosLeft);
	hit_position=TMath::Max(hit_position,fPosRight);
	Double_t scin_corrected_time, postime, negtime;
	if(fTofUsingInvAdc) {
	  timec_pos -= (fPosLeft-hit_position)/
	    fHodoPosInvAdcLinear[index];
	  timec_neg -= (hit_position-fPosRight)/
	    fHodoNegInvAdcLinear[index];
	  scin_corrected_time = 0.5*(timec_pos+timec_neg);
	  if (fCosmicFlag) {
	    postime = timec_pos + (fZpos+(index%2)*fDzpos)/(29.979*fBetaNominal);
	    negtime = timec_neg + (fZpos+(index%2)*fDzpos)/(29.979*fBetaNominal);
	  } else {
	    postime = timec_pos - (fZpos+(index%2)*fDzpos)/(29.979*fBetaNominal);
	    negtime = timec_neg - (fZpos+(index%2)*fDzpos)/(29.979*fBetaNominal);
	  }
	} else {
	  scin_corrected_time = 0.5*(timec_neg+timec_pos);  // add constants for each paddle, 25ns, 25 + zpos, . . . //remove propagation time
	  timec_pos= scin_corrected_time;    
	  timec_neg= scin_corrected_time;	 
	  if (fCosmicFlag) {
	    postime = timec_pos + (fZpos+(index%2)*fDzpos)/(29.979*fBetaNominal);
	    negtime = timec_neg + (fZpos+(index%2)*fDzpos)/(29.979*fBetaNominal);
	  } else {
	    postime = timec_pos - (fZpos+(index%2)*fDzpos)/(29.979*fBetaNominal);
	    negtime = timec_neg - (fZpos+(index%2)*fDzpos)/(29.979*fBetaNominal);
	  }
	}
        ((THcHodoHit*) fHodoHits->At(fNScinHits))->SetPaddleCenter(fPosCenter[index]);
	((THcHodoHit*) fHodoHits->At(fNScinHits))->SetCorrectedTimes(timec_pos,timec_neg,         
								     postime, negtime,
								     scin_corrected_time);
	((THcHodoHit*) fHodoHits->At(fNScinHits))->SetPosADCpeak(adcamp_pos); // need for new TWCOrr
	((THcHodoHit*) fHodoHits->At(fNScinHits))->SetNegADCpeak(adcamp_neg); // need for new TWCOrr
	fGoodPosTdcTimeCorr.at(padnum-1) = timec_pos;
	fGoodNegTdcTimeCorr.at(padnum-1) = timec_neg;
	fGoodPosTdcTimeTOFCorr.at(padnum-1) = postime;
	fGoodNegTdcTimeTOFCorr.at(padnum-1) = negtime;
      } else {
	Double_t timec_pos,timec_neg;
        timec_pos=tdc_pos;
        timec_neg=tdc_neg;
	if(btdcraw_pos) {
	  if(fTofUsingInvAdc) {
	    timec_pos = tdc_pos*fScinTdcToTime
	      - fHodoPosInvAdcOffset[index]
	      - fHodoPosInvAdcAdc[index]/TMath::Sqrt(TMath::Max(20.0*.020,adcint_pos));
	  } else {		// FADC style
	  timec_pos =  tdc_pos*fScinTdcToTime -tw_corr_pos + fHodo_LCoeff[index];
	  }
	}
	if(btdcraw_neg) {
	  if(fTofUsingInvAdc) {
	    timec_neg = tdc_neg*fScinTdcToTime
	      - fHodoNegInvAdcOffset[index]
	      - fHodoNegInvAdcAdc[index]/TMath::Sqrt(TMath::Max(20.0*.020,adcint_neg));
	  } else {		// FADC style
	  timec_neg =  tdc_neg*fScinTdcToTime -tw_corr_neg- 2*fHodoCableFit[index] + fHodo_LCoeff[index];
	  }
	}
        ((THcHodoHit*) fHodoHits->At(fNScinHits))->SetPaddleCenter(fPosCenter[index]);
	((THcHodoHit*) fHodoHits->At(fNScinHits))->SetCorrectedTimes(timec_pos,timec_neg,
								     timec_pos,timec_neg,0.0);
	((THcHodoHit*) fHodoHits->At(fNScinHits))->SetNegADCpeak(adcamp_neg); // needed for new TWCOrr
	((THcHodoHit*) fHodoHits->At(fNScinHits))->SetPosADCpeak(adcamp_pos); // needed for new TWCOrr
	fGoodPosTdcTimeCorr.at(padnum-1) = timec_pos;
	fGoodNegTdcTimeCorr.at(padnum-1) = timec_neg;
	fGoodPosTdcTimeTOFCorr.at(padnum-1) = timec_pos;
	fGoodNegTdcTimeTOFCorr.at(padnum-1) = timec_neg;
      }
      fNScinHits++;		// One or more good time counter
      //
    }
    ihit++;			// Raw hit counter
  }
  //  cout << "THcScintillatorPlane: ihit = " << ihit << endl;

  return(ihit);
}

//_____________________________________________________________________________
Int_t THcScintillatorPlane::AccumulatePedestals(TClonesArray* rawhits, Int_t nexthit)
{
  /*! \brief Extract the data for this plane from raw hit list THcRawHodoHit, accumulating into arrays for calculating pedestals.
   *
   * - Loop through raw data for scintillator plane
   */
  Int_t nrawhits = rawhits->GetLast()+1;
  // cout << "THcScintillatorPlane::AcculatePedestals " << fPlaneNum << " " << nexthit << "/" << nrawhits << endl;

  Int_t ihit = nexthit;
  while(ihit < nrawhits) {
    THcRawHodoHit* hit = (THcRawHodoHit *) rawhits->At(ihit);
    if(hit->fPlane > fPlaneNum) {
      break;
    }
    Int_t element = hit->fCounter - 1; // Should check if in range
    Int_t adcpos = hit->GetRawAdcHitPos().GetPulseIntRaw();
    Int_t adcneg = hit->GetRawAdcHitNeg().GetPulseIntRaw();

    if(adcpos <= fPosPedLimit[element]) {
      fPosPedSum[element] += adcpos;
      fPosPedSum2[element] += adcpos*adcpos;
      fPosPedCount[element]++;
      if(fPosPedCount[element] == fMinPeds/5) {
	fPosPedLimit[element] = 100 + fPosPedSum[element]/fPosPedCount[element];
      }
    }
    if(adcneg <= fNegPedLimit[element]) {
      fNegPedSum[element] += adcneg;
      fNegPedSum2[element] += adcneg*adcneg;
      fNegPedCount[element]++;
      if(fNegPedCount[element] == fMinPeds/5) {
	fNegPedLimit[element] = 100 + fNegPedSum[element]/fNegPedCount[element];
      }
    }
    ihit++;
  }

  fNPedestalEvents++;

  return(ihit);
}

//_____________________________________________________________________________
void THcScintillatorPlane::CalculatePedestals( )
{
  /*! \brief   Calculate pedestals from arrays made in THcScintillatorPlane::AccumulatePedestals
   *
   * - Calculate pedestals from arrays made in THcScintillatorPlane::AccumulatePedestals
   * - In old fortran ENGINE code, a comparison was made between calculated pedestals and the pedestals read in by the FASTBUS modules for zero supression. This is not implemented.
   */
  for(UInt_t i=0; i<fNelem;i++) {

    // Positive tubes
    fPosPed[i] = ((Double_t) fPosPedSum[i]) / TMath::Max(1, fPosPedCount[i]);
    fPosThresh[i] = fPosPed[i] + 15;

    // Negative tubes
    fNegPed[i] = ((Double_t) fNegPedSum[i]) / TMath::Max(1, fNegPedCount[i]);
    fNegThresh[i] = fNegPed[i] + 15;

    //    cout <<"Pedestals "<< i+1 << " " << fPosPed[i] << " " << fNegPed[i] << endl;
  }
  //  cout << " " << endl;

}

//_____________________________________________________________________________
void THcScintillatorPlane::InitializePedestals( )
{
  /*! \brief   called by THcScintillatorPlane::ReadDatabase
   *
   * - Initialize variables used in  THcScintillatorPlane::AccumulatePedestals and THcScintillatorPlane::CalculatePedestals
   * - Minimum number of pedestal events needed for calculation, fMinPeds, hadrcoded to 500
   */
  fNPedestalEvents = 0;
  fMinPeds = 500; 		// In engine, this is set in parameter file
  fPosPedSum = new Int_t [fNelem];
  fPosPedSum2 = new Int_t [fNelem];
  fPosPedLimit = new Int_t [fNelem];
  fPosPedCount = new Int_t [fNelem];
  fNegPedSum = new Int_t [fNelem];
  fNegPedSum2 = new Int_t [fNelem];
  fNegPedLimit = new Int_t [fNelem];
  fNegPedCount = new Int_t [fNelem];

  fPosPed = new Double_t [fNelem];
  fNegPed = new Double_t [fNelem];
  fPosThresh = new Double_t [fNelem];
  fNegThresh = new Double_t [fNelem];
  for(UInt_t i=0;i<fNelem;i++) {
    fPosPedSum[i] = 0;
    fPosPedSum2[i] = 0;
    fPosPedLimit[i] = 1000;	// In engine, this are set in parameter file
    fPosPedCount[i] = 0;
    fNegPedSum[i] = 0;
    fNegPedSum2[i] = 0;
    fNegPedLimit[i] = 1000;	// In engine, this are set in parameter file
    fNegPedCount[i] = 0;
  }
}
//____________________________________________________________________________
ClassImp(THcScintillatorPlane)
////////////////////////////////////////////////////////////////////////////////