Noémie Debroux
email: noemie.debroux@uca.fr

Guillaume Lienemann
email: glienemann@chu-clermontferrand.fr

Benoît Magnin
email: bmagnin@chu-clermontferrand.fr

Carole Le Guyader
email: carole.le-guyader@insa-rouen.fr

Antoine Vacavant
email: antoine.vacavant@uca.fr

A Time-Dependent Joint Segmentation and Registration Model: Application to Longitudinal Registration of Hepatic DCE-MRI Sequences

While segmentation consists in partitioning a given image into meaningful constituents in order to identify relevant structures such as homogeneous regions or edges, registration, given two images, aims at finding an optimal orientationpreserving one-to-one deformation aligning the structures visible in an image into the corresponding ones in the other. Recently, intertwining both tasks into a single framework has proven to yield better results in terms of accuracy -in particular when the images exhibit weak boundary definition -and increase of reliability of the encoded structure matching -since now, not only based on intensity distribution comparison but also on geometrical and topological features -. In line with this idea, we propose going a step further by adding explicitly some dynamics in the modelling, i.e., by making the minimization problem both space and time-dependent so that the correlation between both tasks is achieved through the process, connecting thus the problem to an interpolation one. The shapes to be matched are viewed as Saint Venant-Kirchhoff materials, a special instance of hyperelastic ones, and are implicitly modelled by level-set functions. These are evolved in order to minimize a functional including both a nonlinear-elasticity-based regularizer prescribing the physical nature of the deformation and a term penalizing the shape misalignment, thus promoting structure matching rather than intensity pairing. Theoretical results emphasizing the mathematical soundness of the model are provided, among which the existence of minimizers and the existence of a weak viscosity solution to the related evolution problem. The model is then applied to the longitudinal registration of hepatic dynamic contrast-enhanced MRI sequences and shows good performance. This application has an important impact on the computer-aided

I. INTRODUCTION

Multitask frameworks and especially variational ones have demonstrated significant improvements over sequential approaches. The former cover a large spectrum of image processing problems including combined segmentation and registration models (see references herein below); joint image reconstruction and motion estimation [START_REF] Aviles-Rivero | Compressed sensing plus motion (CS+M): A new perspective for improving undersampled MR image reconstruction[END_REF], [START_REF] Burger | A variational reconstruction method for undersampled dynamic X-ray tomography based on physical motion models[END_REF], [START_REF] Chun | Joint Image Reconstruction and Nonrigid Motion Estimation with a Simple Penalty that Encourages Local Invertibility[END_REF], [START_REF] Tomasi | Shape and motion from image streams under orthography: a factorization method[END_REF]; joint reconstruction and registration for post-acquisition motion correction [START_REF] Corona | Multi-tasking to Correct: Motion-Compensated MRI via Joint Reconstruction and Registration[END_REF] with the goal to reconstruct a single motionfree corrected image and retrieve the physiological dynamics through the deformation maps, joint optical flow estimation with phase field segmentation of the flow field [START_REF] Brune | Detection of intensity and motion edges within optical flow via multidimensional control[END_REF], or joint segmentation/optimal transport models [START_REF] Boink | Combined modelling of optimal transport and segmentation revealing vascular properties[END_REF] (to determine the velocity of blood flow in vascular structures). Sharing representation between tasks and carefully intertwining them allows to reduce error propagation, to create synergies, to compensate for some possible flaws such as image quality impairment, while increasing the accuracy of the outcomes and bridging the gap towards generalization. Joint segmentation and registration models such as [START_REF] Droske | Multiscale joint segmentation and registration of image morphology[END_REF], [START_REF] Han | Mumford-Shah Model for One-to-One Edge Matching[END_REF] (joint phase field approximation and registration), [START_REF] Lord | Simultaneous Registration and Parcellation of Bilateral Hippocampal Surface Pairs for Local Asymmetry Quantification[END_REF] (model based on metric structure comparison), [START_REF] Debroux | A non-local topologypreserving segmentation guided registration model[END_REF], [START_REF] Swierczynski | A level-set approach to joint image segmentation and registration with application to CT lung imaging[END_REF] (level set formulation that merges the piecewise constant Mumford-Shah model with registration principles), [START_REF] Gooya | GLISTR: Glioma Image Segmentation and Registration[END_REF] (grounded in the expectation maximization algorithm), [START_REF] Debroux | A joint segmentation/registration model based on a nonlocal characterization of weighted total variation and nonlocal shape descriptors[END_REF] (based on a nonlocal characterization of weighted-total variation and nonlocal shape descriptors), or [START_REF] An | Lecture Notes in Comput. Sci. 3749, chapter A Variational PDE Based Level Set Method for a Simultaneous Segmentation and Non-Rigid Registration[END_REF], [START_REF] Guyader | A combined segmentation and registration framework with a nonlinear elasticity smoother[END_REF], [START_REF] Ozeré | Joint Segmentation/Registration Model by Shape Alignment via Weighted Total Variation Minimization and Nonlinear Elasticity[END_REF], [START_REF] Rumpf | A nonlinear elastic shape averaging approach[END_REF], [START_REF] Vemuri | Image registration via levelset motion: applications to atlas-based segmentation[END_REF], fall within this framework.

Registration can be seen as the incorporation of prior information -such as topological constraints for instance, since the one-to-one property of the deformation makes the moving shape homeomorphic to the original one -in the segmentation process. By the same token, accurate segmented structures enable one to drive the registration process correctly by transferring edges for instance. The primary scope of the paper is thus to go one step further in this bi-task formalism by explicitly introducing some dynamics in the modelling, i.e., by making the optimization problem both space and timedependent phrased then on a Sobolev space of Banach-spacevalued functions. The model promotes large deformations in a nonlinear-elasticity-based setting, and includes a fidelity measure fostering shape overlaps and relying on nonlocal shape descriptors from the piecewise constant Mumford-Shah model in a spatio-temporal setting.

This model is then applied to the longitudinal computeraided diagnosis and follow-up of hepatic diseases. In clinical routine, Dynamic Contrast Enhancement MRI (DCE-MRI) is a standard for liver cancers such as HCC (Hepato-Cellular Carcinoma) for instance [START_REF] Balogh | Hepatocellular carcinoma: a review[END_REF]. In DCE-MRI four main phases may be visualized: at non-injected phase (phase I), the lesion is difficult to differentiate from healthy tissues; at contrasted arterial phase (phase II), the HCC is enhanced when compared with the surrounding liver tissues; and at the portal and delayed phases (phases III and IV), the signal of the HCC appears slightly lower than the surrounding liver, defined as wash out. In this context, our joint model serves as a registration tool of such sequences and shows good behavior with high Dice coefficient values. As a long-term vision, this will surely improve computerized analysis of DCE-MRI volumes in the follow-up of patients suffering from liver cancers as HCC.

II. MATHEMATICAL MODELLING A. Mathematical background

Let Ω be a connected bounded open subset of R 2 of class C 1 . Let us denote by R : Ω → R the Reference image assumed to be sufficiently smooth and by T : Ω → R the Template image. The shape contained in the Template image is assumed to be modelled by a Lipschitz continuous function Φ 0 (input of the problem) whose zero level line is the shape boundary. Denoting by C the zero level set of Φ 0 and by

w ⊂ Ω the open set it delineates, Φ 0 is chosen such that C = {x ∈ Ω | Φ 0 (x) = 0}, w = {x ∈ Ω | Φ 0 (x) > 0} and Ω \ w = {x ∈ Ω | Φ 0 (x) < 0}.
For theoretical and numerical purposes, we may consider a linear extension operator (see [6, p. 158

]) P : W 1,∞ (Ω) → W 1,∞ (R 2 ) such that for all Φ ∈ W 1,∞ (Ω), (i) P Φ |Ω = Φ, (ii) P Φ L ∞ (R 2 ) ≤ C Φ L ∞ (Ω) and (iii) P Φ W 1,∞ (R 2 ) ≤ C Φ W 1,∞ (Ω)
, with C depending only on Ω. By this extension process, we consider then that Φ 0 ∈ W 1,∞ (R 2 ) to ensure that Φ 0 • ϕ -with ϕ introduced later -is always defined. We recall that, in the general case, if U is an open subset of R N , for 1 ≤ p ≤ +∞, the Sobolev space denoted by W 1,p (U ) consists of the functions in L p (U ) whose partial derivatives up to order 1, in the sense of distributions, can be identified with functions in L p (U ). Let ϕ : Ω → R 2 be the sought deformation. A deformation is a smooth mapping that is orientation-preserving and injective, except possibly on ∂Ω. The deformation gradient is ∇ϕ: Ω → M 2 (R), the set M 2 (R) being the set of all real square matrices of order 2 identified to R 4 . This deformation to be recovered is seen as the optimal solution of a specifically designed cost function, comprising a regularization on ϕ prescribing the nature of the deformation, and a term measuring alignment or how the available data are exploited to drive the registration process. These are depicted hereafter.

B. Deformation regularization

Nonlinear elasticity principles dictate the design of the smoothness on ϕ. The shapes to be matched are viewed as hyperelastic materials (capable of undergoing large deformations), and more precisely as Saint-Venant Kirchhoff ones. The considered regularizer is given, setting ξ = ∇ϕ, by :

QW (ξ) = W SV K (ξ) + µ (det ξ -1) 2 if ξ 2 ≥ α, Ψ(det ξ) if ξ 2 < α, (1)
with

W SV K (F ) = λ 2 (tr E) 2 +µ tr E 2
, the stored energy function of a Saint Venant-Kirchhoff material, λ and µ the Lamé coefficients, E = F T F -I /2 the Green-Saint Venant stress tensor measuring the deviation between ϕ and a rigid deformation, and with the following notation A : B = trA T B, the matrix inner product, and A = √ A : A the related matrix norm (Frobenius norm). Additionally, α = 2 λ+µ λ+2µ and Ψ is the convex mapping defined by Ψ :

s → -µ 2 s 2 + µ (s - 1) 2 + µ(λ + µ) 2(λ + 2µ) :=γ
. Several arguments motivate this choice by comparison to the stored energy function W SV K alone: (i) first, although the stored energy function W SV K is the simplest one that agrees with the generic expression of the stored energy function of an isotropic, homogeneous, hyperelastic material, it lacks a term penalizing the determinant. It thus does not preclude deformations with negative Jacobian determinant; (ii) second, it is not rank-1 convex and consequently not quasiconvex (see [START_REF] Dacorogna | Direct Methods in the Calculus of Variations[END_REF]Chapter 9]), which raises some theoretical issues in terms of existence of minimizers, contrary to QW which exhibits this fine latter property; (iii) thirdly, we see in the expression of QW that when ξ 2 < α, a penalization on the determinant still remains, showing good behavior under compression. By introducing explicitly the time variable t ∈ [0, T ] in the model and the unknown ϕ = ϕ(x, t) being now space and time-dependent, it yields the following regularization on ϕ

Reg(ϕ) = V ∂ϕ ∂t 2 2 dx dt + V QW (∇ x ϕ(x, t)) dx dt, (2)
• 2 denoting the Euclidean norm in R 2 and with V = Ω × [0, T ].

C. Dissimilarity measure

Complying with the idea of promoting structure/shape matching rather than intensity-based rules, we introduce a metric fostering the overlap of shapes as follows. Recall that the shape contained in the Template image is assumed to be modelled by a Lipschitz continuous function Φ 0 whose zero level line is the shape boundary. The idea is to deform the function Φ 0 so that the zero level line of the deformed function, Φ 0 • ϕ, is aligned with the boundary of the shape belonging to the Reference. The trade-off between segmentation and registration -which makes our model a bi-task one -is ensured by processing at the same time the segmentation of the Reference via the level set function Φ (see 4) and the deformation of Φ 0 , and these tasks are connected in a single framework to maximize the overlap of the two resulting shapes. More precisely, the dissimilarity measure is defined by:

Fid(ϕ) = V |H ε (Φ 0 • ϕ(x, t)) -H ε ( Φ(x, t))| 2 dx dt, (3)
H ε denoting a C ∞ -regularization of the one-dimensional Heaviside function and the evolution of Φ being dictated by the following evolution equation, which constitutes a revised version of [START_REF] Guyader | A combined segmentation and registration framework with a nonlinear elasticity smoother[END_REF] to unify local and global (region-based) features to segment the Reference :

                         ∂ Φ ∂t = |∇ Φ| div g(|∇R|) ∇ Φ |∇ Φ| + k (R -c 2 ) 2 -(R -c 1 ) 2 + 4 µ d 2 H( Φ(x) + l) H(l -Φ(x)) Ω x-y,∇ Φ(y) e -x-y 2 2 /d 2 H( Φ(y)+l) H(l-Φ(y)) dy , Φ(x, 0) = Φ 0 (x) , ∂ Φ ∂ n = 0, on ∂Ω , (4)
with

c 1 = Ω H( Φ) R(x) dx Ω H( Φ) dx and c 2 = Ω H(-Φ) R(x) dx Ω H(-Φ) dx
(we dropped the dependency on Φ to lighten the expressions). Φ 0 is naturally taken to be the initial condition of this segmentation process. Function g is an edge-detector function satisfying g(0) = 1, g strictly decreasing and lim r→+∞ g(r) = 0. The first two components of the right hand side of equation ( 4) constitutes a balance between the classical geodesic active contour model ( [START_REF] Caselles | Geodesic Active Contours[END_REF]) and the piecewise constant Mumford-Shah model ( [START_REF] Mumford | Optimal approximation by piecewise smooth functions and associated variational problems[END_REF]) which allows to partition R into two phases w and Ω \ w with respective values c 1 and c 2 . The latter one prescribes some topological constraints ( [START_REF] Guyader | A combined segmentation and registration framework with a nonlinear elasticity smoother[END_REF]) in order for the evolving contour to be homeomorphic to the original one. The overall registration problem thus reads :

J(ϕ) = V ∂ϕ ∂t (x, t) 2 2 dx dt + V [QW (∇ x ϕ(x, t)) (5)
+ ν 2 H ε (Φ 0 • ϕ(x, t)) -H ε ( Φ(x, t)) 2 dx dt, with V = Ω × [0, T ].

D. Theoretical results

We state the main theoretical result related to the existence of minimizers.

Theorem 1 (Existence of minimizers): There exists at least one minimizer of problem [START_REF] Boink | Combined modelling of optimal transport and segmentation revealing vascular properties[END_REF]

in the functional space W = ϕ ∈ L 4 (0, T ; W 1,4 (Ω, R 2 )) | ∂ϕ ∂t ∈ L 2 (0, T ; L 2 (Ω, R 2 )) endowed with the norm ϕ W = ϕ L 4 (0, T ;W 1,4 (Ω,R 2 )) + ∂ϕ ∂t L 2 (0, T ;L 2 (Ω,R 2 )) .
Proof: The proof is divided into the following three steps and relies on Aubin-Lions lemma:

• Coercivity inequality: By taking ϕ(t) = Id ∀t ∈ [0; T ],

we get J(ϕ) < +∞ and the functional is proper. We first derive a coercivity inequality:

J(ϕ) ≥ ∂ϕ ∂t 2 L 2 (]0, T [×Ω,R 2 ) + β 2(c 4 + 1) ϕ 4 L 4 (0, T ;W 1,4 (Ω,R 2 )) -βα 2 meas(Ω) T - µ(λ + 3µ) 2(λ + 2µ) meas(Ω) T + γmeas(Ω) T + k T .
So, the infimum is finite. • Convergence of a minimizing sequence: Let {ϕ k } k∈N be a minimizing sequence. As the functional is proper, there exists φ such that for all k large enough, J(ϕ k ) ≤ J( φ) + 1 < +∞. Then from the coercivity inequality, we can deduce that {ϕ k } k∈N is uniformly bounded in W so we can extract a subsequence still denoted {ϕ k } such that ϕ k φ in W . As W 1,4 (Ω, R 2 ) ⊂ L 2 (Ω, R 2 ) with compact injection, Aubin-Lions lemma states that the embedding of W in L 4 (0, T ;

L 2 (Ω, R 2 )) is compact. As L 4 (0, T ; L 2 (Ω, R 2 )) ⊂ L 2 (]0, T [×Ω, R 2 ) is continuous then the embedding of W into L 2 (]0, T [×Ω, R 2 ) is also compact. We can therefore extract a subsequence of {ϕ k } still denoted {ϕ k } such that ϕ k → φ in L 2 (]0, T [×Ω, R 2 ).
• Weak lower semi-continuity: . L 2 (]0, T [×Ω,R 2 ) is convex and strongly continuous and so it is weakly lower semi-continuous in L

2 (]0, T [×Ω, R 2 ) so ∂ φ ∂t L 2 (]0, T [×Ω,R 2 ) ≤ lim inf k→+∞ ∂ϕ k ∂t L 2 (]0, T [×Ω,R 2 ) . As, ϕ k → φ in L 2 (]0, T [×Ω, R 2 )
, up to a subsequence, one has pointwise convergence of {ϕ k } to φ and the dominated convergence theorem enables us to obtain the weak lower semi-continuity of the data fidelity term. Let {ψ k } k∈N be a sequence that strongly converges to ψ in L 4 (0, T ; W 1,4 (Ω, R 2 )).

Then T 0 ψ(t) -ψ k (t) 4 W 1,4 (Ω,R 2 ) dt -→ k→+∞ 0. By seeing ψ(t) -ψ k (t) 4 W 1,4 (Ω,R 2 )
as a real-valued function depending on t defined on ]0, T [ and by applying the reciprocal of the dominated convergence theorem, we get that ψ(t) -ψ k (t) 4 W 1,4 (Ω,R 2 ) converges to 0 for almost every t ∈]0, T [ up to a subsequence. So for almost every t ∈]0, T [, ψ k (t) strongly converges to ψ(t)

in W 1,4 (Ω, R 2 ) and det ∇ψ k (t) → k→+∞ det ∇ ψ(t) in L 2 (Ω).
From what was done in the stationary case [START_REF] Debroux | A non-local topologypreserving segmentation guided registration model[END_REF], for almost every t ∈]0, T [,

Ω QW (∇ ψ(t)) dx ≤ lim inf k→+∞ Ω QW (∇ψ k (t)) dx.
So by Fatou's lemma, we get:

T 0 Ω QW (∇ ψ(x, t)) dx dt ≤ lim inf k→+∞ T 0 Ω QW (∇ψ k (x, t)) dx dt.
Then it is convex and strongly sequentially lower semicontinuous and so it is weakly lower semi-continuous.

We finally have

T 0 Ω QW (∇ φ(x, t)) dx dt ≤ lim inf k→+∞ T 0 Ω QW (∇ϕ k (x, t)) dx dt,
which concludes the proof.

The well-definedness of Φ is then investigated to ensure that the fidelity term exhibits sufficient regularity and makes sense. Equation ( 4) falls within the framework of viscosity solution theory ( [START_REF] Barles | Existence of weak solutions for general nonlocal and nonlinear second-order parabolic equations[END_REF], [START_REF] Crandall | User's Guide to Viscosity Solutions of Second Order Partial Differential Equations[END_REF]) for equations with a measurable dependence in time (called L 1 -viscosity solution). Under mild assumptions, the following theorem holds. Theorem 2 (Existence of weak solutions of the considered evolution problem): Assuming that g := g(|∇R|), g 1 2 and ∇g are bounded and Lipschitz continuous on R 2 , for fixed c 1 and c 2 , problem (4) admits at least one weak solution.

III. NUMERICAL RESOLUTION

We now aim to apply this theoretical model to the longitudinal study of an image sequence. Let (I n ) N n=0 be the temporal image sequence. The segmentation evolution problem now reads :

                                 ∂ Φ ∂t = |∇ Φ| div g(|∇I n |) ∇ Φ |∇ Φ| + k((I n -c 2,n ) 2 -(I n -c 1,n ) 2 ) + 4 µ d 2 H( Φ(x) + l) H(l -Φ(x)) Ω x-y,∇ Φ(y) e - x-y 2 d 2 H( Φ(y)+l) H(l-Φ(y)) dy, c 1,n = Ω H( Φ)In(x) dx Ω H( Φ) , c 2,n = Ω In(x) H(-Φ) dx Ω H(-Φ) dx , for t ∈ [t n-1 , t n ], N ≥ n ≥ 1, Φ(x, 0) = Φ 0 (x), ∂ Φ ∂ n = 0, on ∂Ω. (6)
For the registration problem, rather than considering a continuum in time which is not realistic, we assume that the problem is sampled in time and drop the regularization in time introduced only for theoretical purposes. We thus solve sequentially the subproblems :

inf ϕi∈Id+W 1,4 0 (Ω,R 2 ) Ω QW (∇ϕ i ) dx + ν 2 Ω H ε (Φ 0 • ϕ 1 • • • • • ϕ i ) -H ε ( Φ(•, t i )) 2 dx,
for i ∈ {1, . . . , N ζ}, ζ being the number of steps saved in the segmentation process between t n-1 and t n . In the end, the overall deformation is given by ϕ

1 • • • • • ϕ N ζ
, and the deformation between the initial frame and the n-th frame by

ϕ 1 • . . . • ϕ nζ .
Remark 1: From a theoretical standpoint, the existence of minimizers for each subproblem is guaranteed: Rellich-Kondrachov's embedding theorem states that weak convergence in W 1,4 (Ω, R 2 ) leads to uniform convergence in Ω, an extension process as before can be applied on all ϕ k , k = 1, • • • , i-1 to ensure the well-definedness of the composition, and the continuous injection

W 1,4 (R 2 , R 2 ) C 0 (R 2 , R 2 )
holds, these three elements combined allowing to handle the fidelity term. In order to deal with the nonlinearity in ∇ϕ, we propose introducing an auxiliary variable V i simulating the Jacobian deformation with a quadratic penalty method as in [START_REF] Debroux | A non-local topologypreserving segmentation guided registration model[END_REF]. The decoupled problem becomes :

inf ϕi∈Id+W 1,2 0 (Ω,R 2 ),Vi∈L 4 (Ω,M2(R)) Ω QW (V i ) dx + ν 2 Ω H ε (Φ 0 • ϕ 1 • • • • • ϕ i ) -H ε ( Φ(•, t i )) 2 dx + γ 2 ∇ϕ i -V i 2 L 2 (Ω,M2(R)) , (7)
for i ∈ {1, . . . , N ζ}. We then use an alternating minimization scheme to solve the problem. We refer the reader to [15, Section 4.3.] for an exhaustive description of the algorithm relying on the derivation of Euler-Lagrange equations solved by an L 2 gradient flow algorithm and an implicit/semi-implicit Euler time stepping scheme.

IV. NUMERICAL SIMULATIONS

The proposed method has been evaluated both qualitatively and quantitatively as a registration process to the longitudinal analysis of 2 DCE-MRI sequences (composed of 4 MRI phases each denoted I, II, III and IV) for a patient suffering from cirrhosis and HCC. DCE-MRI exams were performed on a 1.5T SIGNA™ Artist (General Electric, Milwaukee, WI) with a phased array coil. Manual segmentations of the liver have been elaborated thanks to a Slicer 3D plug-in we are developing [START_REF] Lamy | Preliminary results with a new annotation tool for liver volume and inner vessels from DCE-MRI data[END_REF]. The segmentation step of our method serves here as an interpolation process between the manual segmentations to guide the registration and improve both the matching quality and the topology-preservation property.

To assess the performance of our algorithm, we display the visual outputs in Fig. 1 and use two metrics: the Dice coefficient measuring set agreement with the highest score 1; the mutual information to measure image alignment with larger value meaning better matching. These quantitative measures are reported in Table I along with the max and min value of the deformation Jacobian determinant to appraise the topology preservation. We can see in Fig. 1 that the deformed contour is well aligned with the edge of the liver at all times and the deformation grid does not exhibit overlaps, highlighting the topology preservation property of our algorithm. This is particularly visible in the zoom-in views of the complex topology exhibited by the liver around the vena cava: the segmentation delineates well the concavity without violating topology preservation. It thus shows the capability of our model to handle large deformations on complex shapes and intensity variations between frames. This is further confirmed by the analysis of the Dice coefficients between the deformed initial frame (T • ϕ) and the corresponding temporal frame (R) and their mutual information. Indeed, the latter always significantly increases between the non-deformed frame T and the deformed one T • ϕ, while the Dice coefficient is greater than 98.7% at all times. Favoring shape matching rather than intensity pairing has therefore a positive impact on the registration quality when intensity changes are involved between the temporal frames. Furthermore, we can see in Table I that the determinant remains positive at all times which translates into the ability of our model to generate physically meaningful deformations along the longitudinal analysis. This experiment is thus a proof of concept that our method can be used for the longitudinal analysis of complex organs such as the liver including intensity variations.

V. DISCUSSION AND FUTURE WORKS

In this article, we have proposed a novel joint registrationsegmentation model that minimizes a functional by integrating both space and time-dependent terms. Our contribution enables one to generate large deformations thanks to a nonlinear elasticity regularization process and promotes structure matching rather than intensity-based similarity measuring.

We have successfully exploited our model in a computerized diagnosis and follow-up process by computing the registration of liver shapes within DCE-MRI sequences. By construction, these MRI volumes have variable intensity ranges (due to contrast agent absorption), which would have penalized standard approaches that consider more pixel values than geometrical or topological structures. In this context, we can handle concavities that are generally considered as strong topological constraints for registration. Our experiments show that our model can offer accurate performance by means of excellent Dice coefficient and mutual information measures.

As future works, we first plan to increase the capability of our model by developing a fully 3D approach rather than a slice-based one, and by studying multi-modal joint registration-segmentation through the context of personalized MRI/CT alignment. Moreover, we would like to tackle the problem of internal tissue registration by studying tumors throughout dynamic sequences. In this case, we can extend our framework by considering two level sets: one for the liver volume, and another one for the tumor. We can also produce an approximate registration of the tumors by applying the deformations obtained from the whole liver to internal tissues. Finally, we can exploit further our model by calculating a more precise segmentation of the liver and internal lesions within medical 3D volumes.

Fig. 1 .

 1 Fig. 1. Visual assessment of the longitudinal analysis performed on the sequence 19-02-2013 Phase I to IV followed by 08-08-2013 Phase I to IV. For each time step, the image with the initial contour along with the image with the deformed contour and the deformation grid are shown. Zoom-in views are provided below to highlight the preservation of topology in complex shapes such as concavities, and the quality of the alignment between the deformed contour and the image edges.

TABLE I QUANTITATIVE

 I EVALUATION OF THE REGISTRATION IN THE LONGITUDINAL STUDY COMPRISING THE FOUR STEPS OF TWO DCE MRI EXAMS WITH TWO METRICS: DICE COEFFICIENT, MUTUAL INFORMATION AND MIN/MAX OF THE DEFORMATION JACOBIAN DETERMINANT.

This work was co-funded by the French Agence Nationale de la Recherche (grant ANR-18-CE45-0018, project R-Vessel-X, http://tgi.ip.uca.fr/r-vesselx), by the European Union with the European Regional Development Fund (ERDF, 18P03390/18E01750/18P02733) and by the Haute-Normandie Regional Council via the M2SINUM project.

follow-up of patients suffering from liver cancers.


Robert J Hargreaves
email: robert.hargreaves@cfa.harvard.edu

Iouli E Gordon

Michael Rey

Andrei V Nikitin

Vladimir G Tyuterev

Roman V Kochanov

Laurence S Rothman

& Lyons

Waldmann, I. P A Tsiaras

Typeset using

Keywords:

) ApJ, 847, 105 (provided at seperate temperatures in the TheoReTS information system), are now combined with HITRAN2016 methane data to produce a single line list suitable for high-temperature line-byline calculations up to 2000 K. An effective-temperature interpolation model was created in order to represent continuum-like features at any temperature of interest. This model is advantageous to previously-used approaches that employ so-called "super-lines", which are suitable only at a given temperature and require separate line lists for different temperatures. The resultant HITEMP line list contains ∼32 million lines and is significantly more flexible than alternative line lists of methane, while accuracy required for astrophysical or combustion applications is retained. Comparisons against experimental observations of methane absorption at high temperatures have been used to demonstrate the accuracy of the new work. The line list includes both strong lines and quasi-continuum features and is provided in the common user-friendly HITRAN/HITEMP format, making it the most practical methane line list for radiative transfer modeling at high-temperature conditions.

INTRODUCTION

On Earth, atmospheric methane (CH 4 ) is a prominent greenhouse gas that has seen a steady increase over the last decade (Fletcher & Schaefer 2019). Terrestrial CH 4 has both natural and anthropogenic sources, with atmospheric monitoring of CH 4 typically achieved using infrared spectral observations (Jacob et al. 2016). CH 4 is also the main constituent of natural gas, and plays a central role in combustion. At high temperatures, CH 4 spectra can be used for diagnostics of hydrocarbon combustion processes throughout the infrared (Nagali et al. (Cushing et al. 2011;Kirkpatrick et al. 2012) and is also present in the atmospheres of the giant planets (Irwin et al. 2005;Mueller-Wodarg et al. 2008) and Titan (Karkoschka 1994;Atreya et al. 2006).

Since the detection of 51 Pegasi b (Mayor & Queloz 1995), there are now in excess of 4000 known exoplanets. Studies of transiting exoplanets have been able to probe the atmospheres of a small number of these objects (Tsiaras et al. 2018), with observations of water vapor (Grillmair et al. 2008) and carbon monoxide absorption (Konopacky et al. 2013). Models predict CH 4 to be more abundant than carbon monoxide below ∼1300 K (Burrows & Sharp 1999), yet observations of CH 4 have only been reported in the spectra of five exoplanets to date: HD 189733b (Swain et al. 2008), HD 209458b (Swain et al. 2009), XO-1b (Tinetti et al. 2010), HR 8799b (Barman et al. 2015) and 51 Eridani b (Macintosh et al. 2015).

Many exoplanet observations have used instruments

with low resolving powers (Brogi & Line 2019), where R = λ/∆λ 200, which can limit the capability to identify individual molecular species. However, recent spectroscopic techniques such as cross-correlation (Snellen et al. 2014) and Doppler tomography (Watson et al. 2019) are able to take advantage of high resolution instruments (R ∼ 25, 000 -100, 000) to definitely confirm detections of H 2 O (Birkby et al. 2017), CO (Snellen et al. 2010), TiO (Nugroho et al. 2017), as well as neutral and ionized atoms (Hoeijmakers et al. 2018), from exoplanet transit spectra. These methods have also highlighted the need for line lists to be both accurate and complete at high resolutions (Hoeijmakers et al. 2015).

The pressing need for improvements to line lists for planetary spectroscopy (including CH 4 ) have been emphasized in a number of review papers (Tinetti et al. 2013;Bernath 2014;Fortney et al. 2016;Tennyson & Yurchenko 2017;Fortney et al. 2019). These improvements are essential to make the most of measurements from the forthcoming Atmospheric Remote-sensing Infrared Exoplanet Large-survey (ARIEL) mission (Tinetti et al. 2018), which is dedicated to exoplanet observations. Furthermore, the James Webb Space Telescope will provide a significant advancement in the capability to characterize exoplanet atmospheres using moderate resolution (R 3500) spectroscopy (Greene et al. 2016).

Methane spectroscopy

The polyad nature of CH 4 is a consequence of all four vibrational modes having the relationship ν 1 ≈ ν 3 ≈ 2ν 2 ≈ 2ν 4 ≈ 3000cm -1 . Each polyad is identified by P n , where n = 2(v 1 +v 3 )+v 2 +v 4 (with v i equal to the number of quanta of each mode), but named according to the number of vibrational bands within each polyad. For 106 example, the second polyad P 2 contains 5 vibrational 107 bands (ν 1 , ν 3 , 2ν 2 , 2ν 4 , ν 2 +ν 4 ), and is therefore referred 108 to as the pentad (Boudon et al. 2006). Due to the tetra- as TDS (Tyuterev et al. 1994), STDS pion 1998) and MeCaSDa (Ba et al. 2013), have been and lower-temperature spectra (Nikitin et al. 2017a, 137 2018; Rodina et al. 2019;Nikitin et al. 2019). Many 138 of these studies, as well as HITRAN2016 updates, have 139 already benefited from supplementary information for 140 the resonance interaction parameters within vibrational 141 polyads. These are derived from an ab initio poten-142 tial energy surface that made analyses of experimen-143 tal spectra more consistent and reliable, as described 144 in Tyuterev et al. (2013). However this was only done 145 for cold bands and for relatively low polyads up to 146 ∼7300 cm -1 . The difficulty of extending assignments 147 is strongly exacerbated at higher temperatures. For 148 this reason, a number of high-temperature laboratory 149 measurements have been made of CH 4 in both emission 150 (Nassar & Bernath 2003;Thiévin et al. 2008;Hargreaves 151 et al. 2012;Amyay et al. 2018a,b;Georges et al. 2019) 152 and absorption (Alrefae et al. 2014;Hargreaves et al. 153 2015; Ghysels et al. 2018;Wong et al. 2019).

154

On the theoretical side, the hot bands and high-J (Rey et al. 2016). These works demonstrated that ab initio line lists of CH 4 could approach the accuracy required for high-temperatures, but the inclusion of billions of transitions made the resulting full line-by-line lists impractical for typical applications. When comparing these line lists, Hargreaves et al. (2015) recommended the separation of strong and continuum-like features. Indeed, it was shown by Rey et al. (2014a) that it is necessary to account for approximately 1 million rovibrational transitions per 1 cm -1 for CH 4 opacity calculations at 2000 K. To make online computations of the absorption cross-section faster, it was suggested to model the quasi-continuum formed by the contributions of huge amounts of very weak lines using so called "super-lines", as originally implemented in the TheoReTS database (Rey et al. 2016). Super-lines represent integrated intensity contributions from tiny transitions on a pre-defined grid of small wavenumber and temperature intervals.

Updated state-of-the-art ab initio line lists have since been published, ExoMol '34to10' (Yurchenko et al. 2017) and Rey et al. (2017) (referred to here as RNT2017), both of them using the super-line approach for the compression of relatively weak absorption/emission features complemented with lists of medium and strong lines.

To obtain the full CH 4 spectrum, both the strong and super-line components are required. In each case, these line lists still require a large quantity of strong lines to cover the temperature range of calculations. Furthermore, a separate super-line component is provided at each temperature, which makes them difficult to integrate into existing radiative transfer codes and significantly less flexible than a standard line list.

The HITRAN and HITEMP databases

The HITRAN database contains detailed spectroscopic line-by-line parameters of 49 molecules with many of their isotopologues (along with absorption crosssections for almost 300 molecules, collision-induced absorption spectra for many collisional pairs, and aerosol properties). HITRAN2016 (Gordon et al. 2017) is the most recent version of the database, and is freely available at HITRANonline 1 . Recent efforts have been undertaken to expand the use of HITRAN towards planetary atmospheres, with the inclusion of additional broadening species (Wilzewski et al. 2016;Tan et al. 2019). However, the CH 4 line list in HITRAN2016 is un- (cm -1 ), between two rovibronic states is given as

233 S ij (T ) = A ij 8πcν 2 ij g I a Q(T ) exp -c 2 E T 1 -exp -c 2 ν ij T ,
(1) where A ij (s -1 ) is the Einstein coefficient for sponta-234 neous emission, g is the upper state statistical weight, 235 E (cm -1 ) is the lower-state energy, Q(T ) is the total 236 internal partition sum, I a is the natural terrestrial iso-

242

The units 3 used throughout HITRAN editions do not 243 2 One should note that isotopic abundance is dependent upon the environment and HITRAN is consistent with specific terrestrial values given by De Biévre et al. (1984). For applications that do not assume these isotopic mixtures (e.g., exoplanetry atmospheres), this weighting should be renormalized by the user. 3 Line positions in HITRAN and HITEMP are provided in reciprocal centimeter (cm -1 ) and denoted ν (thereby dropping the tilde that is the official designation of wavenumber, ν), and pressure in atm (atmosphere). Intensity is traditionally expressed as cm -1 /(molecule cm -2 ) rather than simplifying to the equivalent cm molecule -1 .

strictly adhere to the SI system for both historical and application-specific reasons.

The HITRAN Application Programming Interface, HAPI (Kochanov et al. 2016) The RNT2017 line list was created in three steps. The first was to provide over 150 billion transitions (with a lower-state rovibrational energy cutoff of 33,000 cm -1 ) from first-principles quantum mechanical variational calculations using the molecular potential energy surface of Nikitin et al. (2011Nikitin et al. ( , 2016)). The line intensities were calculated from the purely ab initio dipole moment surfaces of Nikitin et al. (2017b). The resulting line list ranges from 0-13,400 cm -1 (i.e., > 746 nm) with a maximum temperature of 3000 K.

To improve the accuracy of the ab initio line positions, a second step applied empirical corrections for 3.7 million of the strongest transitions. This involves ∼100,000 Table 1. Summary of the individual 12 CH4 line lists used in this work from Rey et al. (2017). At each temperature, the number of lines (N RNT (T )) and intensity sums (ΣS RNT (T )) are given for the total line list, along with the strong and super-line components.

energy
ΣS RNTstr (T ) b , ΣS RNTsup (T ) b , ΣS RNTtot (T ) b , T νmax a N RNTstr (T ) ×10 -17 N RNTsup (T ) ×10 -18 N RNTtot (T ) ×10 -17 (K) (cm -1 ) (cm -1 /(molecule cm -2 )) (cm -1 /(molecule cm -2 )) (cm -1 /(molecule cm -2 )) 300
S ij (T ) S ij (T 0 ) = Q(T 0 ) Q(T ) exp c 2 E T 0 - c 2 E T 1 -exp(-c 2 ν ij /T ) 1 -exp(-c 2 ν ij /T 0 )
(2) where T 0 = 296 K for the HITRAN and HITEMP 2) can be rearranged as

505 sity I 1 at T 1 , if I strong (ν 1 , T 1 ) > I 1 > I cut (ν 1 , T 1 ) then I 1
537 ln S ij (T )Q(T )R(T 0 ) S ij (T 0 )Q(T 0 )R(T ) = c 2 E T 0 - c 2 E T (3)
where Table 2. Intensity sums for the total 12 CH4 line list from this work (ΣS TWtot (T )), along with the strong and effective components. For comparison with Tab. (1), these intensity sums have been calculated at the same temperatures.

R(T ) = 1 -exp(-c 2 ν ij /T ). Thus, a plot of 538 ln[S ij (T )Q(T )R(T 0 )/S ij (T 0 )Q(T 0 )R(T )] against -c 2 /
ΣS

TWstr (T ) a , ΣS TWeff (T ) a , ΣS TWtot (T ) a , ΣS TWtot (T ) / T ×10 -17 ×10 -17 ×10 -17 ΣS RNTtot (T ) (K) (cm -1 /(molecule cm -2 )) (cm -1 /(molecule cm -2 )) (cm -1 /(molecule cm -2 )) (%) 300

ROOM TEMPERATURE COMPARISONS

The resulting HITEMP line list from this work can be compared to high-resolution absorption cross-sections of HITRAN2016 at 296 K, calculated using the HAPI routines (Kochanov et al. 2016). These comparisons can test the accuracy of the known line positions as well as the validity of unassigned features. (i.e., R ∼ 20, 000 for the pentad and R ∼ 28, 000 for the icosad). It can be expected that this uncertainty will increase for higher temperatures. A diagram for the completeness of the 12 CH4 HITEMP line list. Each colored region identifies the wavenumber range of completeness (up to the limits indicated) when used at high temperature. For T < 800 K, the HITEMP line list is complete between ν = 0-13,400 cm -1 , whereas at 2000 K, it is complete between ν = 0-6600 cm -1 . Use of the 12 CH4 HITEMP line list at temperatures/wavenumbers within the grey region will lead to errors as a consequence of incompleteness and is not recommended.

eters of CH 4 currently available for use at high temper- vations of transiting exoplanets are sensitive to these window regions (in addition to the band centers) as light passes through the limb of the planet's atmosphere. The contribution to these windows at high temperatures is caused by the continuum lines, and sensitive measurements of the CH 4 continuum-like absorption would be valuable for validation. Application of this work to simulations/retrievals of brown dwarfs and exoplanets will assist in the validation of these regions at the highest temperatures.

The completeness of RNT2017 is clearly described in Rey et al. (2017), where efforts were made to keep each polyad complete up to the temperatures and spectral limits given in Tab. (1). Our current work has been based on the RNT2017 line lists and therefore retains its completeness. The temperature-dependant wavenumber limits are also applicable to this work (and therefore the CH 4 HITEMP line list), which have been illustrated in Fig. (11). For temperatures of spectral regions outside of these bounds, a scaling factor was recommended for the RNT2017 super-line lists to extrapolate the absorption and compensate for the lack of hot lines in these regions. However, scaling of the effective lines is not recommended since discontinuities can appear at the wavenumber limits, as shown in Fig . (3) at 12,600 cm -1 . These discontinuities are a consequence of the effective lines being retrieved from fewer temperatures. For example, a super-line grid point <6600 cm -1 can be populated at all temperatures, whereas one >12,600 cm -1 is only populated up to 800 K, which reduces the quality of the effective line fits. For users who require the full spectroscopic detail of the complete RNT2017 line lists, we refer to the original work (Rey et al. 2017) available from the TheoReTS data system (Rey et al. 2016).

To create an accurate and complete line list of CH 4 requires the calculation of billions of transitions. This makes them difficult to use in their entirety, and RNT2017 (and ExoMol 34to10) have attempted to mitigate this effect by retaining a relatively small number of strong lines that account for the structure of polyad bands, while compressing the remaining calculated transitions into super-lines. However, these super-lines are temperature dependent and lack flexibility due to the loss of unique line information. For RNT2017, this significantly increases the total number of lines required to cover the full temperature range because a separate list is required at each temperature. In addition, this In this present study, we adopt a different strategy for spectral data compression at high temperatures. To this end we have combined the most accurate theoretical calculations for CH 4 to date (i.e., RNT2017), with HITRAN2016 parameters. Technically, this task is non-trivial due to separation between "strong" and "weak" transitions. In the TheoReTS database the strong line lists are temperature specific, accounting for the Boltzmann population of lower levels. For the present HITEMP line list, we have combined all strong lines from RNT2017 at various temperatures into a single line list. That is to say, all transitions which could result in sufficiently sharp features at temperatures up to 2000 K are included. In addition, we have determined "effective lines" that can model the quasi-continuum features. These effective lines should not be confused with the previously described super-lines, even though they have been determined from the super-line lists. The main difference is that we have attributed an effective lower-state energy to each of the effective lines (i.e., it does not correspond to true quantum state or individual ro-vibrational transition), which enables them to be converted to any temperature using the standard formulae.

109

  hedral symmetry of the CH 4 molecule, the degenerate 110 overtone and combination vibration states involved in 111 successive polyads are split into sub-levels, which com-112 plicates ro-vibrational band patterns for analyses. Early 113 versions of spectroscopic databases specifically devel-114 oped for CH 4 and other high-symmetry molecules, such 115

  (Gordon et al. 2017) de-120 tails the most accurate collection of line parameters for 121 CH 4 , with a primary focus towards the modeling of the 122 terrestrial atmosphere. This is also the focus of the 123 GEISA (Jacquinet-Husson et al. 2016), MeCaSDa (Ba 124 et al. 2013) and GOSAT (Nikitin et al. 2015b) databases. 125 These linelists, which are based on experimental mea-126 surements and/or empirical fits of laboratory spectra, 127 suffer from incompleteness issues for high-temperature 128 conditions because of insufficient information on exper-129 imentally measured and assigned transitions. They are 130 therefore unsuitable for astrophysical applications with 131 a large range of temperatures. 132 Assigning individual transitions becomes a significant 133 challenge in dense spectra with numerous blended fea-134 tures, as is the case for CH 4 . Since HITRAN2016, there 135 has been steady progress in assigning room-temperature 136

  155 transitions have been included in global variational CH 4 156 line lists : '10to10' (Yurchenko & Tennyson 2014), as part of the ExoMol project (Tennyson et al. 2016), and the Rey et al. (2014a) line list (referred to here as RNT2014) as a part of TheoReTS project

1

  https://hitran.org suitable for spectroscopy at high temperatures due to is-206 sues of incompleteness. This is a consequence of the ab-207 sence of many vibrational hot bands, high ro-vibrational 208 transitions or any other extremely weak transitions (at 209 terrestrial temperatures), due to their negligible effect 210 in terrestrial atmospheric applications. 211 The HITEMP database (Rothman et al. 2010) was 212 established specifically to model gas-phase spectra in 213 high-temperature applications, and can be thought of as 214 a "sister" to HITRAN (with data also provided through 215 HITRANonline). One substantial difference between 216 HITRAN and HITEMP is the number of transitions 217 included for each molecular line list, a consequence of 218 the inclusion of numerous vibrational hot bands, high 219 ro-vibrational transitions and overtones. This differ-220 ence is most apparent for H 2 O, where there are cur-221 rently ∼800 times the number of lines in HITEMP2010 222 when compared to HITRAN2016. Typically, these 223 additional transitions constitute numerous lines (often 224 millions) from ab initio or semi-empirical calculations, 225 which are then combined with accurate parameters from 226 HITRAN. The HITEMP database has been undergoing 227 a large scale update (Li et al. 2015; Hargreaves et al.228 2019) and, prior to this work, included seven molecules: 229 H 2 O, CO 2 , N 2 O, CO, NO, NO 2 , and OH. 230 For HITRAN and HITEMP, the temperature-231 dependent spectral line intensity of a transition, ν ij 232

  237 topic abundance 2 , and c 2 = hc/k = 1.4387770 cm K, 238 the second radiation constant. To remain consistent, 239 the spectroscopic parameters in HITRAN and HITEMP 240 are provided at a reference temperature of 296 K and 241 the line intensities are scaled to terrestrial abundances.

  levels extracted from analyses of experimental laboratory room-temperature spectra. No empirical corrections were applied to line intensities, which were computed from an ab initio dipole moment surface using a variational method. A third and final step follows the recommendation of Hargreaves et al. (2015) to separate the empirically corrected line lists into two components: "strong" and "super" lines. To obtain the full CH 4 spectrum at each temperature, both the strong and super-line lists are required. The number of lines in each subsequent line list (at each temperature) is shown in Tab. (1). Full details are described by Rey et al. (2017) with only important points explained here. From the billions of transitions that are computed, an intensity cutoff function, I cut (ν, T ), is used to exclude the weakest transitions that have a negligible contribution to the total opacity at each temperature. The cutoff function has the approximate structure of an extremely low-resolution CH 4 spectrum and is dependent on the wavenumber and temperature. To separate between strong and super-lines at each temperature, a temperature-dependant scale factor (α strong (T )) is applied to the cutoff functions such that I strong (ν, T ) = α strong (T )I cut (ν, T ). All transitions that have an intensity I(ν, T ) > I strong (ν, T ) are retained for the strong line lists. These strong lines are necessary for accurate simulation of sharp features in absorption/emission spectra. Transitions that have an intensity I strong (ν, T ) > I(ν, T ) > I cut (ν, T ) are compressed into so-called super-lines (Rey et al. 2016). These superlines are provided on a 0.005 cm -1 grid and account for billions of weak transitions. The compression of the full line list at each temperature reduces the number of lines necessary for line-by-line calculations and increases the efficiency of radiative transfer calculations. However, the downside of this compression means that the parameters of individual contributing transitions are not stored (e.g., ν, I, E , J ). It is also worth noting that

Figure 1 .

 1 Figure 1. The intensities and positions of strong and superlines from RNT2017 (Rey et al. 2017) at 800 K. The intensity cutoff, I cut (ν, 800 K), and strong line threshold, I strong (ν, 800K), are given as the dashed lines. For reference, each polyad region has been indicated.

  486 line lists. Consequently, all intensities of the RNT2017 487 strong line lists were converted to 296 K, then merged 488 into a global list of ∼27 million unique transitions. 489 The challenge of the second step is to convert the 490 super-line lists into "effective" lines that can be used 491 in line-by-line radiative transfer calculations. These are 492 much more flexible than temperature-specific line lists, 493 cross sections or k-correlation tables and make the fi-494 nal HITEMP line list more practical. However, the 495 RNT2017 strong line lists are provided at separate tem-496 peratures, meaning it is possible for a strong line at T 1 497 to be compressed into a super-line at T 2 . Hence, it is 498 also necessary to remove the contribution of the global 499 lines from each super-line list to avoid double counting 500 of individual transitions. 501 The global line list is calculated at all temperatures 502 given in Tab. (1), and the same temperature-dependent 503 thresholds from RNT2017 (I strong (ν, T ) and I cut (ν, T )) 504 are applied. Considering a transition at ν 1 with inten-

  506is part of the super-line list at T 1 . The line intensity I 1 507 will be included as part of the super-line intensity of the 508 nearest 0.005 cm -1 grid point to ν 1 . The super-line lists 509 are then reprocessed to remove the global strong line 510 contributions. In a small number of cases, the strong 511 line intensity at T 1 was greater than the correspond-512 ing super-line intensity at T 1 . This issue arises because 513 empirical corrections to the RNT2017 strong line lists 514 could not be disentangled from the empirical correc-515 tions applied to constituent transitions of each super-516 line, before they were compressed and the line informa-517 tion lost. It was deemed necessary to remove the line 518 intensity from the super-line lists, even when this inten-519 sity had to be removed from a neighboring super-line (to 520 avoid double counting of the strong line intensity). This 521 error is a consequence of attempts to reconstruct the 522 original RNT2017 line list (with 150 billion transitions) 523 prior to compression and can be completely avoided by 524 working with the original line list prior to compression. 525 We strongly recommend that for future investigations, 526 all line lists be retained, prior to the compression into 527 super-lines. 528 The reprocessed super-line lists are used to produce ef-529 fective lines that account for the continuum-like absorp-530 tion of CH 4 . These effective lines must have an effective 531 lower-state energy (allowing conversion of intensities be-532 tween temperatures) and can then be included with the 533 global line list above. From the intensity ratio of a line 534 as given in Eqn. (2), it is possible to determine the E of 535 a transition by comparing the line intensity at different 536 temperatures. Eqn. (

Figure 2 .

 2 Figure2. Effective lower-state energies (E eff ) have been calculated from the reprocessed super-lines ofRey et al. (2017). A sample grid point is shown for the pentad (a), octad (b), tetradecad regions (c), and between the icosad and triacontad regions (d). On the left panels, the reprocessed super-line intensity ratios are plotted for the sample grid points (ν in cm -1 ), using Eqn. (3). The retrieved values of E eff are provided (in cm -1 ) for a single line fit (dashed blue) and dual line fit, where the cold and hot component fits are solid green and red lines, respectively. The right panels display the reprocessed super-line intensities as a function of temperature for the same grid points, with the shaded region highlighting an upper/lower bound of a factor of two. In each case, the retrieved values of E eff have been used to calculate the intensity contribution from the single line fit (dashed blue) and dual line fits (green and red) at each temperature, with the combined dual line fit given as a dashed orange line.

Figure 3 .

 3 Figure 3. Comparisons of the RNT2017 lines lists against the more flexible line list from this work at (a) 500 K, (b) 1000 K, (c) 1500 K, and (d) 2000 K.In each panel, the shaded region indicates the spectral region that is beyond the RNT2017 line lists bounds at each temperature and are therefore not considered complete. These cross sections have been calculated using HAPI(Kochanov et al. 2016).

7615

  These tests have been performed with line wings set to 25 cm -1 for each line; the time may vary depending on the width of line wings taken into account.

a

  See Fig. (4) for number of lines per isotopologue. b Dependant on the temperature coverage given in Tab. (1). c Dependant on the spectral range given in Tab. (1).

Figure 4 .

 4 Figure 4. The number and coverage of all CH4 isotopologues included as part of this work towards a line list for HITEMP, with line intensities plotted at 296 K. The strong and effective lines of 12 CH4 have been indicated.

Fig. ( 5 Figure 5 .

 55 Fig. (5) details these comparisons for a 10 cm -1 por-783

Figure 6 .

 6 Figure 6. Absorption cross sections (in cm 2 /molecule) from HITRAN2016 (black) calculated at 296 K for 100 Torr of 12 CH4 and convolved to a resolution of 1.0 cm -1 , compared to this work (dashed red) and ExoMol 34to10 (blue). Each panel contains an overview of the icosad (a), triacontad (b), tetracontad (c) and pentacontakaipentad (d).

932Figure 7 .

 7 Figure 7. (a) An experimental transmission spectrum at 296 K from Hargreaves et al. (2015) near 3060 cm -1 , compared to this work. (b) An experimental absorption cross section at 295 K from Wong et al. (2019) near 6070 cm -1 , compared to this work. In both panels, the contribution of each isotopologue has been indicated, with only those isotopologues that have lines within the corresponding spectral regions (of this work) shown. light the quality of the line list from this work. Figs. (5) 937

974CH 4

 4 also begins to dissociate towards higher temper-975 atures, making it a challenge to record high-resolution 976 spectra without also observing contamination species. 977 Many more measurements are required throughout the 978 infrared in order to further increase assignments (and 979 thereby improve theoretical calculations), but also to 980 further validate methane line lists. 981 The comparisons presented in this work demonstrate 982 that this line list contains the most accurate line param-

Figure 8 .

 8 Figure 8. Comparisons with transmission spectra of Hargreaves et al. (2015) for the pentad region at 1173 K, with 60 Torr of 12 CH4 and a path length of 50 cm. Panel (a) displays the full band at a resolution of 0.015 cm -1 , with (b) showing a zoomed in feature. Panel (c) details the same as (a) but with a convolved resolution of 0.15 cm -1 , the same zoomed feature is given in (d). In all panels the observations are given in black, this work red, and ExoMol in blue.

Figure 9 .

 9 Figure 9. Comparisons with transmission spectra of Hargreaves et al. (2015) for the octad region at 1173 K, with 60 Torr of 12 CH4 and a path length of 50 cm. Panel (a) displays the full band at a resolution of 0.015 cm -1 , with (b) showing a zoomed in feature. Panel (c) details the same as (a) but with a convolved resolution of 0.15 cm -1 , the same zoomed feature is given in (d). In all panels the observations are given in black, this work red, and ExoMol in blue.

Figure 10 .

 10 Figure 10. Empirical absorption cross sections of CH4 from Wong et al. (2019) at 600 K, compared to calculations using this work and ExoMol 34to10 for the tetracontad (a), icosad (b), and triacontad (c) regions. In each case, the cross sections have been calculated for 100 Torr of CH4 and convolved to 2 cm -1 . The lower panels in each plot display the residuals.

  Figure 11.A diagram for the completeness of the 12 CH4 HITEMP line list. Each colored region identifies the wavenumber range of completeness (up to the limits indicated) when used at high temperature. For T < 800 K, the HITEMP line list is complete between ν = 0-13,400 cm -1 , whereas at 2000 K, it is complete between ν = 0-6600 cm -1 . Use of the 12 CH4 HITEMP line list at temperatures/wavenumbers within the grey region will lead to errors as a consequence of incompleteness and is not recommended.

Figure 12 .

 12 Figure 12. Empirical absorption cross sections of CH4 from Wong et al. (2019) at 1000 K, compared to calculations using this work, RNT2017 and ExoMol 34to10 for the triacontad region. In each case, the cross sections have been calculated for 100 Torr of CH4 and convolved to 2 cm -1 .

  makes the lines lists difficult to use and are often converted into k-correlation tables to speed up atmospheric calculations (Grimm & Heng 2015). However, these approaches do not provide the flexibility and practicality of a single line list. The line list created for this work has been able to 1107 combine the accuracy of HITRAN2016 with the com-1108 pleteness of RNT2017 to form a single line list of CH 4 1109 for high-temperature applications. This line list uses 1110 the familiar HITRAN/HITEMP format, making it com-1111 patible with existing radiative-transfer software. Fur-1112 thermore, the second generation of HAPI (Kochanov 1113 et al. 2016) is able to perform much faster radiative-1114 transfer calculations by using "just-in-time" compiled 1115 codes while still using the line-by-line Voigt profile cal-1116 culation with no interpolation. With this approach, 1117 the CH 4 line list consisting of ∼32 million lines shown 1118 in Fig. (4) can be processed in approximately 450 sec-1119 onds on a 12 core 2.6 GHz CPU. These speed improve-1120 ments make radiative-transfer calculations with the CH 4 1121 HITEMP line list practical for all users. 1122 Consequently, this work is currently the most accu-1123 rate, and practical line list of CH 4 for high-temperature 1124 applications. This work is freely available through the 1125 HITEMP section 6 of HITRANonline and is provided as 1126 a compressed ascii file along with a supporting Python 1127 tool to extract the line list for a specified spectral range.

  combined the separate 12 CH 4 line lists 1130 of RNT2017 with HITRAN2016 to provide the most ac-1131 curate line list of CH 4 for high-temperature applications. 1132 This work encompasses the 0-13,400 cm -1 spectral re-1133 gion and is sufficiently complete to be used in line-by-1134 line calculations up to 2000 K. As a result, this work has 1135 been included as part of the HITEMP database and is 1136 available via HITRANonline. 1137 To avoid incorrect conclusions during inter-1138 comparisons with previous work, it is necessary to briefly 1139 summarize the difference in terminology and distinctions 1140 with respect to each other. The data compression strat-1141 egy applied in the ab initio born databases of RNT2017 1142 (Rey et al. 2017) and then in ExoMol 34to10 (Yurchenko 1143 et al. 2017) consisted of a summation of weak line con-1144 tributions within small wavenumber intervals, with the 1145 integrated features defined as super-lines. Typically, 1146 these super-lines are provided on a regular grid and 1147 include numerous transitions with a large variety of 1148 lower-state energy levels. Strong transitions were then 1149 provided as separate lists. The super-line approach per-1150 mitted one to speed-up simulations of quasi-continuum 1151 absorption/emission cross-sections while the complete-1152 ness was maintained. However, this strategy leads to 1153 the loss of information on individual lower-state en-ergies. Thus a direct extrapolation of super-lines to other temperatures using standard conversion formulae is not possible. It was therefore necessary for RNT2017 and 34to10 to provide super-line lists for a range of temperatures.

  The primary advantage of this strategy is the effective lines can be included alongside the strong lines to form a single, practical line list of CH 4 . This strategy has enabled the production of a line list containing ∼32 million lines, provided in the HITRAN/HITEMP format and that is compatible with existing radiative-transfer software. This work is able to reproduce the RNT2017 1189 intensities up to 2000 K, but it must be stressed that 1190 using this work outside of the recommended tempera-1191 ture/wavenumber ranges may lead to issues of complete-1192 ness. The effective lines have to be used with caution 1193 as they are convenient for radiative-transfer simulations 1194 but should not be used for calculating energy levels as 1195 they do not represent transitions between real quantum 1196 states. 1197 The HITEMP line list of CH 4 produced for this work 1198 has been validated against available high-temperature 1199 measurements of CH 4 polyads up to 9000 cm -1 1200 (>1.11 µm), the triacontad region. Comparisons to al-1201 ternative CH 4 line lists demonstrate that this work is 1202 the most accurate at reproducing observations at high 1203 temperatures. 1204 The majority of the lines in our CH 4 line list are due 1205 to the principal isotopologue of methane; however the 1206 line list has been supplemented with lines of 13 CH 4 , 1207 12 CH 3 D and 13 CH 3 D taken from HITRAN2016. In the 1208 future, it will be productive to include improved high-1209 temperature line lists for these isotopologues, although 1210 they are not expected to be dominant absorption fea-1211 tures in CH 4 spectra.

  1212

  , is available via HITRANonline and is provided for users to work with the HITRAN and HITEMP line lists. The line-by-

		294	in laboratory experiments, with weaker features being
		295	provided from calculated line lists such as MeCaSDa (Ba
		296	et al. 2013). Campargue et al. (2012) provide ∼2500
		297	assignments between 6230-7920 cm -1 . For unassigned
		298	lines in this region, E has been determined for approx-
		299	imately half of these lines from spectra at 80 and 300 K,
	line nature and consistency between the HITRAN and	300	and remaining lines contain an estimated E . Between
	HITEMP databases mean that they are extremely flex-	301	7920-10,450 cm -1 , empirical line positions and intensi-
	ible when modeling a variety of environments. The	302	ties are provided without assignments and with a con-
	HITRAN and HITEMP parameters undergo rigorous	303	stant E (Brown 2005; Béguier et al. 2015a,b). Finally,
	validations against observations (Olsen et al. 2019; Har-	304	limited lower rotational assignments are given for lines
	greaves et al. 2019) and are regularly used in radia-	305	between 10,920-11,502 cm -1 (Benner et al. 2012).
	tive transfer codes such as LBLRTM (Clough et al.	306	For all spectral ranges, line-shape parameters have
	2005), NEMESIS (Irwin et al. 2008), the Reference	307	been provided from appropriate empirical observations.
	Forward Model (Dudhia 2017), RADIS (Pannier &	308	When these were unavailable, line-shape parameters
	Laux 2019) and the Planetary Spectrum Generator (Vil-	309	have been calculated using the algorithms described by
	lanueva et al. 2018).	310	Brown et al. (2013) and Lyulin et al. (2009).
	This article describes the addition of CH 4 to the	311	The main issue for the modeling of CH 4 absorp-
	HITEMP database, bringing the total number of	312	tion/emission at elevated temperature is to account for
	HITEMP molecules to eight. The aim of this line list is	313	the rapidly increasing contributions of hot bands, in
	to be accurate and complete, but at same time practi-	314	which a huge amount of excited rovibrational levels for
	cal (in terms of time required to calculate opacities) for	315	high-energy polyads (Tyuterev et al. 2013; Nikitin et al.
	high-temperature applications.	316	2015a; Rey et al. 2017) are involved. As mentioned,
		317	HITRAN2016 is unsuitable for high-temperature appli-
	2. LINE LISTS COMPARED IN THIS WORK	318	cations due to lack of completeness for hot bands and
	Over the last decade, there has been a significant in-	319	high-J transitions, but also because the assignment de-
	crease in the capability of theoretical calculations for	320	ficiencies and limited knowledge of lower-state energies,
	CH 4 spectroscopy at high temperatures (Rey et al. 2014a; Yurchenko & Tennyson 2014; Rey et al. 2017; Yurchenko et al. 2018), which coincides with the re-	321 322 323	E , for large spectral regions introduce errors at temper-atures beyond room-temperature. This is particularly true for the portion of the line list beyond 6230 cm -1
	quirement for sufficiently accurate high-temperature line	324	(i.e., < 1.3 µm).
	lists in order to characterize brown dwarfs and exoplan-
	ets (Tennyson & Yurchenko 2017; Fortney et al. 2019).
	This article broadly describes the three state-of-the-art	325	2.2. RNT2017 and TheoReTS calculated data
	line lists of CH 4 that have been used (and compared) in	326	For this study we use RNT2017, the latest high-
	this work.	327	temperature theoretical line list for 12 CH 4 constructed
	2.1. HITRAN2016	328 329	by Rey et al. (2017) and provided as part of the Reims-Tomsk collaboration via the TheoReTS data system
	In HITRAN2016 (Gordon et al. 2017), CH 4 (molecule	330	(Rey et al. 2016). RNT2017 contains significant im-
	6) contains parameters for four isotopologues: 12 CH 4 ,	331	provements with respect to the previous RNT2014 (Rey
	13 CH 4 , 12 CH 3 D and 13 CH 3 D. Line parameters are pro-	332	et al. 2014a) work, for which a good general agreement
	vided at 296 K and intensities are scaled for natu-	333	with experimental spectra up to 1200 K has been ob-
	ral abundances (0.988274, 0.011103, 6.15751×10 -4 and	334	served by Hargreaves et al. (2015) for the pentad (P 2 )
	6.91785×10 -6 , respectively). The partition function	335	and octad (P 3 ) regions (2.0-3.8 µm) . RNT2017 has re-
	from Gamache et al. (2017) is recommended when using	336	cently been validated against experimental observations
	HITRAN2016, and is also provided at HITRANonline.	337	up to 1000 K for the tetradecad (P 4 ), icosad (P 5 ) and
	For 12 CH 4 there are 313,943 transitions up to	338	triacontad (P 6 ) regions (1.11-1.85 µm) by Wong et al.
	11,502 cm -1 (P 8 ). Below 6230 cm -1 , there are both	339	(2019) at resolutions of 0.02, 0.2 and 2.0 cm -1 . In addi-
	upper-state and lower-state assignments for vibrational	340	tion, the region near 1.7 µm has also been validated to
	and rotational quanta for almost all transitions, however	341	accurate (±0.002 cm -1 ) observations at 1000 K by Ghy-
	there are only limited assignments beyond 6230 cm -1 .	342	sels et al. (2018) along with comparisons to MeCaSDa,
	The majority of assigned transitions have been validated		HITRAN2016 and ExoMol 10to10.

  -6 cm -1 /atm. In the context of

	670
	671	3.3. Comparison to RNT2017
	672	Fig. (3) contains a comparison of calculated absorp-	722	plete line lists. On a line-by-line basis, some of the
	673	tion cross sections using the total line list from this	723	strongest features will have an uncertainty much less
	674	work against line lists of RNT2017 at four tempera-	724	than this value, whereas above 6300 cm -1 the majority
	675 676 677 678 679	tures (combining strong and super-line components). These cross sections are calculated on a fine 0.001 cm -1 grid for the full 0 to 13,400 cm -1 spectral range with 100 Torr of 12 CH 4 , with calculations performed using HAPI (Kochanov et al. 2016). To make the comparisons	630 631 632 725 726 727 728 729	high-temperature applications, there is a large room for improvement for these line-shape parameters. For in-of hot lines have not been catalogued and the positions are not well known. For the strong line intensities, a cautious 20% uncertainty is suggested, although for strong lines within polyads up to P 5 the intensities can be much more accurate. The uncertainly of the effective
	680	appropriate, the same broadening and temperature de-	730	stance, the HITRAN default format allows only temper-line intensities and positions has not been assigned since
		pendence has been applied to the RNT2017 line lists as	731	the contributing lines are not observed and difficult to
			732	quantify. However the comparisons above indicate that
			733	these effective lines fall within the uncertainties of the
				strong lines.
			740
			741	3.4. HITRAN2016 replacements and isotopologue
			742	inclusions
				Where possible, owing to their high reliability,
			654
	702		655	appropriate calculations. Broadening parameters for
	703 704 705 706 707 708 709 710 711 712	tative of the accuracy of the current work because it is heavily weighted to the intensity sum of the dyad and pentad regions, but these intensity sums can be consid-ered an indicator of the uncertainty of the effective line intensities at each temperature. In total, this work requires fewer lines than the RNT2017 line lists, yet is significantly more flexible and able to reproduce the RNT2017 absorption. Cross-sections calculated for this work have been done so using the second generation of HAPI (Kochanov et al. 2016), which is available online 4 . The updates to HAPI mean	656 657 658 659 660 661 662 663 664 752 753 754 755 756 757 758 759 760	other gases will be added to the database in the near future as a response to the increasing amount of rele-vant experimental and theoretical studies. For instance, Gharib-Nezhad et al. (2019) recently measured broad-ening of CH 4 lines by H 2 over an extended range of temperatures. Finally, the HITRAN database has re-cently introduced advanced line-shape profiles (Wcis lo et al. 2016), due to the flexibility offered by the rela-tional database structure. These advanced line shapes can decrease residuals in terrestrial atmospheric spectra above 6230cm -1 , where assignments are limited or lack-ing, have very few HITRAN2016 replacements. In total, the 12 CH 4 line list from this work contains 81,245 lines replaced by HITRAN2016, which amounts to approxi-mately 50% of the HITRAN2016 12 CH 4 line list. The corresponding theoretical line lists of methane isotopologues is currently not sufficient for high-temperature applications, although progress is being made towards assignments (Rey et al. 2014b, 2015;

633

ature dependence for γ air , and using this temperature 634 dependence for γ self is only an approximate solution. 635 Furthermore, recent works that study the line shape 636 effects over a broad range of temperatures (Gamache & 637 Vispoel 2018; Stolarczyk et al. 2019) propose the use 638 of a double power law as opposed to a power law with 639 a single exponent. With that being said, Vispoel & 640 Lepère (2019) recently studied CH 4 lines broadened by 641 N 2 but did not observe a large discrepancy between a 642 single power law and double power law up to 700 K. 643 Another consideration for line broadening of CH 4 is by 644 "planetary" gases, including CO 2 , H 2 , He and H 2 O. As 645 previously discussed, HITRAN provides line broadening 646 by CO 2 , H 2 , He and H 2 O (Wilzewski et al. 2016; Tan 647 et al. 2019) for several gases. But for CH 4 , broaden-648 ing by H 2 O is the only additional perturber currently 649 available (Tan et al. 2019). To obtain water-broadened 650 parameters, Tan et al. (2019) recommend multiplying 651 γ air by a single scaling factor of 1.36 and multiplying 652 n air by a factor of 1.26. These factors can be applied 653 to the HITEMP line list from this work when doing 665 to the sub-percent level. While HITEMP line lists could 666 also benefit from their inclusion with respect to high-667 resolution combustion measurements, the main target of this work is astrophysical applications where such 669 accuracy on the line-shape parameters is not required. 681 was used for this work. In all cases, this work is able to 682 reproduce the absorption features seen when using the 683 RNT2017 line lists. For the lowest temperatures, the 684 first three window regions display differences when used 685 at high resolution: a consequence of a slight overesti-686 mation of the effective line strengths of the cold lines 687 at 296 K. However, these differences contribute an ex-688 tremely small amount to the total absorption. Tab. (2) 689 includes intensity sums from the single line list of this 690 work (ΣS TW (T )) calculated at the same temperatures 691 as RNT2017. These intensity sums have been compared 692 to those of RNT2017, given in Tab. (1). It can be seen 693 that the total intensity for the corresponding wavenum-694 ber limits agree with RNT2017 to within 2% up to 695 1100 K. This difference increases to a maximum of 6% 696 at 1500 K, before reducing towards 2000 K. It should be 697 noted that the strong-line sums (or effective/super-line 698 sums) are not directly comparable between this work 699 and RNT2017 because of differences between the num-700 ber of lines, and therefore intensity, included in each sum 701 (see Sect. (3.1)). This metric is not necessarily represen-734 For the remaining line parameters, an uncertainly 735 20% is given for γ self , γ air and n, and between ±0.01-736 0.001 cm -1 for δ. 737 To remain consistent with HITRAN, the strong and 738 effective line intensities have been scaled by 0.98827 to 739 account for the natural terrestrial abundance of 12 CH 4 . 743 HITRAN2016 parameters for 12 CH 4 have been used 744 in place of the RNT2017 values in this work (TW). 745 Matches were identified based on the criteria ν TW = 746 ν HIT2016 ± 0.01 cm -1 , E TW = E HIT2016 ± 0.01 cm -1 , 747 S TW = S HIT2016 ± 20 % and a consistent J between 748 line lists. The wavenumber criteria were relaxed to 749 ±1.0 cm -1 for transitions greater than 10,000 cm -1 due 750 to a reduced accuracy of the RNT2017 data. A re-751 quirement for matching J and E means that regions

Table 3 .

 3 A summary of the CH4 line list for HITEMP.

	Item	Details
	Isotopologues included a	12 CH4, 13 CH4, 12 CH3D,
		13 CH3D
	Total number of lines	31,880,412
	Proportion of effective lines 16.0%
	νmax(T ) b	13,400 cm -1 (>746 nm)
	Tmax(ν) c	2000 K
	T0	Parameters are provided
		at 296 K
	S(T0)	Intensities are scaled for
		natural abundances and
		given at T0 (see text)

https://hitran.org/hitemp/

ACKNOWLEDGMENTS

We would like to acknowledge J. A. Karns for contributions made towards Python data input routines. Update of the HITRAN and HITEMP databases was supported through the NASA Aura and PDART grants NNX17AI78G and NNX16AG51G, respectively. Support from the French ANR e-PYTHEAS project and from the D. Mendeleev program of Tomsk State University are also acknowledged.

that absorption cross sections calculated from a line list 714 4 https://github.com/hitranonline/hapi2 of ∼32 million can be processed in approximately 450 715 seconds on a 12 core 2.6 GHz CPU 5 .

716

For the strong lines that have been empirically cor-717 rected by RNT2017, an uncertainty of between ±0.01-718 0.001 cm -1 is obtained, which increases for lines with no 719 empirical corrections. This uncertainty corresponds to 720 the value reported by Rey et al. (2017) for the comsection calculated using the ExoMol 34to10 line list for 786 the same conditions is also included. In each of these 787 polyads regions, this work is able to replicate all features remarkably well with similar results seen across


Kok-Giap Haw

Simona Moldovan

Lingxue Tang

Zhengxing Qin

Qianrong Fang

Shilun Qiu

Valentin Valtchev
email: valentin.valtchev@ensicaen.fr

A sponge-like small pore zeolite with great accessibility to the micropores

Small pore zeolites, due to pore size constraint, have limited applications. Post-modification on the materials is therefore important to widen their use. A CHA-type zeolite (SSZ-13) with unique sponge-like structure is obtained by fluoride leaching. The development of macroporosity started from the crystal surfaces and continued progressively into the crystal by prolonging treatment. Nitrogen physisorption measurements showed an increase in micropore volume and specific surface area as a consequence of the dissolution of low crystalline part of the zeolite. The impact of etching on the accessiblity through the pore network, evaluated by means of breakthrough experiments of CO2/N2 and CO2/CH4 binary mixtures, showed an improved accessibility thanks to the interconnected macropores which shortens the diffusion pathlength. The set of experimental data shows the sponge-like SSZ-13 crystals retaining the intrinsic zeolitic properties but having an improve accessibility and crystallinity.

Introduction

SSZ-13 is an aluminosilicate zeolite with CHA-type topology. 1 The zeolite has attracted considerable attention thanks to its important application in selective catalytic reduction of nitrogen oxides (NO x ), [START_REF] Beale | Recent advances in automotive catalysis for NOx emission control by small-pore microporous materials[END_REF] and potential applications in methanol to olefins process (MTO), [START_REF] Yang | Recent Progress in Methanol-to-Olefins (MTO) Catalysts[END_REF] methylation of ethylene to propylene (ETP), [START_REF] Dai | Verifying the mechanism of the ethene-topropene conversion on zeolite H-SSZ-13[END_REF][START_REF] Jun | Conversion of Y into SSZ-13 zeolites and ethylene-topropylene reactions over the obtained SSZ-13 zeolites[END_REF] and CO 2 storage. [START_REF] Hudson | Unconventional, highly selective CO 2 adsorption in zeolite SSZ-13[END_REF][START_REF] Pham | Carbon dioxide and nitrogen adsorption on cation exchanged SSZ-13 zeolites[END_REF] The 3-dimensional pore system in SSZ-13 consists of large cavities (7.3 x 12 Å) that are connected by a small 8MR window (3.8 x 3.8 Å). These small pore openings allow only small molecules to diffuse through the pores with an excellent product shape selectivity. Thus the small pore zeolite, i.e. SSZ-13, is particularly appropriate for the reactions such as NH 3 -SCR and MTO in respect to the medium (ZSM-5) and large (Mordenite) pore zeolites.

8,9

However, it is very often that during a catalytic reaction the benefit of small pores in SSZ-13, i.e. the shape selectivity, is at the expense of its accessibility, thereby limiting the physical transport of reactants to the active sites. [START_REF] Tao | Mesopore-Modified Zeolites: Preparation, Characterization, and Applications[END_REF] Gao et al. found that the diffusion limitations play a significant role in the Cu-SSZ-13's catalytic performance. [START_REF] Gao | Structure-activity relationships in NH 3 -SCR over Cu-SSZ-13 as probed by reaction kinetics and EPR studies[END_REF] It is even worse in a micron size zeolite catalyst since only a minor fraction of the zeolites' active volume is being used in catalysis. Hierarchical zeolites which couple secondary network of meso-and macropores to the micropores could potentially solve the problem by enhancing the accessibility and thereby catalytic effectivity of the active sites in the zeolite crystal. It is thus important to enhance their mass transport properties. Secondary pores, i.e. meso-or macro-pores, in zeolites can be introduced either by post-synthesis methods or in-situ during crystal growth. By introducing secondary pores to these microporous materials, the transport of molecules through the catalytically active sites is enhanced and their catalytic performance is improved. [START_REF] Groen | Alkaline-mediated mesoporous mordenite zeolites for acidcatalyzed conversions[END_REF][START_REF] Mokrzycki | Properties of Desilicated ZSM-5, ZSM-12, MCM-22 and ZSM-12/MCM-41 Derivatives in Isomerization of α-Pinene[END_REF] Soft and hard templating routes are commonly used for in-situ introduction of secondary porosity. Depend on the target application, the secondary pores as well as the physicochemical properties of zeolites, could be customized during post-synthesis modification either by steaming and acid or alkaline leaching. Neutral pH etchant is employed when milder leaching is required. Neutral etching has the advantages of being environment friendly, easy to control and safe to handle. One of the most widely used neutral etchant is ammonium fluoride. 14 It was shown recently that post-treatment in ammonium fluoride solution allows the unbiased extraction of T-atoms from the pore system of a large pore FAU-type zeolite and generated additional microporous space thanks to the opening of sodalite cages. [START_REF] Qin | Opening the Cages of Faujasite-Type Zeolite[END_REF] The ammonium fluoride treatment was also successfully used to generate mesopores in medium (ZSM-5) and large (Mordenite) pore pentasil type zeolites.

Synthesis of SSZ-13

SSZ-13 was synthesized by interzeolite conversion of zeolite Beta in the presence of N,N,N-trimethyl-1-adamantammonium (TMAD) as a structure directing agent, reported elsewhere.

19

The molar composition of the initial gel was 0.236TMAD:0.5KOH:30H 2 O. Zeolite Beta was used as silicon and aluminum source and no additional Si or Al was added to the starting mixture. In a typical experiment, 1.7 g zeolite beta was added to a solution containing 0.60 g KOH, 1.94 g TMAD, 14.82 g H 2 O and stirred for 2 h at room temperature. The final molar composition of the gel after the addition of zeolite beta was 0.236TMAD:0.5KOH:1SiO 2 :0.

Characterizations

All materials were characterized by X-ray diffraction (XRD) using a PANalytical B.V. Empyrean powder diffractometer with Cu K α radiation (λ = 1.5418 Å) operating at 40 kV and 40 mA. The samples were scanned in the 2θ range of 4 -40° with a step size of 0.02°. Scanning electron microscopy (SEM) micrographs were taken on a JEOL JSM7400F microscope operated at 15 kV and JEOL JSM-7900F microscope operated at 1 kV. Before measurement, the sample was placed on a sample holder, coated with isopropyl alcohol-based carbon conductive adhesive. Transmission electron microscopy (TEM) analysis was performed on a JEOL JEM-2100F. Prior to the measurement, a diluted colloidal suspension of sample was sonicated for 5 min and then dropped on a carbon-filmcovered 300-mesh copper electron microscope grids and dried. Tomography study was performed using a JEOL JEM-2100F TEM equipped with a high-resolution objective lens pole piece at 200 kV. The tomography series were acquired between tilting angles of ±70° with a 2° Saxton scheme and the subsequent series alignments were performed with the IMOD software using Au nanoparticles (5-7 nm) as fiducial markers.

To resolve details at maximum resolution we used the SIRT algorithm implemented in fast software running on multicore computers, Tomo3D. The size of the TEM projections used for the reconstruction was 2k x 2k pixels. Nitrogen physisorption measurements were performed on a Micromeritics ASAP 2040 Prior to the measurements, the adsorbent was purged with a pure He (10 ml/min) through the column at 25 °C with a residence time of 5 h. Breakthrough experiments were performed at 25 °C with the flow of feed gas at 0.1 MPa of pressure and 8 ml/min of flow rate.

Results and discussion

Homemade rectangular SSZ-13 crystals with dimension between 2 to 5 μm are used as the precursor. The parent SSZ-13 exhibits high X-ray crystallinity without traces of concurrent phases (Fig. 1). An increase in crystallinity is observed on the samples after treatment. Such increase in the crystallinity is attributed to the removal of the amorphous and/or low crystalline part of the sample. The crystallinity of the samples is retained even after 90 min of treatment showing that the zeolite structure is stable during prolonged treatment. SEM micrograph of the parent SSZ-13 shows individual rectangular shape crystals with 2-5 micron in diameters and welldeveloped crystal faces (Fig. 2a). SEM inspection of the treated samples shows changes in surface morphology (Fig. 2b). Rectangular pores with straight edges are developed on the crystal faces. These pores are attributed to the dissolution of the nano crystalline domains. [START_REF] Qin | The mosaic structure of zeolites[END_REF] The interface between the coherent crystalline domains and between crystals represents a zone with high concentrations of structural defects, which is more vulnerable to chemical attack. Thus a mosaic structure of rectangular shaped pores is formed upon unbiased dissolution with the NH 4 F etching solution. An increase in the etching time leads to the formation of larger and deeper pores with retaining rectangular shape since the dissolution rate is similar in different crystallographic directions (Fig. S1 andS2). The size of these secondary pores varies between 50 and 500 nm. The final stage of this process is desegmentation of the zeolite crystals into much smaller particles (Fig. S1e,f). The treatment in the presence of hydrogen peroxide (Fig. 2f) is more intense than the sole use of NH 4 F (Fig. 2c). We interpret the higher concentration of pores as aggressive etching in the presence of H 2 O 2 with the potential to generate secondary pores from the smallest surface imperfections. The pores penetrate deeply in the bulk of the crystals, giving an open structure with spongelike appearance which differs from the smooth crystal surface of the pristine SSZ-13, as illustrated under SEM and TEM (Fig. 3). A careful inspection on the TEM tomography showed that these macropores are well connected (Fig. S3). A good connectivity between macropores is an important criterion in zeolite catalyst thanks to its exceptional accessibility. Nitrogen adsorption-desorption analysis of SSZ-13 samples showed Type I isotherms with a sharp uptake at a low pressure of P/P 0 < 0.05, which is characteristic of microporous materials (Fig. S4). The S BET and V mic of the NH 4 F treated samples are in good agreement to the XRD analyses (Table S1). Micropore volume is an important characteristic of zeolite-type materials which relates to their crystallinity. For instance, an increase in the S BET and V mic in the treated sample is attributed to the elimination of low crystalline phase from SSZ-13 which increases their crystallinity. The absence of hysteresis loop This journal is © The Royal Society of Chemistry 20xx

Please do not adjust margins Please do not adjust margins (Fig. S4) indicates the pores generated are macropores with pore diameter > 50 nm, which is in good agreement to TEM and SEM observations. This conclusion is supported by negligible increase of the total pore volume, V T which is due to the fact that N 2 adsorption analysis is not appropriate for studying macroporous materials. The Si/Al ratio of the parent SSZ-13 sample is 7.8 according to the ICP analysis (Table S1). Since NH 4 F treatment is non-chemically selective the elemental compositions of the treated samples (i.e. Si/Al ratio) are retained.

In comparison, both HCl and NaOH treated samples showed significant drop in crystallinity as shown in XRD patterns and nitrogen physisorption (Fig. S5 andS6). Both caustic and acid treatments are effective in desilication and dealumination as shown in EDX analysis (Table S2) but not effective in generating sponge-like SSZ-13 as shown in nitrogen physisorption (Table S2), and SEM investigation (Fig. S7).

The yield and weight loss of SSZ-13 after NH 4 F, H 2 O 2 , NaOH and HCl treatments at different time intervals are summarized in Table S3. The synergy between the two etchant makes the process of zeolite dissolution much more efficient. This will allow also to decrease the amount of NH 4 F employed and thus to make the treatment more economic and environmentally benign. Our interpretation is that hydrogen peroxide, as a powerful oxidant, prepare the surface T-atoms towards reaction with NH 4 F. The same reaction has been used for silicon etching in semiconductor industry. 14 The mechanism involving hydrogen peroxide and NH 4 F was reported elsewhere. 14 The performance of the sponge-like materials was analyzed and compared by means of breakthrough studies. The breakthrough curves from the separation between CO 2 /CH 4 (Fig. 5) and CO 2 /N 2 (Fig. S8) binary mixtures. All of the breakthrough curves intensities are higher than the loaded concentration of 25%, which is due to the strong adsorption of CO 2 when it replaces N 2 and CH 4 . The accessibility of the sponge-like material is better than the pristine SSZ-13 due to the presence of interconnected macropores which shortens the diffusion pathlength. The increase in accessibility is at the expense of separation performance. This result reveals the advantage of sponge-like structure which provides rapid access and effective use of intrinsic zeolite micropores.

Highest

Conclusions

In summary, we have synthesized SSZ-13 material in spongelike structure with great accessibility properties without compromising their intrinsic properties, which would bring substantial practical advantages. The sponge-like structure will bring benefits of rapid access and effective use of zeolite micropores and increasing the effective surface areas. This work demonstrated a simple and efficient route to develop sponge-like porous structure on a small pore zeolite, which could be extended to other zeolite materials.

032Al 2 O 3 :

 3 30H 2 O. After homogenization, the mixture was transferred to Teflonlined stainless-steel autoclave. The hydrothermal synthesis was carried out at 150 o C for 24 h. After crystallization step, the solid product was collected by vacuum filtration and washed thoroughly with deionized water until a near neutral pH and dried overnight at 60 o C. The as-synthesized SSZ-13 sample was calcined at 600 °C (at 100 °C•h -1 ramp) in air flow for 10 h. Preparation of sponge-like SSZ-13 NH 4 F solution (40 wt%, aq.) was prepared by dissolving 4.0 g of NH 4 F in 3.0 g of H 2 O and 3.0 g of H 2 O 2 . 0.3 g of SSZ-13 was then added into the H 2 O 2 -NH 4 F solution. The zeolite treatment was performed at 20 o C under ultrasonic radiation for t = 15, 30, 60, 90, 120, and 180 min respectively. The heat from the radiation was compensated by adding ice into the ultrasonic bath and the temperature was maintained at 20 ± 2 o C. The samples were rinsed thoroughly with deionized water under vacuum filtration and dried overnight at 60 o C. The effect of H 2 O 2 -assisted etching was compared with the sole NH 4 F treatment under the same conditions. The samples obtained were denoted as H 2 O 2 -NH 4 F USt and NH 4 F USt, where t is the time in minutes.For comparison purposes, etching of SSZ-13 in an acidic etchant HCl (0.5 M), and an alkaline etchant NaOH (0.5 M), was performed. 0.3 g of SSZ-13 was added into 10 g of the corresponding acid or alkaline solutions. Etching was performed at 20 o C under ultrasonic (US) radiation for t = 30 and 60 min. The heat from the radiation was compensated by addition of ice into the ultrasonic bath and the temperature was maintained at 20 ± 2 o C. The samples were rinsed thoroughly with deionized water under vacuum filtration and dried overnight at 60 o C. The samples obtained were denoted as HCl USt and NaOH USt respectively, where t is the time in minutes.

2 g - 1 )

 21 Please do not adjust marginsPlease do not adjust margins surface area analyzer. The calcined samples were analyzed after degassing at 300 °C. The specific surface areas (S BET / m were obtained using the Brunauer-Emmett-Teller (BET) equation. The microporous volume (V micro / cm 3 g -1 ) was obtained from t-plot based on the Harkins-Jura equation. The pore size distribution and total pore volume (V T / cm 3 g -1 ) was obtained from the desorption branch using the Barrett-Joyner-Halenda (BJH) algorithm assuming cylindrical pores. NMR measurements were done with 4-OD mm zirconia rotors and with a spinning speed of 12 kHz. Magic angle spinning nuclear magnetic resonance (MAS NMR) 27 Al and 29 Si spectra were recorded on a Bruker Avance 500 spectrometer operating at 130.3 MHz, with a π/12 pulse and a recycle delay of 1 s. Al(NO 3 ) 3 1M was used as a reference for 27 Al. Elemental analysis of the sample was measured by Energy Dispersive X-Ray Fluorescence Spectrometer (EDX) operated under TEM (JEM-2100F). Breakthrough experiment Na-exchanged SSZ-13 samples were used in breakthrough experiments. The samples were ion exchanged with 0.3 M NaNO 3 solution (solid:liquid = 1:100, weight ratio) twice. Breakthrough experiments were performed with binary mixtures of CH 4 /CO 2 /He (25/25/50) and N 2 /CO 2 /He (25/25/50), respectively. The gas mixtures were equilibrated at 25 o C with a residence time of 5 h before the measurement. A 4.32 mm diameter stainless steel adsorption column with a length of 10 cm packed with approximately 0.8 g of sample.

Fig. 1

 1 Fig. 1 XRD patterns of parent and sponge-like SSZ-13 samples.

Fig. 2

 2 Fig. 2 SEM micrographs of parent SSZ-13 (a), general and magnified view of NH4F US60 (b, c) and H2O2-NH4F US60 (d, e, f). Inset showed the region magnified in c and f.

Fig. 3

 3 Fig. 3 SEM (a,b,c) and TEM (d,e,f) micrographs showing the mosaic structure (a,d,e) and sponge-like structure (b,c,f) of sample H2O2-NH4F US90.

4 .

 4 weight loss is recorded on the NH 4 F treated sample due to the formation of sponge-like structure. In sole use of hydrogen peroxide, 8% of weight loss is recorded. Since no change of morphology is detected (vide supra SEM, TEM and N 2 physisorption), such weight loss might attribute to the extraction of non-crystalline part of the sample. The same observations applied for NaOH and HCl treated samples since no changes in morphologies are observed. 27 Al and 29 Si MAS NMR analysis provides valuable information for the short range order in zeolite frameworks. The NMR spectra of the parent and the treated SSZ-13 samples are presented in Fig. There is only single well-resolved signal at 48 ppm in the 27 Al MAS NMR spectra (Fig. 4a) of the parent and treated samples attributed to tetrahedrally coordinated

Fig. 4 27

 4 Fig. 4 27 Al (a) and 29 Si (b) MAS NMR spectra of the as-synthesized parent and sponge-like SSZ-13 samples.

Fig. 5

 5 Fig. 5 Breakthrough curves of pristine and sponge-like SSZ-13 (H2O2-NH4F US60) in the presence of CH4/CO2 gas mixtures at 25 o C.

  , which then reacts with NH 4 F and forms soluble hexafluorosilicates. The exposed silicon is consecutively oxidized to silicon oxide intermediate and reacts with NH 4 F. Thus, the etching proceeds spontaneously.

Briefly, the reaction begins with the oxidizing agent (i.e. H 2 O 2 ) reacts with silicon and forms silicon oxide Please do not adjust margins Please do not adjust margins intermediate

This journal is © The Royal Society of Chemistry 20xx J. Name., 2013, 00, 1-3 | 3

This journal is © The Royal Society of Chemistry 20xx J. Name., 2013, 00, 1-3 | 5

Acknowledgements

V.V. and Q.F. thank the National Natural Science Foundation of China (21571079, 21621001, 21390394, 21571076 and 21571078) for the financial support. S.M. and V.V. acknowledge the financial support from the CARNOT project ESP 3DNANOZET 5281.

Notes and references

Conflicts of interest

There are no conflicts to declare.