Page MenuHomeHEPForge

Implementation of higher-spin K-matrix
Needs ReviewPublic

Authored by johndan on Nov 26 2020, 12:36 PM.

Details

Reviewers
tlatham
jback
Summary

This differential contains commits intended to allow use of K-matrix implementation in non S-wave scenarios.
(it builds on D39)

Test Plan

Trial fits to a B+ 3-body decay

Diff Detail

Repository
rLAURA laura
Branch
johndan-KMatrix_D0K_higherSpin
Lint
No Linters Available
Unit
No Unit Test Coverage
Build Status
Buildable 119
Build 119: arc lint + arc unit

Event Timeline

There are a very large number of changes, so older changes are hidden. Show Older Changes
johndan retitled this revision from Allow floating of parameters in the k-matrix to Implementation of higher-spin K-matrix.Nov 26 2020, 4:19 PM
johndan edited the summary of this revision. (Show Details)
johndan updated this revision to Diff 145.Nov 26 2020, 4:58 PM

Fix typo in LauKMatrixPropagator.hh

johndan updated this revision to Diff 146.Nov 27 2020, 2:07 PM
  • Typo in LauKMatrixPropagator.cc
  • Updates to K-matrix propagator to use refined barrier factor
  • Clean up indentation
  • Parallel updates to KMatrixPropagator.hh
  • Cleanup indentation
  • Revert "Allow slope of exponential non-resonant shape to vary positive (for low-stats fits)"
johndan updated this revision to Diff 150.Nov 30 2020, 1:46 PM
  • Correct name of s0prod LauParameter, thus distinguishing it from s0Scatt
  • Need to name scattering parameter for it to be floated
jback requested changes to this revision.Dec 17 2020, 10:42 PM

I've had a first look at the code and have provided some suggestions for changes that I think are needed.

inc/LauKMatrixPropagator.hh
65

We should add another parameter in the constructor to define the spin L value (the last argument). By default L = 0, and for specific channels it could be set appropriately by the code. This integer should probably be another enumeration, something like:

enum KMatrixSpin {SWave, PWave, DWave, MaxSpin};

467

I think it will be much better to define another enumeration that specifies the spin of the K-matrix, irrespective of the channels. That way we can use non-zero spins for other modes such as pi pi etc.

enum KMatrixSpin {SWave, PWave, DWave, MaxSpin};

I would change D0K_Pwave to D0K etc. and then we specify the spin via the constructor. Would we never use S-wave for D0K?

src/LauKMatrixPropagator.cc
236

We should reuse the gamma matrix product a few lines earlier (set it as another matrix), since we need to reduce matrix multiplications as much as possible!

874

I think we can remove the if statement and just use:

gamma = this->calcGamma(iChannel, s, L_);

where L_ is the spin integer from the constructor. This code will then be general for any channel.

909

We could use TMath::Sqrt( fabs(sqrtTerm) ) to allow cases below threshold.

913

We need (a/R^2) instead of (1/R) here. We should set default values of a = 0, 1 and 3 (which could be an optional parameter that can be varied) for L = 0, 1, and 2, respectively.

Some analyses don't include the denominator term, so we could make this optional (perhaps via a boolean passed to the function) to allow systematic comparisons, but will be enabled by default. Don't have to do this now.

932

Check spin enumeration integer again.

979

"D0K" and "Dstar0K" channel enumerations

This revision now requires changes to proceed.Dec 17 2020, 10:42 PM
jback added inline comments.Dec 18 2020, 1:56 AM
inc/LauKMatrixPropagator.hh
65

Instead of defining the KMatrixSpin enum (mentioned also below) we just need to set a spin integer L_.

johndan marked 8 inline comments as done.Tue, Jan 5, 3:21 PM

Thank you for the comments @jback . I've implemented them in most cases, or left a reply in the few where further input from you is needed.

General question, @tlatham : should I cherry-pick commits to bring this differential in line with the master (e.g. c38fa2a34f4a772754d1e8f610def10747d95104 - LauFormulaPar fix)?

src/LauKMatrixPropagator.cc
874

I agree that calcGamma should be general for any channel, so no need to check the phasespace index here.

But shouldn't we retain the if statement in order to avoid calling "calcGamma" (and it's internal "if" statements) if the spin is zero?

913

Absolutely. By 'parameter' and 'varied', do you mean a constant (that can be set by the user, but is fixed during the fit) or a LauParameter?

Concerning the denominator, I've included a function to ignore this if needed.

932

Yes. For L==2 should the function be simply "3 cos^2 - 1"? Or are other normalisation factors required?

johndan updated this revision to Diff 169.Tue, Jan 5, 3:22 PM
  • D*K system should be in P-wave, not just D0K !
  • Changes after round of code-review with @jback
  • Remove whitespace and make minor bugfixes found during compilation tests
jback requested changes to this revision.Wed, Jan 6, 4:09 PM

Thanks @johndan, the changes you have made so far look good and makes the code more general. I have replied to your questions, and have also suggested that we calculate the barrier radii-squared terms at initialisation. Not sure what we can do for avoiding hard-coding the D* mass at initialisation (regarding access to LauResonanceMaker info).

src/LauKMatrixPropagator.cc
874

Yes, OK, this will avoid unnecessary calculations for spin zero.

906

We could store "gamAInvRadSq_[iCh] = a_/(radii_[iCh]*radii_[iCh])" at initialisation to avoid repeating these calculations every time calcGamma is called.

913

Fixed during the fit, but can be initialised to a different constant value for systematic checks.

932

Let us use "3cos^2 - 1".

This revision now requires changes to proceed.Wed, Jan 6, 4:09 PM
In D41#1074, @johndan wrote:

General question, @tlatham : should I cherry-pick commits to bring this differential in line with the master (e.g. c38fa2a34f4a772754d1e8f610def10747d95104 - LauFormulaPar fix)?

I think you should be able to rebase the local branch (where you're working on these changes) on master. But I also don't think it's strictly necessary since I think this will be done when "landing" the revision anyway.

What might be more crucial is when D39 gets landed to then do a rebase on master since there are changes in that revision that are also here IIRC. But let's cross that bridge when we come to it.

tlatham added inline comments.Wed, Jan 6, 4:55 PM
src/LauKMatrixPropagator.cc
78–79

You can use something like:

LauResonanceMaker::get().getResInfo("D*0")->getMass()->value()

Not very pretty but it should work I think.

johndan marked 9 inline comments as done.Thu, Jan 7, 1:04 PM

Thanks @jback and @tlatham, I've handled all the remaining comments. There are some deeper issues to resolve (see the Laura++ mattermost channel) before we land this. D39 needs to come first, anyway, and that one is ready to go.

johndan updated this revision to Diff 174.Thu, Jan 7, 1:07 PM
  • Correct displaced doxygen comment
  • Simplify setting of denominator parameter in barrier factor from file
  • Cache the a/R^2 calculation - thanks John!
  • Neater access to D*0 mass information from LauResonanceMaker - thanks Tom!
  • Add definition of angular momentum term for spin-2 (max spin) to K-matrix
johndan updated this revision to Diff 175.Thu, Jan 7, 1:07 PM
  • Correct displaced doxygen comment
  • Simplify setting of denominator parameter in barrier factor from file
  • Cache the a/R^2 calculation - thanks John!
  • Neater access to D*0 mass information from LauResonanceMaker - thanks Tom!
  • Add definition of angular momentum term for spin-2 (max spin) to K-matrix
johndan updated this revision to Diff 178.Fri, Jan 8, 9:04 PM
  • Remove done TODO comment and make const several function arguments
  • Add function to LauKMatrixPropagator that returns K-matrix spin, and make various function arguments const
  • Add barrier factors to numerator of K-matrix amplitude in ProdPole and ProdSVP, which need to know the spin of the K-matrix using the new getL() function
  • Fix typo in ProdPole.cc and ProdSVP.cc
jback accepted this revision.Mon, Jan 11, 2:22 PM
This comment was removed by jback.
This revision is now accepted and ready to land.Mon, Jan 11, 2:22 PM
jback requested changes to this revision.Mon, Jan 11, 6:55 PM

From further discussion elsewhere, it looks like the code needs further edits to allow channel-dependent L values (for a given total angular momentum J).

This revision now requires changes to proceed.Mon, Jan 11, 6:55 PM
johndan updated this revision to Diff 179.Mon, Jan 11, 8:13 PM
  • Changes to allow different orbital angular momenta for different channels inside the K-matrix
johndan updated this revision to Diff 180.Mon, Jan 11, 8:13 PM
  • Changes to allow different orbital angular momenta for different channels inside the K-matrix

I've included the changes to tell the K-matrix about the orbital angular momentum of each channel, as discussed on Mattermost. I've also made 3 comments in the code for you, @jback. Thanks!

src/LauKMatrixProdPole.cc
87

@jback, here I assume the 0th channel is the one for which the amplitude is requested. This assumption has precedence in your code, right?

src/LauKMatrixPropagator.cc
51

@jback I do not include the J^P for the K-matrix, I simply removed L. I think it's reasonable to think that the user has to be trusted to get that right as it's equally important that they check the validity of the coupled channels in this regard and we don't have the information about the J^P of the coupled channel particles that we would need in order to check for them.

98

@jback I was thinking about this in passing: this will give a L=2 barrier-factor-squared denominator of the form q^4 + 6q^2/R^2 + 9/R^4 which is, I think, not what you should have here. I think we are wanting the second term to be 3q^2R^2 aren't we?

jback added inline comments.Tue, Jan 12, 5:40 PM
src/LauKMatrixProdPole.cc
87

Well, you could have the case where the first channel has L = 0, but others have L > 0. So we should check if the propagator has any non-zero L channels, perhaps by returning a cached boolean that keeps track of this condition at initialisation.

Do we need to include the resBWFactor? I thought this was already in the "D" matrices (denominator terms involving q, a and R).

src/LauKMatrixPropagator.cc
51

OK, sounds reasonable.

98

No, the denominator factor uses L/2 for the power not L, so we get q^2 + 3/R^2 as expected.

johndan updated this revision to Diff 184.Thu, Jan 14, 8:05 PM
  • Correctly initialise a_ vector. Remove calculation of angular factor from K-matrix propagator. Remove any need to update the propagator with the LauKinematics object
  • Remove angular parts of K-matrix header, and add extra function to indicate the final-state channel of the various coupled channels
  • Update ProdPole and SVP to use resAmp (and leave amplitude to LauAbsResonance)
  • Update LauAbsResonance to apply spin to constructed objects that don't use the 'normal' constructor
  • Changes to LauResonanceMaker to allow retrieval of parent Blatt-Weisskopf factor without using resInfo that is not defined for K-matrix
  • Changes to LauBlattWeisskopfFactor to accommodate tangential use by K-matrix to get parent factor without access to a resInfo
johndan updated this revision to Diff 185.Thu, Jan 14, 8:05 PM
  • Correctly initialise a_ vector. Remove calculation of angular factor from K-matrix propagator. Remove any need to update the propagator with the LauKinematics object
  • Remove angular parts of K-matrix header, and add extra function to indicate the final-state channel of the various coupled channels
  • Update ProdPole and SVP to use resAmp (and leave amplitude to LauAbsResonance)
  • Update LauAbsResonance to apply spin to constructed objects that don't use the 'normal' constructor
  • Changes to LauResonanceMaker to allow retrieval of parent Blatt-Weisskopf factor without using resInfo that is not defined for K-matrix
  • Changes to LauBlattWeisskopfFactor to accommodate tangential use by K-matrix to get parent factor without access to a resInfo
johndan updated this revision to Diff 189.Fri, Jan 15, 11:21 AM
  • Add K-matrix name to gCoupling parameters to avoid clashes in models with >1 K-matrix
  • Resolve merge conflict
  • Allow other parameters of the SVP part to be floated
  • Fix merge conflict
johndan marked 6 inline comments as done.Fri, Jan 15, 11:46 AM

Thanks, @jback. A few answers to your responses, with just one that doesn't (yet) make sense to me...

src/LauKMatrixProdPole.cc
87

On the first point, I don't think it matters if other channels have L>0 does it? Why do we need to track it?

On the second, I think you're spot on - the resBWFactor has been dropped as you say.

src/LauKMatrixPropagator.cc
98

Hmmm not sure. If a = 3 then calcGamma gives:
pow(q,L) / pow( q*q + gamAInvRadSq_[iCh] , L_[iCh]/2. );
where gamAInvRadSq_ = a/(R*R)

i.e. the return is:

L = 2 : sqrt[ (qR)^2L / ( (qR)^2 + 3 )^L ]
L = 2 : sqrt[ (qR)^4 / ( (qR)^2 + 3 )^2 ]
L = 2 : sqrt[ (qR)^4 / ( (qR)^4 + 6 (qR)^2 + 9 ) ]

But looking at equation 98 in the Laura++ paper I would expect the coefficient of (qR)^2 in the denominator to be 3 not 6. What am I missing?

johndan updated this revision to Diff 190.Fri, Jan 15, 3:04 PM
johndan marked an inline comment as done.
  • Remove duplicate par_mSq0_
jback added inline comments.Fri, Jan 15, 6:00 PM
src/LauKMatrixProdPole.cc
87

Yes, we don't need to track this since we have changed how the angular dependence is included.

src/LauKMatrixPropagator.cc
98

You are not missing anything; it's just that we can't exactly reproduce the z^4 + 3z^2 + 9 factor using a perfect square product between two D terms without introducing complex terms. Remember, this barrier term is an ad-hoc "best guess" at what the momentum dependence should be. The important thing is that we have a quadratic dependence on z^2, and we can change "a" and "R" at initialisation for systematic variations.

johndan marked 5 inline comments as done.Thu, Jan 21, 12:41 PM
johndan added inline comments.
src/LauKMatrixPropagator.cc
98

OK! Agreed, the dependence is what matters.

johndan updated this revision to Diff 197.Thu, Jan 21, 12:54 PM
johndan marked an inline comment as done.
  • Use default values at declaration in LauAbsResonance to significantly simplify the constructors