Project

General

Profile

Tracker Software Review: Meeting 7

September 8, 2017: 1400 BST

Agenda

  1. Introduction - K. Long
  2. Actions from previous meetings
    • PR:
      1. AD: Handle errors correctly in MINUIT fits
      2. AD: Optimize Chi2 cut for “MINUIT” variant
      3. AD, CH: Check and implement MCS errors
      4. AD: Plot chisq per dof for 'perfect' events
      5. AD, KL: Define multi-track requirements
      6. AD: Compare data/MC efficiencies for 8681 Done
      7. AD: Show distributions with sign per station per tracker Done
      8. AD: Show individual residual distributions to understand details - is the station 1 distribution a binning artifact or a constraint from the fit? Done
    • Track fit:
      1. Repeat study that was done to investigate whether field integration (e.g. using Runge Kutta) is required to compensate for magnetic field non-uniformity
      2. Make study of magnetic-field alignment with a view to establishing whether the present algorithm is sufficiently insensitive to reasonable assumptions of maximum field-misalignment
      3. Study p-value: split up contributions to p-value shape using MC and study shape in data, e.g. is there an error in the resolution per plane, handling of MCS, handing of energy loss, b field..
    • Tracker, MC:
      1. MU/EO: Understand reasons for shape of hits-per-station histogram - low-level
      2. PK: Check dead channel & noise handling in MC.
      3. Tracker group: Summarize, show what has been studied so far for systematics (alignment, scattering, energy loss, field uniformity, etc).
  1. Report on actions: Pattern recognition: A. Dobbs
  2. Report on actions: Track fit: C. Hunt
  3. Next Steps for review
  4. AoB

Dial-in information

http://mice.iit.edu/phonebridge.html


Attendance: KL, CR, PK, AD, CH, MU, DR

Notes from Meeting 5: Aug 16, 2017:

  • AD: Compare data/MC efficiencies for 8681
    • done. 8681, 2.8.5 MC reprocessed with 2.9.1-PatRec,Kalman -- results consistent with 2.9.1 reco
    • Illustrative plots added
  • AD: Optimize Chi2 cut for “MINUIT” variant
    • pushing cuts lower gives higher efficiencies, though purity may be an issue. Noted that critical parameter is efficiency
    • Q: how loose should the LSQ cuts be? AD points to loosest cut in table from previous meeting.
  • AD: Show residual distributions with sign per station per tracker
    • sent to email list, added to meeting wiki
  • AD: Plot chisq per dof for 'perfect’ events
    • stands
  • Define multi-track requirements
    • stands. Need discussion and inputs on how to define
  • Adding errors:
    • for circle fit: suggestion is to define a mean error and apply to every station rather than station-by-station errors.
    • Q: how to define mean error? based on the geometric mean of the errors?
    • Q: does taking the chi^2/mean_error^2 bias the pr fit?
    • errors on the resulting circle fit parameters will feed into the longitudinal fit
  • Discussion about aperture cut
    • CR notes that aperture cut should be standardized by tracker group. AD notes that the cut in PR is loose -- 150.
    • PK suggests a flag to indicate if track passed or failed an aperture cut
    • All agree that fiducial cut should be standardized,
    • Q: does the current aperture cut throw out any potential tracks? If it does, cut should be loosened, but consensus seems to be that tracking should not do any "analysis-style" cuts but rather optimize track-finding.
  • Effect of magnetic field non-uniformity, alignment on track fit -- CH had issues reading mag field map from Joe Langland's program. CH will check against comsol map for 7469.

Code snippet of pattern recognition residual analysis

bool AnalyserTrackerMCPRResiduals::analyse(MAUS::ReconEvent* const aReconEvent,
                                           MAUS::MCEvent* const aMCEvent) {
  if (!aReconEvent || !aMCEvent) {
    return false;
  }

  // Find a hit from a muon in the tracker references plane for each tracker
  MAUS::SciFiHit* tku_ref_hit = nullptr;
  MAUS::SciFiHit* tkd_ref_hit = nullptr;
  for (auto&& hitref : *(aMCEvent->GetSciFiHits())) {
    MAUS::SciFiHit* hit = &hitref;
    int pid = hit->GetParticleId();
    if (pid != 13 && pid != -13)
      continue;
    if (hit->GetChannelId()->GetStationNumber() == 1 && \
        hit->GetChannelId()->GetPlaneNumber() == 0) {
      if (hit->GetChannelId()->GetTrackerNumber() == 0) {
        tku_ref_hit = hit;
      } else if (hit->GetChannelId()->GetTrackerNumber() == 1) {
        tkd_ref_hit = hit;
      }
    }
    if (tku_ref_hit && tkd_ref_hit)
      break;
  }

  if (!tku_ref_hit || !tkd_ref_hit) {
    return false;
  }

  // Access the SciFi Event
  MAUS::SciFiEvent* sfevt = aReconEvent->GetSciFiEvent();
  if (!sfevt) {
    return false;
  }

  // Loop over all the scifi tracks in this event
  int nTkURecTracks = 0;
  int nTkDRecTracks = 0;
  MAUS::SciFiHelicalPRTrack* tku_trk = nullptr;
  MAUS::SciFiHelicalPRTrack* tkd_trk = nullptr;
  for (auto trk : sfevt->helicalprtracks()) {
    if (trk->get_tracker() == 0) {
      ++nTkURecTracks;
      tku_trk = trk;
    } else if (trk->get_tracker() == 1) {
      ++nTkDRecTracks;
      tkd_trk = trk;
    }
  }

  // Require 1 recon track in each tracker
  if (nTkURecTracks != 1 || nTkDRecTracks != 1) {
    return false;
  }

  if (tku_trk->get_spacepoints_pointers().size() != 5 ||
      tkd_trk->get_spacepoints_pointers().size() != 5) {
    return false;
  }

  // At this point we should an event that is good to analyse
  // Calculate the recon postion at the reference station, and the recon momentum (average)
  double tku_x = 0.0;
  double tku_y = 0.0;
  for (auto sp : tku_trk->get_spacepoints_pointers()) {
    if (sp->get_station() == 1) {
      tku_x = sp->get_position().x();
      tku_y = sp->get_position().y();
    }
  }

  double tkd_x = 0.0;
  double tkd_y = 0.0;
  for (auto sp : tkd_trk->get_spacepoints_pointers()) {
    if (sp->get_station() == 1) {
      tkd_x = sp->get_position().x();
      tkd_y = sp->get_position().y();
    }
  }

  double tku_dx = tku_x + tku_ref_hit->GetPosition().x();
  double tku_dy = tku_y - tku_ref_hit->GetPosition().y();

  double tku_ptmc = sqrt(tku_ref_hit->GetMomentum().x()*tku_ref_hit->GetMomentum().x() +
                         tku_ref_hit->GetMomentum().y()*tku_ref_hit->GetMomentum().y());
  double tku_ptrec = 0.3*mBfield*tku_trk->get_R(); // 0.3 comes from mom being in MeV and rad in mm
  double tku_dpt = tku_ptrec - tku_ptmc;

  double tku_pzrec = tku_ptrec / tku_trk->get_dsdz();
  double tku_dpz =  tku_pzrec + tku_ref_hit->GetMomentum().z();

  double tkd_dx = tkd_x - tkd_ref_hit->GetPosition().x();
  double tkd_dy = tkd_y - tkd_ref_hit->GetPosition().y();

  double tkd_ptmc = sqrt(tkd_ref_hit->GetMomentum().x()*tkd_ref_hit->GetMomentum().x() +
                         tkd_ref_hit->GetMomentum().y()*tkd_ref_hit->GetMomentum().y());
  double tkd_ptrec = 0.3*mBfield*tkd_trk->get_R(); // 0.3 comes from mom being in MeV and rad in mm
  double tkd_dpt = tkd_ptrec - tkd_ptmc;

  double tkd_pzrec = tkd_ptrec / tkd_trk->get_dsdz();
  double tkd_dpz = tkd_pzrec + tkd_ref_hit->GetMomentum().z();

  // Fill the histograms
  mHTkUMCPositionX->Fill(tku_ref_hit->GetPosition().x());
  mHTkUMCPositionY->Fill(tku_ref_hit->GetPosition().y());
  mHTkUMCMomentumT->Fill(tku_ptmc);
  mHTkUMCMomentumZ->Fill(tku_ref_hit->GetMomentum().z());

  mHTkURecPositionX->Fill(tku_x);
  mHTkURecPositionY->Fill(tku_y);
  mHTkURecMomentumT->Fill(tku_ptrec);
  mHTkURecMomentumZ->Fill(tku_pzrec);

  mHTkUPositionResidualsX->Fill(tku_dx);
  mHTkUPositionResidualsY->Fill(tku_dy);
  mHTkUMomentumResidualsT->Fill(tku_dpt);
  mHTkUMomentumResidualsZ->Fill(tku_dpz);

  mHTkDMCPositionX->Fill(tkd_ref_hit->GetPosition().x());
  mHTkDMCPositionY->Fill(tkd_ref_hit->GetPosition().y());
  mHTkDMCMomentumT->Fill(tkd_ptmc);
  mHTkDMCMomentumZ->Fill(tkd_ref_hit->GetMomentum().z());

  mHTkDRecPositionX->Fill(tkd_x);
  mHTkDRecPositionY->Fill(tkd_y);
  mHTkDRecMomentumT->Fill(tkd_ptrec);
  mHTkDRecMomentumZ->Fill(tkd_pzrec);

  mHTkDPositionResidualsX->Fill(tkd_dx);
  mHTkDPositionResidualsY->Fill(tkd_dy);
  mHTkDMomentumResidualsT->Fill(tkd_dpt);
  mHTkDMomentumResidualsZ->Fill(tkd_dpz);

  return true;
}