New Insights on Low Energy $\pi N$ Scattering Amplitudes: Comprehensive Analyses at $\mathcal{O}(p^3)$ Level

A production representation of partial-wave $S$ matrix is utilized to construct low-energy elastic pion-nucleon scattering amplitudes from cuts and poles on complex Riemann sheets. Among them, the contribution of left-hand cuts is estimated using the $\mathcal{O}(p^3)$ results obtained in covariant baryon chiral perturbation theory within the extended-on-nass-shell scheme. By fitting to data on partial-wave phase shifts, it is indicated that the existences of hidden poles in $S_{11}$ and $P_{11}$ channels, as conjectured in our previous paper~\citep{Wang:2017agd}, are firmly established. Specifically, the pole mass of the $S_{11}$ hidden resonance is determined to be $(895\pm81)-(164\pm23)i$ MeV, whereas, the virtual pole in the $P_{11}$ channel locates at $(966\pm18)$ MeV. It is found that analyses at the $\mathcal{O}(p^3)$ level improves significantly the fit quality, comparing with the previous $\mathcal{O}(p^2)$ one. Quantitative studies with cautious physical discussions are also conducted for the other $S$- and $P$-wave channels.


Introduction
Pion-nucleon (πN ) scattering, as one of the fundamental processes of low-energy quantum chromodynamics (QCD), has been extensively studied since the middle of the last century, see e.g. Ref. [1]. Up to now, there exists a wealth of long-time accumulated experimental data on differential cross section and polarizations [2]. To decode the inherent physics, efforts have been devoted to performing partial wave analysis (PWA) of the πN scattering amplitudes by groups: Karlsruhe [3,4], Matsinos [5], GWU [6], etc. Pole information about excitations, e.g. the Roper resonance N * (1440) [7], can be extracted from PWA, usually relying on various phenomenological models used in the analyses. However, a modern study of πN scatterings is no longer restricted to the extraction of the nucleon excitations, but aims at a precision description of the πN amplitudes, not only in the physical region but also in the subthreshold region, and offers reliable inputs for any other related physics, for instance, πN sigma term serving as a crucial input for the interpretation of dark-matter searches [8,9]. As such, it is necessary to perform investigations in a model-independent way so as to gain an overall understanding of the underlying dynamics of πN physics.
At low energies, baryon chiral perturbation theory (BChPT) is one of such model-independent methods in the study of πN scatterings. As an expansion in powers of external momenta and light-quark masses, BChPT relies on a hierarchy of the contributions (Feynman diagrams) known as power counting. The presence of baryon fields as explicit degrees of freedom introduces new scales that do not vanish in the chiral limit, leading to the occurrence of power counting breaking (PCB) terms [10] in loop diagrams. In the past thirty years, to remedy this issue various approaches have been proposed: e.g., the heavy baryon (HB) formalism [11,12], the infrared regularization (IR) prescription [13,14], and the extended-on-mass-shell (EOMS) scheme [15][16][17] (see also Ref. [18] for extension beyond low-energy region). Consequently, analyses on low-energy as well as resonance region of πN scatterings have been done in the various approaches, see, for example, Refs. [19][20][21][22][23][24][25][26][27][28][29]. In particular, results obtained in the EOMS scheme has proven to be more suitable for amplitude analysis in the sense that, compared to the other schemes, the analytic structure of the calculated physical quantities is properly kept [30][31][32]. Though great achievements have been gained using BChPT with EOMS scheme, one should still keep in mind that those results only fulfill the unitarity perturbatively. To restore exact unitarity, a plenty of unitarization methods have been developed, but with shortcomings like violation of analyticity, breaking of crossing symmetry, etc.. Thus, more rigorous model-independent tools are required.
A rigorous manner is to construct partial-wave πN amplitudes under the guidance of the axiomatic principles of S-matrix theory, such as unitarity, crossing symmetry, analyticity and so on. Attempts in this direction have been performed in Refs [33][34][35][36] by imposing dispersion relations. Likewise, the recent Roy-Steiner (RS) equation analysis in Refs [37,38] is an alternative and has gained great achievements in the description of πN physics. More recently, in Ref. [39], a fresh look at low energy πN scatterings has been established with the help of the Peking University (PKU) representation [40][41][42][43], a production parametrization of partial wave S-matrix on the whole complex plane for two-body elastic scattering amplitudes. In this paper, we intend to carry out a more comprehensive study of the partial wave πN scattering amplitudes using PKU representation.
The advantage of the use of the PKU representation is two-fold. On the one hand, it is suitable for pole analysis. The PKU representation separates partial waves into various terms contributing either from poles or branch cuts. The corresponding phase shifts extracted from PKU representation are sensitive to subthreshold poles, enabling one to determine pole positions rather accurately. Furthermore, each phase shift contribution has a definite sign, which makes possible the disentanglement of hidden poles from a background. On the other hand, it respects causality honestly. In the PKU representation, the pole contributions are regarded either as hidden poles or as known poles fixed by experiments, while the cut contribution is estimated from perturbative BChPT amplitudes and uncertainties from such an estimation is known to be severely suppressed. Importantly, the philosophy of the PKU representation is not to directly unitarize the amplitude itself, instead, it unitarizes the left-hand (and inelastic) cuts of the perturbative amplitude and hence hazardous spurious poles, violating causality, can be avoided [44]. In addition, the consistency of that representation with crossing symmetry is already examined in Refs. [45,46]. In Ref. [39], the application of the PKU representation to pion-nucleon scatterings has led to very interesting findings: a resonance below πN threshold (denoted as crazy resonance for short) in the S 11 channel is discovered for the first time, and a companionate virtual state of the nucleon in the P 11 channel is found. Those findings are drawn with the left-hand cuts (l.h.c.s) estimated by tree-level BChPT amplitudes. That is, only the kinematical cut (−∞, (M N − m π ) 2 ] and the segment cut [(M 2 N −m 2 π ) 2 /M 2 N , M 2 N +2m 2 π ] due to the u-channel nucleon exchange are taken into account, where M N and m π are physical masses of the nucleon and the pion, respectively. In fact, as demonstrated in Refs [47,48], the full structure of the l.h.c.s for πN scattering contains a circular cut as well. Thus the motivation of this work boils down to two aspects as follows. Firstly, we intend to consider the whole structure of the l.h.c.s and verify that the S 11 hidden pole still exists in O(p 3 ) analysis. Meanwhile, the O(p 3 ) calculation also provides cross-check of the O(p 2 ) result in the P 11 channel. Secondly, since the O(p 3 ) calculation is expected to determine the contributions from left-hand cuts more precisely, we try to scrutinize the minor discrepancies in the other channels, which were not able to be fitted well in O(p 2 ) case. Note that the one-loop BChPT amplitudes are renormalized within the EOMS scheme mentioned above. The utility of EOMS-renormalized amplitudes guarantees the correct analyticity behaviour as required by the foundation of PKU representation, i.e., the principles of S-matrix theory.
In Sect. 2 basic formalisms of BChPT, partial wave projection and the PKU representation are introduced. Sect. 3.1 collects the numerical results of the known contributions used in the PKU-representation analysis at O(p 3 ) level, then in Sect. 3.2 the two hidden states, pointed out in Ref. [39], are examined and, meanwhile, the other four Sand P -wave channels are also studied. The scattering length and effective range parameters for all the Sand P -waves are estimated in Sect. 3.4. Finally, conclusions and outlook are made in Sect. 4. Appendices. A.1 and A.2 exhibit the calculations of the chiral amplitudes for tree diagrams and one-loop diagrams respectively, while Appendix. A.3 shows the procedure of renormalization and EOMS subtraction. In Appendix. B the uncertainties of the contributions from l.h.c.s are discussed.
2 Theoretical framework 2.1 Formal aspects of πN scattering amplitudes The isospin structure of the πN amplitude can be decomposed as where χ i and χ f correspond to the isospinors of initial and final nucleon states, respectively. The amplitudes with isospins I = 1 2 , 3 2 can be obtained by As for the Lorentz structure, for an isospin index I ∈ { 1 2 , 3 2 } or I ∈ {+, −}, where s = W 2 ≡ (p + q) 2 , t ≡ (p − p) 2 are Mandelstam variables, and q (p) and q (p ) are the 4-momenta of initial and final state pions (nucleons), respectively. A and B are scalar functions of s and t.
One can substitute the nucleon spinors u(p, s) andū(p , s ) in Eq. (4) by helicity eigenstates in the centre of mass frame to obtain the following helicity amplitudes the first and second subscripts refer to the helicities of the initial and final nucleon respectively, and subscripts "±" are shorthands for helicity ±1/2; z s = cos θ with θ the scattering angle. The partial wave amplitudes for total angular momentum J can be written as where d J is the standard Wigner d-function. Finally one can get the six Sand P -wave amplitudes (in L 2I 2J convention) as follows: Throughout this work, the partial-wave label L 2I,2J is always suppressed if no confusion is caused.
The branch-cut structure of the partial-wave pion-nucleon scattering amplitudes is generally discussed in Refs. [47,48] with details, which can be shown schematically in Fig. 1. It is worth stressing that the circular cut stems from the t channel continuum [4m 2 π , 4M 2 N ], which is absent at O(p 2 ) level but starts to appear at O(p 3 ) level. Figure 1: The branch cuts of the partial wave πN elastic scattering S matrix on the s plane (schematic drawing). The cuts are represented by thick lines and circle. The right hand cut corresponds to the line [s R , +∞), and all the others are left hand cuts. The following abbreviations are used: The gray disk indicates the BChPT valid region characterized by the convergence radius r χ . The intersection of the BChPT convergence circle and the on-axis left-hand cuts is denoted by s c .

Perturbative BChPT description
The A and B scalar functions in Eq. (4) can be calculated perturbatively by invoking BChPT. In the framework of covariant BChPT of SU (2) case, the nucleon field is notated as N , and the pion as u-field where F is the pion decay constant in the chiral limit and τ stands for Pauli matrices in the flavour space. Relevant pieces of chiral effective Lagrangians for a calculation up to the leading one-loop order, i.e., O(p 3 ) level, are shown as follows [49]: Here M and g are the mass of the nucleon and the axial current coupling constant in the chiral limit, while c i s and d j s are low-energy coupling constants (LECs); · · · means matrix tracing in the flavour space. Furthermore, the relevant building blocks are [49] D µ = ∂ µ + Γ µ , with l µ , r µ , s, p being external sources and B 0 a constant related to quark condensation. Finally, the interaction vertex of four pions is also relevant to the O(p 3 ) calculation, and the following pion-pion interaction Lagrangian [50] in pure meson ChPT is needed With the above Lagrangians, it is readily to derive the πN scattering amplitudes up to O(p 3 ), which are necessary for estimating the contribution of the l.h.c.s in the PKU representation to be discussed in the next subsection. Though the leading one-loop BChPT amplitudes can be found in Refs [26,32], for the sake of easy reference, here we demonstrate in Appendix A the details of the leading one-loop calculation within the framework of BChPT using EOMS scheme.
The tree and loop amplitudes are collected in Appendices. A.1 and A.2 respectively. As pointed out in Ref. [10], the loop amplitudes contain both ultraviolet (UV) divergences and PCB terms, i.e., the terms in lower chiral order than the naively expected order. For the UV divergence, the dimensional regularization method is employed to pick out the relevant pieces (see Appendix. A.2). We use on-shell renormalization scheme to handle the self-energy of the nucleon, which maintains the correct analytical behavior at the point s = M 2 N , see Appendix. A.3.1. On the other hand, to remedy the power counting breaking problem, we use EOMS scheme to separate the PCB polynomials, which are then absorbed by the LECs appearing in the chiral effective Lagrangians. Those procedures are shown in Appendix. A.3.2.
The analytical expressions for all the involved loop integrals are compiled in Appendix A.4. Those expressions are obtained by means of dispersion relations with the spectral functions calculated using Cutkosky rule. Therefore, in principle, the BChPT amplitudes shown in this paper are calculable on the whole complex Mandelstam plane.

Unitary PKU representation
The PKU representation of the partial wave πN elastic scattering S matrix can be written as where S b , S v and S r represent the individual contributions of bound states, virtual states and resonances respectively (see Ref. [39] for their explicit expressions). The exponential term amounts to the background that carries the information of l.h.c.s and right-hand inelastic cut (r.h.i.c.) above inelastic thresholds. The kinematic factor ρ(s) is given by with s L = (M N − m π ) 2 and s R = (M N + m π ) 2 . The function f (s) satisfies a dispersion relation, where L and R abbreviate l.h.c.s and r.h.i.c., respectively. The discontinuities of f (s) along the various cuts can be deduced from the partial-wave BChPT amplitudes T (L 2I,2J ) [c.f. Eq. (7)] through the following relations: discf (s ) = disc ln S(s ) 2iρ(s ) , Specifically, the partial-wave πN elastic scattering amplitudes T (s ) are obtained with Eq. (7) by inserting the BChPT amplitudes calculated up to O(p 3 ) level in the previous subsection. Taking into account the realistic l.h.c.s and r.h.i.c., as depicted in Figure 1, Eq. (16) is expanded to the form where S in and S out are the S matrices along the circular cut from inside and outside, respectively, and 0 < η < 1 is the inelasticity along the r.h.i.c.. Furthermore, three cut-off parameters, s c , θ c and Λ R , have been taken for the kinematical cut (−∞, (M N − m π ) 2 ], the circular cut and the r.h.i.c. [(2m π + M N ) 2 , ∞), in order. The first term of Eq. (18) is negative definite and dominates the contribution to f (s). The other three terms are quite small numerically. The second term originates from the u-channel nucleon exchange which can be approximated as contact interaction at low energies. The third piece only contributes at O(p 3 ) level, while the fourth term is positive definite but has been shown to be of insignificant impact [39]. Eventually, the resulting function f (s) provides a definitely negative contribution to the partial-wave phase shift, given that the cut-off parameters are chosen inside reasonable domains.
The fact that f (s) contributes negatively is of crucial importance in exposing hidden poles on the second Riemann sheet. Apart from f (s), the other source of negative contribution is thoroughly owing to bound states, which are clearly seen experimentally and hence their contribution can be fixed a priori. As a result, the disentanglement of the hidden poles with positive contributions from a negative definite background can be easily achieved, which will be explored in the next section.
3 Partial wave phase shifts and pole analyses 3

.1 Pole-hunting strategy
The production presentation of the S matrix results in additive phase shifts with definite signs that are implied by the contributing sources. In general, a bound state contributes negatively, while a pole (virtual or resonance state) on the second Riemann sheets gives positive contribution. The sum of the cut contributions, identical to ρ(s)f (s), is negative. Therefore, it is rather convenient to carry out exploration of hidden poles by using partial wave phase shifts, denoted by δ(s) henceforth.
For easy explanation of the pole-hunting strategy, the S-matrix in Eq. (14) can be rewritten as In above, the term S h (s) accounts for the contribution from all possible hidden poles, e.g. virtual states, resonances below the threshold or with extremely large widths, or shadow poles. They cannot be observed directly by experiments, but their existence (and pole positions) can be examined through PKU representation once the last two terms in the first equality are known.
To that end, in what follows, we will specify the known-pole and cut contributions by discussing the necessary inputs. The known-pole contribution can be fixed with their corresponding experimental information. The known poles contain bound states and above-threshold resonances, both of which are well determined experimentally. In our case, the Sand P -wave known poles below 2 GeV are complied in Table 1. The bound state in the P 11 channel is the nucleon itself. For the known resonances, some of them are above the inelastic ππN threshold and hence are actually located on the third and higher Riemann sheets. Nevertheless, under narrow width approximation, their corresponding shadow poles on the second Riemann sheet may be estimated by which contribute to S matrix in the same way as the standard resonances on 2nd Riemann sheet, like ∆(1232). Table 1: Known poles in the two S-wave and four P -wave channels.
The cut contribution can be estimated from BChPT amplitudes as demonstrated in subsection 2.3. There are three free cut-off parameters, i.e. s c , θ c and Λ R in Eq. (18). It is emphasized that the cut-off parameters for left-hand cuts, i.e. s c and θ c , are in principle unknown parameters. However, to get a first glimpse of the physics, one can assign the values of s c and θ c in accordance with the validity region of BChPT. Similar to the previous work, we assume that the results of l.h.c.s from the perturbative BChPT calculation remain correct till meeting the N * (1440) shadow pole, above which complicated coupled-channel dynamics takes place. The distance between the shadow pole and the center point of chiral expansion is denoted as r χ . The boundary of BChPT valid region should be a circle centering at s χ = M 2 N + m 2 π with the radius r χ . The intersections of this circle with the l.h.c.s give s c = −0.08 GeV 2 and θ c = 1.18 radian, see Figure 1. The two parameters are actually related by throughout this work: θ c will always be obtained from s c , which sometimes takes values other than the chiral estimation −0.08 GeV 2 , by using the above relation. As for Λ R , following Ref. [39], it is set to Λ R = 4 GeV. In addition, the experimental data of inelasticity η taken from [6] is used in the calculation of the r.h.i.c. contribution. In above we discussed the validity region of BChPT. However, in no way one should expect that the contribution beyond the validity region of BChPT should be negligible. The two parameters s c and θ c are actually taken as free parameters in the fit. Their deviation from the boundary of validity region signals the contribution from high energy region 1 .
Furthermore, to compute the cut contribution, the masses, the nucleon axial charge and the pion decay constant in the BChPT amplitudes take the following values: while the relevant LECs, c i and d j in Eqs. (10) and (11), are set to the fit values of Fit I in Table. I of Ref. [27] 2 . Other recent determinations, e.g. Refs. [28,51], yield very similar LECs values, and it has already been illuminated in Ref. [39] that the change of the LECs has inconsiderable impact on the l.h.c.s. Moreover, the values of the LECs are determined in Ref. [27] by performing fit to the phase shifts generated by the recent RS equation analysis of the πN scatterings in Ref. [38], where both central values and error bars are given. To be consistent, for the pole analyses we will employ the recent RS partial wave phase shifts as well, rather than the ones from PWA by GWU group without errors [6] which are used in Ref. [39].
With the above preparations, we are now in the position to perform phase shift analyses, and explore possible hidden poles with the equipment of the PKU representation.

Phase shift analyses
In Ref. [39], among the six Sand P -wave channels, the S 11 and P 11 waves have been paid special attention to in the study at O(p 2 ) level. It is important and necessary to verify if the hidden poles, discovered therein for these two channels, still keep their presence in an O(p 3 ) analysis. The O(p 3 ) analysis will take into account the full structure of l.h.c.s, as mentioned in the Introduction, and hence is adequate to verify the very existence of those poles. Besides, higher-order contributions in the chiral expansion will not introduce any other new branch-cut structures, and they should have little impact on the O(p 3 ) results obtained here. One refers to Appendix. B to see the uncertainties caused by numerical details such as the cut-off parameters and the unitarization methods.
In what follows, it will be shown that the O(p 3 ) results indeed provide further and stronger evidences indicating that the hidden poles in S 11 and P 11 channels definitely exist. In particular, the extra hidden poles in S 11 and P 11 waves do play an essential role in establishing a meaningful numerical fit -without them the chi-squares would be unacceptably large. On the other hand, we will also carry out careful investigations for the other Sand P -wave channels, which were less studied in Ref. [39]. It can be observed that the S 11 channel suffers from a large discrepancy between the RS data and the sum of all the known contributions. Unfortunately, this disagreement cannot be rescued by fine-tuning the involved parameters in the numerical computation. Instead, a hidden pole contribution, as a new dynamical structure, must be included so as to compensate such discrepancy, as suggested in the O(p 2 ) study carried out in Ref. [39]. In the right panel of Fig. 2, the known contributions are displayed order by order. Compared to the O(p 2 ) case, our thorough O(p 3 ) calculation now enlarges the discrepancy 3 , which strongly indicates the necessity of the hidden pole. Furthermore, the fit result is unique, at least for two poles. For example, if we add one resonance together with one virtual state, the virtual state will move automatically to the pseudo-threshold and the remaining resonance results in the same location as the one resonance scene. In addition, if only two virtual states are taken into consideration, they will automatically collide on the real axis and turn into one resonance. This property of the stability of pole positions is an imperceptible advantage of the present method.  To proceed, the position of the hidden pole can be determined by a fit to the RS phase shift data. The fit results in a crazy resonance below threshold, which is quite stable against the variation of the cut-off parameters. The hidden pole positions with various cut-off parameters are listed in Table. 2 4 . The fit results are plotted in Fig. 3 and the RS data are well described due to the inclusion of the hidden pole 5 .
The pole position of this crazy resonance is given by The numbers in the first brackets are systematical errors responsible for the variation of the cutoff parameter s c (see Table 2), while the ones in the second brackets are the averaged statistical uncertainties from fitting 6 . By adding the two uncertainties in quadrature, our final reported result is It is compatible with the determination reported by the O(p 2 ) analyses [39], i.e. 861(53) − 130(75)i MeV 7 .

P 11 channel
Likewise in the P 11 channel, the RS data is far from being described by the known contributions as shown in the left panel of Fig. 4. Furthermore, it can be seen from the right panel of Fig. 4 that the discrepancy in P 11 channel becomes larger when the chiral order increases. Hence, the same prescription proposed in Ref. [39] should be employed to resolve the inconsistence problem in the P 11 channel.
Following Ref. [39], two virtual states are added to the S matrix, and one stays in the near-threshold region as the major hidden contribution, while the other is nearly absorbed by the pseudo-threshold (M N − m π ) 2 . 8 The locations of the near-threshold virtual state with different cut-off parameters are listed in Table. 3, and the corresponding phase shifts are shown in Fig. 5. Note that a P -wave threshold constraint, i.e. δ(k) = O(k 3 ) for 3-momentum k, has been considered in the fits. It should be pointed out that no satisfied fit can be achieved with s c = −0.08 GeV 2 . Nevertheless, if |s c | becomes a bit larger, e.g. s c = −1.00 GeV 2 , the fit quality is significantly improved. However, this is not true for the O(p 2 ) fit in Ref. [39]. Therein a rather large |s c | like s c = −9.00 GeV 2 has to be adopted. This improvement should be owing to the full consideration of the various l.h.c.s, such as the circular cut. As a consequence, the low-energy fit is now no longer sensitive to the physics in high energy region. Finally, based on the results in Table. 3, we report the location for the extra near-threshold virtual pole: Note that similar to what is discussed in Ref. [39], the O(p 3 ) perturbative S matrix gives a first-sheet zero at √ s = 973 MeV, which is compatible with the fit result Eq. (25). In principle, one may also perform fit with the cut-off parameter s c released, meanwhile, θ c is determined by Eq. 21. However, due to the existence of many local-minimums, such fit does not converge and the resultant statistical errors are untrustworthy any more.
The physical mechanism of the P 11 virtual state is already discussed in Ref. [39]: it is a partner of the nucleon pole, which is different from the S 11 hidden pole, i.e. (likely) from potential mechanism. Case in S 31 channel is quite interesting. As the left panel of Fig. 6 shows, firstly one finds some disagreements between the known contributions and the data, but those disagreements are minor: at least the known terms give the correct sign of the phase shift, and it seems that fine-tuning of the cut-off parameter s c can make the curve match the first two data points near threshold, which is totally different to the cases in S 11 and P 11 channels. Unexpectedly, situation becomes awful once s c is released as a fitting parameter. The RS data cannot be fitted well due to the mismatching of line shape -the data curve is too "straight" to be described by the current formula. Besides, this line-shape mismatching can not be cured by means of truncating the chiral series, since O(p i ) (i = 1, 2, 3) calculation all give similar line shapes, see the right panel of Fig. 6.
This observation may indicate that the evaluation of l.h.c.s is not accurate enough, which has various possible origins: some additional minor cuts like the left-hand cuts from tor u-channel resonance exchange, improper shadow pole positions of the known resonances due to the narrow width approximation (see Eq. (20)), uncertainties from LECs, or some quite distant hidden contributions. It is impossible to take all the factors into consideration and do an all-embracing analyses.
However if one insists on getting a good fit to the data, the easiest way is to employ a pole term as an extra background, similar to what have been done in the P 11 channel. Here one virtual state is put into the S matrix, and the results are shown in Table. 4. A good fit needs a larger |s c | value, e.g. when |s c | > 9 GeV 2 one gets χ 2 /d.o.f < 3. At the same time, the line shape becomes better when |s c | is large enough, see Fig. 7. It should be emphasized here the O(p 2 ) fit fails to give a good chi-square when |s c | ≤ 25 GeV 2 .   It is interesting to find that the location of the virtual state serving as background is not so far away from the threshold (although it lies on the left of the u-channel nucleon cut), and seems stable, which implies that the extra corrections of l.h.c.s in this channel is rather important.

P 31 channel
Similar to the S 31 channel, there exist some minor discrepancy between the known contributions and the data in P 31 channel, see the left panel of Fig. 8. The discrepancy remains at each chiral order, as can be seen from the right panel of Fig. 8.
The s c parameter can be tuned to be 0.14 GeV 2 to satisfy the P -wave threshold constraint, i.e. δ(k) = O(k 3 ), but the tuning of s c does not work since it gives a disastrous fit except at threshold. This drives us to add extra poles again to simulate the missing corrections of the l.h.c.s. The fit with extra poles are shown in Table. 5, and plotted in Fig. 9.
Only when s c is around −9 GeV 2 a good fit can be achieved. Furthermore, if one uses the l.h.c.s from O(p 2 ) level and perform the same fit, |s c | should be rather large e.g. |s c | > 16 GeV 2 , to get a result with the quality as good as s c = −1 GeV 2 here.

P 13 channel
The situation in P 13 channel is similar to that in P 31 channel: the left panel of Fig. 10 is the comparison between known terms and the data, and the right panel is the known contributions at different chiral orders.
The fit results with an extra resonance as background are shown in Table. 6 and Fig. 11, and the fit quality is worse than the P 31 channel. However, O(p 2 ) estimations can never give a result with χ 2 /d.o.f < 10 for any s c value.    Figure 11: The fit with one extra resonance as background in P 13 channel.

P 33 channel
The discrepancy between known contributions and data in P 33 channel is the smallest among the six channels. The comparison between the known terms and the data is shown in the left panel of Fig. 12. The higher the chiral order of the BChPT amplitudes used in the estimation of l.h.c.s is, the better the description of the phase shift becomes, see the right panel of Fig. 12.
In this channel the P -wave threshold constraint can be satisfied by tuning s c to be −0.18 GeV 2 . This is a parameter-free calculation. It is amazing that the data are automatically described with a quality χ 2 /d.o.f = 3.39, see Fig. 13. One can also use a different strategy by varying the pole position of ∆(1232): its location is only shifted slightly (from 1.211 − 0.050i GeV to 1.208 − 0.047i GeV), and the χ 2 /d.o.f is reduced to 1.78. 3.3 Further remarks on the S 31 , P 31 and P 13 channels From the above discussions, it is found that, to obtain good description of the data in S 31 , P 31 and P 13 channels, one still needs some extra terms, which are for the moment simulated by single pole terms. Nevertheless, one should be cautious to regard them as physically existing poles. The following criterion might be helpful to distinguish a real pole from a fake one.
• Large discrepancies between the known contributions and the data. As already mentioned in previous O(p 2 ) analyses, the discrepancies in S 11 and P 11 channels are crucial since even if the left-hand cut terms are switched off, the discrepancies are still significant, which forces us to add new poles. On the contrary, the discrepancies in other channels are minor, because at least one can fine-tune the cut-off parameter to make the scattering lengths from the known terms conform with the data. To be conservative, we prefer the missing contributions in those channels to be corrections of the background term a priori.
• Stability of the pole positions. It is obvious that the S 11 hidden pole stays stable and fits the data well for all the cut-off values we chose. On the contrary , the poles in S 31 , P 31 and P 13 channels only fit the data when the cut-off parameters are large, and some of them are not so stable (e.g. in P 31 channel the poles change from virtual states to resonances).
• Alternative physical mechanism in support of the poles. As pointed out in previous O(p 2 ) analyses, the P 11 pole can be regarded as the shadow pole of the nucleon, and the S 11 pole may be of potential dynamical nature. However we fail, at least at the moment, to find a possible physical interpretation of the poles in S 31 , P 31 and P 13 channels. Actually we have tried to use square-well potential as a toy model to fit the S 31 data but it failed: the potential would become repulsive, leaving no nearby poles, meanwhile the line-shape cannot be well fitted.
The poles in S 31 , P 31 and P 13 channels do not meet all the three criterion raised above, hence we can not claim those extra poles are truly existing in reality. It is also worth noticing that the extra contributions in S 31 , P 31 and P 13 channels are indeed much smaller than the S 11 and P 11 hidden poles. For S 11 , P 31 and P 13 channels, the comparison of the phase shifts contributing from the extra resonances is shown in Fig. 14. For P 11 and S 31 channels where virtual states are added, the extra P 11 state contributes to the phase shift about two times as large as the S 31 one.

Threshold parameters
All the Sand P -wave phase shifts have been analyzed by means of the PKU representation with proper inputs taking from O(p 3 ) BChPT. For some channels, especially P 31 and P 13 waves, distant resonances have to be incorporated artificially as background to simulate any other underlying contributions. Nonetheless, our resultant predictions are in good agreement with the data in the low-energy region close to threshold, which therefore can reliably be applied to the extraction of threshold parameters, namely the coefficients of effective range expansion (ERE). As usual, the ERE is defined by where L is the angular momentum, k is modulus of the 3-momentum in the center of mass frame, and δ is the phase shift. The coefficient a L in Eq. 26 is usually called scattering length (or volume for P -waves), while r L is termed as effective range. For L = 0 they are both in dimension of length, characterizing the effective size of the target and the range of the interacting potential as well known in potential scattering theory, respectively. For a better extraction of threshold parameters, the fit results with χ 2 /d.o.f. > 5 are excluded. For the scattering lengths (volumes) a L and the parameters r L ,  Table 7: Sand P -wave scattering lengths (in the unit of 10 −3 m −1 π ) and volumes (in the unit of 10 −3 m −3 π ).
Channel r L from PKU analyses  Table 8: Sand P -wave effective ranges (in the unit of m −1 π for Swaves and m π for P -waves).
Furthermore, the non-relativistic limit of PKU representation (i.e., Ning Hu representation [52]) is where k is the 3-momentum in the center of mass frame, p n denotes the pole position on the complex k plane, and the exponential term is the background term with R being a constant parameter. According to Ref. [53], R may be positive definite and approximately equal to the potential range. Here one can also investigate the R parameters using the results from PKU representation. The values of R parameters in each Sand P -wave channels (calculated under the best fit of the cut-off parameter s c ) are exhibited in Table. 9. Note that in πN scatterings the potential can be regarded as the effect of two-pion exchanges in t channel, thus R should be Channel S 11 S 31 P 11 P 31 P 13 P 33 R parameter 0.56 0.47 0.52 0.38 0.38 0.30 Table 9: Sand P -wave R parameters from Eq. (27) (in the unit of m −1 π ).
roughly 1/(2m π ), hence the results in Table. 9 turn out to be reasonable. Similarly, application of the simplified version of PKU representation to nucleon-nucleon scatterings can be found in Ref. [54].

Conclusions and outlook
In summary, this paper, as a follow-up of Ref. [39], applies PKU representation to the analyses of elastic πN scattering processes in Sand P -waves, with the left-hand cut contributions calculated from the results obtained in covariant baryon chiral perturbation theory up to O(p 3 ) level. It is found generally that the fit quality has been greatly improved compared with the previous O(p 2 ) results. The existence of the hidden states in S 11 and P 11 channels are firmly established irrespective of the chiral order as well as other numerical details.
The other four channels are also investigated quantitatively: a good description of the data in P 33 wave can be achieved with a cut-off parameter determined by the threshold constraint, while in the other three waves some finer structures are needed that are temporarily simulated by extra poles. As byproducts, the parameters of effective range expansion for all the channels are calculated.
Although the analyses of low-energy pion-nucleon scatterings based on PKU representation with left-hand cuts simply from perturbative calculation are completed, the present work could still be improved in future. The major weakness of the present work, though not vital in drawing the current conclusions as we believe, comes from the uncertainty when calculating distant lefthand cuts. Hence through our work, we only cautiously draw physical conclusions extracted from discrepancies at qualitative level, and avoid any results drawn only depending on minor discrepancies between data and the cuts. It is worth stressing that our major conclusions are actually established based on the negative definiteness of the background contribution, which is supported by quantum mechanical scattering theory. The S 11 and P 11 poles turn out to be very stable and hence trustworthy. There have been successful examples for this approach: the f 0 (500) pole and the K * (800) pole positions in Refs. [42,43] are found in very good agreements with the results from Roy-like equation analyses(see Refs. [55,56]). One may improve the calculation by invoking resonance plus Reggeon exchange model to estimate the distant lefthand cuts. Constrains like crossing symmetry may also be incorporated. Finally, the physical meaning of the newly found hidden poles are still open for discussions.  • At O(p 1 ): • At O(p 2 ):

A.1 Tree-level
• At O(p 3 ) (Born diagram): • At O(p 3 ) (contact diagram): Note that Eq. (30) can be absorbed in Eq. (28) by redefinition of the g coupling. From the A, B functions above one can readily get the analytical expressions of the amplitudes in each channel after partial wave projection, as done in Appendix. B of Ref. [39].

A.2.1 Diagrams
According to the power counting rule of BChPT, for a loop diagram with L independent loops, I π internal pion lines, I N internal nucleon lines, and V n vertices of O(p n ), its chiral order Consequently, the diagrams, which are classified into three groups and displayed in Figs  To be specific, the first group (G 1 ) contains the diagrams of contact type, which do not have companions related by crossing symmetry; the second (G 2 ) and third groups (G 3 ) are Born diagrams, which have corresponding u-channel crossing diagrams but are not shown explicitly; moreover, the time-reversal counterparts of the G 2 diagrams are not displayed either. The sum of those diagrams is • Three point functions

A.2.3 Results of the one-loop amplitudes
For brevity, the abbreviations of the loop functions in Table. 10 will be used from now on. All the A, B functions corresponding to the diagrams in Figs. 18, 19 and 20 are shown diagram by diagram in the following.

A.3.1 UV cancellation
The UV divergent pieces in loop functions are dealt with in dimensional regularization and expressed in terms of R as defined in Eq. (77), see Appendix. A.4. The MS − 1 subtraction scheme is employed in this paper, i.e. the counter terms should completely cancel the divergent pieces proportional to R . However, the diagram (3n) in Fig. 20 needs to be treated specifically, since MS − 1 scheme itself can not erase the second order pole at s = M 2 in Eq. (51), which is a disaster of analyticity properties especially when investigating the u-channel nucleon cut. To remedy such a problem, one should regard the loop in the internal nucleon propagator in diagram (3n) as the renormalized self-energy Σ r ( / p), which satisfies the on-shell renormalization conditions Such conditions can be fulfilled by choosing suitable values of the nucleon wave function renormalization constant Z N and the mass counter term δM . To be specific, the bare parameters can be split as follows where the LECs with superscript "0" denote the bare quantities and the ones with "r" stand for renormalized quantities. The corresponding coefficients γ are and in addition, where Table. 10 for the abbreviations of the involved loop functions. Note that except for Z N and δM , the other values in Eq. (55) are the same as in Ref. [26,28]. Lastly, the renormalized result of diagram (3n) reads where 1. About C t 2. About C t2m 3. About C s • Four point functions Now only the expressions of those scalar loop functions are needed. One point and two point functions are simple when the dimension D is set to 4: where with γ E the Euler constant; and B 0 (p 2 ; m 2 , M 2 ) = 1 − R + ln In Eqs. 76 and 78 we set ν = M . However, three and four point functions do not have such compact analytical expressions. The most commonly used representation of those functions in textbooks is the Feynman representation, which expresses the loop functions by integrals with respect to Feynman parameters. However in our calculation it is found that the Mandelstam spectral representation is more preferable numerically. In the Mandelstam spectral representation the loop functions are calculated by means of dispersion integrals: At O(p 3 ) level there are two traditional schemes to unitarize the amplitude, i.e. K-matrix approach in Eq. (82), and Padé approximation in Eq. (83) where T i represents the perturbation amplitude at O(p i ) level in each channel; ReT 3 stands for the real part of O(p 3 ) amplitude in physical region, and generally should be regarded as T 3 −iρT 2 1 . The two methods both give manifestly unitary T , i.e. ImT = ρ|T | 2 . The contributions of l.h.c.s given by perturbation calculation and the two unitarization methods are shown in Fig. 22. It is found that the uncertainties in P 31 and P 13 channels may be larger since the two unitarization methods both give positive results in these two channels.
In the end we emphasize that, as claimed in the preceding paper [39], conventional unitarization methods may contain spurious poles that violate causality, thus they are not trustworthy in analyzing pole structures.