ANTI-CORRELATED NOISE FILTER

Abstract

An imaging system (100) includes an anti-correlated noise filter (120), which jointly filters noise from a first portion (116) and a second portion (118), and the first portion (116) and the second portion (118) include anti-correlated noise.

Claims

1. An imaging system comprising: an anti-correlated noise filter configured to jointly filter noise from a first portion and a second portion, and the first portion and the second portion include anti-correlated noise.

2. The imaging system according to claim 1, wherein the first portion and the second portion include spectral CT data from a basis decomposition of at least one of: projection data of an object or subject; or image data of an object or subject.

3. The imaging system according to claim 1, wherein the anti-correlated noise filter jointly filtered noise is suppressed based on at least one of: a weighted difference between initially combined data of the first portion and the second portions, and a sum of a filtered first portion and a filtered second portion; and the filtered first portion and the filtered second portion selected to minimize noise in a spectral monochromatic image which includes a weighted combination of the filtered first portion and the filtered second portion.

4. The imaging system according to claim 1, wherein the anti-correlated noise filter is further configured to filter noise according to a function defined by: $\left(\hat{p},\hat{s}\right)=\arg .Math..Math.{\min }_{\left(p,s\right)}.Math.R\left(p\right)+R\left(s\right)+\frac{1}{2}.Math.\int {{\lambda }^u\left(p+s-u^0\right)}^2+\frac{1}{2}.Math.\int {{\lambda }^p\left(p-p^0\right)}^2+\frac{1}{2}.Math.\int {{\lambda }^s\left(s-s^0\right)}^2,$ where R(p) and R(s) are roughness penalties for p and s, respectively, u.sup.0 is an image volume where the correlated noise maximally cancels out with the initially decomposed portions, p.sup.0 and s.sup.0, e.g., u.sup.0=p.sup.0+s.sup.0, p and s are the filtered image volumes, and λ.sup.u, λ.sup.p and λ.sup.s are weights.

5. The imaging system according to claim 4, wherein the function is implemented by: $\begin{array}{c}{s}_{i,j,k}^{n+1}=\frac{{\lambda }_{i,j,k}^u\left(u_{i,j,k}^0-p_{i,j,k}^n\right)+{\lambda }_{i,j,k}^s.Math.s_{i,j,k}^0+{\mathrm{\delta \Sigma }}_D.Math.{\sigma }_D^n.Math.s_D^n}{{\lambda }_{i,j,k}^u+{\lambda }_{i,j,k}^s+{\mathrm{\delta \Sigma }}_D.Math.{\sigma }_D^n}.Math..Math.\mathrm{and}\\ p_{i,j,k}^{n+1}=\frac{{\lambda }_{i,j,k}^u\left(u_{i,j,k}^0-s_{i,j,k}^n\right)+{\lambda }_{i,j,k}^p.Math.s_{i,j,k}^0+{\mathrm{\delta \Sigma }}_D.Math.{\varphi }_D^n.Math.s_D^n}{{\lambda }_{i,j,k}^u+{\lambda }_{i,j,k}^p+{\mathrm{\delta \Sigma }}_D.Math.{\varphi }_D^n}\end{array}$ where λ.sup.u, λ.sup.p, λ.sup.s, σ.sub.D are weights, p.sup.0 and s.sup.0 are the first and second decomposed portions, p.sup.n and s.sup.n are current values of the n.sup.th iteration of p.sup.0 and s.sup.0, p.sup.n+1 and s.sup.n+1 are a next iteration filter first and second portion, D includes a set of each orthogonal three dimensional direction {E(ast), W(est), S(outh), N(orth), U(p), and (d)O(wn)}, and i,j,k represent a current voxel in the image volume or a position in projection space volume, and δ is a Huber parameter.

6. The imaging system according to claim 1, wherein the anti-correlated noise filter is further configured to filter noise according to the function defined by: $\left(\hat{s},\hat{p}\right)=\arg .Math..Math.\underset{\left(s,p\right)}{\min }.Math.\alpha .Math..Math.R\left(s\right)+\left(1-\alpha \right).Math.R\left(p\right)$ subject to the constraints that (s.t.) 1. s and p are obtained by removing negatively correlated estimated noise from s.sup.0 and p.sup.0, respectively; 2. {circumflex over (m)} monochromatic image is unchanged; and 3. image frequencies outside band frequencies are unchanged, where R(p) and R(s) are roughness penalties or regularization terms for p and s, respectively, {circumflex over (m)} is an energy level parameter in keV unit, and α is an algorithm control parameter.

7. The imaging system according to claim 6, wherein the function is implemented by detecting a spectral virtual monochromatic image, {circumflex over (m)}, in which the anti-correlated noise is minimized, and the generating a new s and p based on the detected monochromatic image, and {circumflex over (m)} is defined by: $\hat{m}=\arg .Math..Math.\underset{m}{\min }.Math.\frac{R\left(c_s\left(m\right).Math.s+c_P\left(m\right).Math.p\right)}{c_s\left(m\right).Math.R\left(s\right)+c_p\left(m\right).Math.R\left(p\right)},$ where c.sub.s(m) and c.sub.p (m) are the coefficients of the first decomposed portion, s and the second decomposed portion, p, respectively, that enable the algorithm to obtain the monochromic image {circumflex over (m)} for energy {circumflex over (m)} keV.

8. The imaging system according to claim 6, wherein the function is implemented by detecting a spectral virtual monochromatic image, {circumflex over (m)} by defining a selection region of the combined spectral data using a predetermined threshold value, a local standard deviation calculated for a neighborhood of size ne for the combined spectral data, a set, q, of locations is created of an r smallest local standard deviations located in the selection region and {circumflex over (m)} is defined by: $\hat{m}=\arg .Math..Math.\underset{m}{\min }.Math.\Sigma .Math..Math.\mathrm{localstddev}\left(c_s\left(m\right).Math.s+c_p\left(m\right).Math.p,q,\mathrm{ne}\right),$ where the local standard deviation is calculated only over the set q and ne specify the neighborhood of the local standard deviation.

9. The imaging system according to claim 1, wherein the anti-correlated noise filter is further configured to filter noise according to the function defined by:
{circumflex over (p)}=p.sup.0Â and ŝ=ŝ=s.sup.0+Â, where $\hat{A}={\hat{A}}_{L_1}+\mathrm{ScaleUp}\left({\hat{A}}_{L_2}^d-\mathrm{ScaleDown}\left({\hat{A}}_{L_1},d\right),d\right),.Math.{\hat{A}}_{L_1}=\underset{A}{\mathrm{argmin}.Math.}.Math.R\left(s^0+A\right)+R\left(p^0-A\right)+\int {\lambda }_1.Math.h_{\delta }\left(A\right),.Math.{\hat{A}}_{L_2}^d=\underset{A}{\mathrm{argmin}.Math.}.Math.R\left(s_d^0+A\right)+R\left(p_d^0-A\right)+\int {\lambda }_2.Math.A^2,$ s.sub.d.sup.0=ScaleDown (s.sup.0, d), p.sub.d.sup.0=ScaleDown (p.sup.0, d), d is a scale parameter, R(.Math.) is a roughness penalty or regularization term, λ.sub.1 and λ.sub.2 are weights, Â is the estimated anti-correlated noise image, A is a prior estimate of the anti-correlated noise image, h.sub.δ(A)=δ.sup.2(√{square root over (1+(A/δ).sup.2)}−1)is the pseudo-Huber penalty function and δ is the pseudo-Huber parameter.

10. The imaging system according to claim 1, wherein the anti-correlated noise filter is further configured to filter noise according to the function defined by:
{circumflex over (p)}=p.sup.0−Â and ŝ=s.sup.0+Â, where $\hat{A}={\hat{A}}_{L_1}+\mathrm{ScaleUp}\left({\hat{A}}_{L_2}^d,d\right),.Math.{\hat{A}}_{L_1}=\underset{A}{\mathrm{argmin}.Math.}.Math.R\left(s^0+A\right)+R\left(p^0-A\right)+\int \frac{{\lambda }_1}{n}.Math.h_{\delta }\left(A\right),\mathrm{and}.Math..Math.{\hat{A}}_{L_2}^d=\underset{A}{\mathrm{argmin}.Math.}.Math.R\left(s_d^0+A\right)+R\left(p_d^0-A\right)+\int \left(\frac{{\lambda }_2}{n^2}+{\lambda }_3.Math.F\right).Math.A^2,$ where s.sub.d.sup.0=ScaleDown (s.sup.0+Â.sub.L.sub.1, d), p.sub.d.sup.0=ScaleDown (p.sup.0−Â.sub.L.sub.1, d), d is a scale parameter, R(.Math.) is a roughness penalty or regularization term, λ.sub.1, λ.sub.2 and λ.sub.3 are weights, Â is the estimated anti-correlated noise image, A is a prior estimate of the anti-correlated noise image, h.sub.δ(A)=δ.sup.2(√{square root over (1+(A/δ).sup.2)}−1) is the pseudo-Huber penalty function, δ is the pseudo-Huber parameter, n is an estimated noise map, and $F=\mathrm{ScaleDown}\left(\frac{\sigma \left(s_d^0+{\hat{A}}_{L_1}\right).Math.\sigma \left(p_d^0-{\hat{A}}_{L_1}\right)}{1+{\sigma \left({\hat{A}}_{L_1}\right)}^2}\right),$ where σ(x) is the local standard deviation of the image .sup.x .

11. The imaging system according to claim 1, wherein the first portion and the second portion are basis pairs and include at least one of: a photoelectric absorption component and a Compton-scatter component; a water component and an Iodine component; a water component and a Calcium component; or an acetal homopolymer resin component and a tin components.

12. The imaging system according to claim 1, wherein the anti-correlated noise filter is further configured to iteratively filter noise from the first portion and the second portion until a stopping criteria is reached.

13. The imaging system according to claim 1, wherein the anti-correlated noise filter is further configured to filter separately the first portion and the second portion with a Structure Propagation (SP) filter prior to jointly filtering anti-correlated noise from the SP filtered first portion and the SP filtered second portion.

14. A method of filtering image data, comprising: jointly filtering noise from a first portion and a second portion, and the first portion and the second portion include anti-correlated noise.

15. The method according to claim 14, wherein the first portion and the second portion are formed from a basis decomposition of spectral CT data, which includes at least one of: projection data of an object or subject; or imaging data of an object or subject.

16. The method according to claim 14, wherein jointly filtering includes at least one of: weighting a difference between initially combined data of the first portion and the second portions, and a sum of a filtered first portion and a filtered second portion; and selecting the filtered first portion and the filtered second portion to minimize noise in a spectral monochromatic image which includes a weighted combination of the filtered first portion and the filtered second portion.

17. The method according to claim 14, wherein filtering noise filtered according to the function defined by: $\left(\hat{p},\hat{s}\right)={\mathrm{argmin}}_{\left(p,s\right)}.Math.R\left(p\right)+R\left(s\right)+\frac{1}{2}.Math.\int {\lambda }^u\left(p+s-u^0\right)+\frac{1}{2}.Math.\int {{\lambda }^p\left(p-p^0\right)}^2+\frac{1}{2}.Math.\int {{\lambda }^s\left(s-s^0\right)}^2,$ where R(p) and R(s) are roughness penalties for p and s, respectively, u.sup.0 is an image volume where the correlated noise maximally cancels out with the initially decomposed portions, p.sup.0 and s.sup.0, e.g., u.sup.0=p.sup.0+s.sup.0, p and s are the filtered image volumes, and λ.sup.u, λ.sup.p and λ.sup.2 are weights.

18. The method according to claim 14, wherein the function is implemented by: $s_{i,j,k}^{n+1}=\frac{{\lambda }_{i,j,k}^u\left(u_{i,j,k}^0-p_{i,j,k}^n\right)+{\lambda }_{i,j,k}^s.Math.s_{i,j,k}^0+\delta .Math.{.Math.}_D.Math.{\sigma }_D^n.Math.s_D^n}{{\lambda }_{i,j,k}^u+{\lambda }_{i,j,k}^s+\delta .Math.{.Math.}_D.Math.{\sigma }_D^n}$ $\mathrm{and}$ $p_{i,j,k}^{n+1}=\frac{{\lambda }_{i,j,k}^u\left(u_{i,j,k}^0-s_{i,j,k}^n\right)+{\lambda }_{i,j,k.Math.}^p.Math.s_{i,j,k}^0+\delta .Math..Math.{.Math.}_D.Math.{\varphi }_D^n.Math.s_D^n}{{\lambda }_{i,j,k}^u+{\lambda }_{i,j,k}^p+\delta .Math..Math.{.Math.}_D.Math.{\varphi }_D^n}$ where λ.sup.u, λ.sup.p, λ.sup.s, σ.sub.D are weights, p.sup.0and s.sup.0 are the first and second decomposed portions, p.sup.n and s.sup.n are current values of the n.sup.th iteration of p.sup.0 and s.sup.0, p.sup.n+1 and s.sup.n+1 are a next iteration filter first and second portion, D includes a set of each orthogonal three dimensional direction {E(ast), W(est), S(outh), N(orth), U(p), and (d)O(wn)}, and i,j,k represent a current pixel, and δ is a Huber parameter.

19. The method according to claim 14, wherein filtering noise is filtered according to the function defined by: $\left(\hat{s},\hat{p}\right)=\arg .Math..Math.\underset{\left(s,p\right)}{\min }.Math.\alpha .Math..Math.R\left(s\right)+\left(1-\alpha \right).Math.R\left(p\right)$ subject to the constraints that (s.t) 1. s and p are obtained by removing negatively correlated estimated noise from s.sup.0 and p.sup.0, respectively; 2. {circumflex over (m)} monochromatic image is unchanged; and 3. image frequencies outside band frequencies are unchanged, where R(p) and R(s) are roughness penalties or regularization terms for p and s, respectively, {circumflex over (m)} is an energy level parameter in keV unit, and α is an algorithm control parameter.

20. The method according to claim 19, wherein the function is implemented by detecting a spectral virtual monochromatic image, {circumflex over (m)}, in which the anti-correlated noise is minimized, and the generating a new s and p based on the detected monochromatic image, and {circumflex over (m)} is defined by: $\hat{m}=\arg .Math..Math.\underset{m}{\min }.Math.\frac{R\left(c_s\left(m\right).Math.s+c_p\left(m\right).Math.p\right)}{c_s\left(m\right).Math.R\left(s\right)+c_p\left(m\right).Math.R\left(p\right)},$ where c.sub.s(m) and c.sub.p (m) are the coefficients of the first decomposed portion, s and the second decomposed portion, p, respectively, that enable the algorithm to obtain the monochromic image {circumflex over (m)} for energy {circumflex over (m)} M in keV.

21. The imaging system according to claim 19, wherein the function is implemented by detecting a spectral virtual monochromatic image, {circumflex over (m)} by defining a selection region of the combined spectral data using a predetermined threshold value, such as −200 HU. A local standard deviation is calculated for a neighborhood of size ne for the combined spectral data. A set, q of locations is created of the r smallest local standard deviations located in the selection region and an example of {circumflex over (m)} is defined by: $\hat{m}=\arg .Math..Math.\underset{m}{\min }.Math..Math.\mathrm{localstddev}\left(c_s\left(m\right).Math.s+c_p\left(m\right),p,q,\mathrm{ne}\right),$ where the local standard deviation is calculated only over the set q and ne specify the neighbourhood of the local standard deviation.

22. The method according to claim 14, wherein filtering noise is filtered according to the function defined by:
{circumflex over (p)}=p.sup.0−Â and ŝ=s.sup.0+Â, where $\hat{A}={\hat{A}}_{L_1}+\mathrm{ScaleUp}\left({\hat{A}}_{L_2}^d-\mathrm{ScaleDown}\left({\hat{A}}_{L_1},d\right),d\right),.Math.{\hat{A}}_{L_1}=\underset{A}{\mathrm{argmin}}.Math.R\left(s^0+A\right)+R\left(p^0-A\right)+\int {\lambda }_1.Math.h_{\delta }\left(A\right),.Math.{\hat{A}}_{L_2}^d=\underset{A}{\mathrm{argmin}}.Math.R\left(s_d^0+A\right)+R\left(p_d^0-A\right)+\int {\lambda }_2.Math.A^2,$ where s.sub.d.sup.0=ScaleDown (s.sup.0, d), p.sub.d.sup.0=ScaleDown (p.sup.0, d), d is a scale parameter, R(.Math.) is a roughness penalty or regularization term, λ.sub.1 and λ.sub.2 are weights, Â is the estimated anti-correlated noise image, A is a prior estimate of the anti-correlated noise image, h.sub.δ(A)=δ.sup.2(√{square root over (1+(A/δ).sup.2)}−1) is the pseudo-Huber penalty function and δ is the pseudo-Huber parameter.

23. The method according to claim 14, wherein filtering noise is filtered according to the function defined by:
{circumflex over (p)}=p.sup.0−Â and ŝ=s.sup.0+Â, where $\hat{A}={\hat{A}}_{L_1}+\mathrm{ScaleUp}\left({\hat{A}}_{L_2}^d,d\right),.Math.{\hat{A}}_{L_1}=\underset{A}{\mathrm{argmin}}.Math.R\left(s^0+A\right)+R\left(p^0-A\right)+\int {\lambda }_1.Math.h_{\delta }\left(A\right),\mathrm{and}$ ${\hat{A}}_{L_2}^d=\underset{A}{\mathrm{argmin}}.Math.R\left(s_d^0+A\right)+R\left(p_d^0-A\right)+\int \left(\frac{{\lambda }_2}{n^2}+{\lambda }_3.Math.F\right).Math.A^2,$ where s.sub.d.sup.0=ScaleDown (s.sup.0+Â.sub.L.sub.1, d), p.sub.d.sup.0=ScaleDown (p.sup.0−Â.sub.L.sub.1, d), d is a scale parameter, R(.Math.) is a roughness penalty or regularization term, λ.sub.1, λ.sub.2 and λ.sub.3 are weights, Â is the estimated anti-correlated noise image, A is a prior estimate of the anti-correlated noise image, h.sub.δ(A)=δ.sup.2(√{square root over (1+(A/δ).sup.2)}−1) is the pseudo-Huber penalty function, δ is the pseudo-Huber parameter, n is an estimated noise map, and $F=\mathrm{ScaleDown}\left(\frac{\sigma \left(s_d^0+{\hat{A}}_{L_1}\right).Math.\sigma \left(p_d^0-{\hat{A}}_{L_1}\right)}{1+{\sigma \left({\hat{A}}_{L_1}\right)}^2}\right),$ where σ(.sup.x) is the local standard deviation of the image .sup.x .

24. The method according to claim 14, wherein the first portion and the second portion are formed from basis decomposition of CT spectral imaging data and decomposed into basis pairs which include at least one of: a photoelectric absorption component and a Compton-scatter component; a water component and an Iodine component; a water component and a Calcium component; or an acetal homopolymer resin component and a tin components.

25. The method according to claim 14, wherein filtering further includes: filtering separately the first portion and the second portion using a Structure Propagation (SP) filter prior to jointly filtering anti-correlated noise from the SP filtered first portion and the SP filtered second portion.

26. A non-transitory computer readable storage medium encoded with computer readable instructions, which, when executed by a processor, causes the processor to: jointly filter noise from a first portion and a second portion, and the first portion and the second portion include anti-correlated noise, and the filter operates iteratively according to at least one of the following functions: $\left(\hat{p},\hat{s}\right)={\mathrm{argmin}}_{\left(p,s\right)}.Math.R\left(p\right)+R\left(s\right)+\frac{1}{2}.Math.\int {{\lambda }^u\left(p+s-u^0\right)}^2+\frac{1}{2}.Math.\int {{\lambda }^p\left(p-p^0\right)}^2+\frac{1}{2}.Math.\int {{\lambda }^s\left(s-s^0\right)}^2,$ where R(p) and R(s) are roughness penalties for p and s, respectively, u.sup.0 is an image volume where the correlated noise maximally cancels out with the initially decomposed portions, p.sup.0 and s.sup.0, e.g., u.sup.0=p.sup.0+s.sup.0, p and s are the filtered image volumes, and λ.sup.u, λ.sup.p and λ.sup.s are weights; or $\left(\hat{s},\hat{p}\right)=\arg .Math..Math.\underset{\left(s,p\right)}{\min }.Math.\alpha .Math..Math.R\left(s\right)+\left(1-\alpha \right).Math.R\left(p\right)$ subject to the constraints that (s.t.) 1. s and p are obtained by removing negatively correlated estimated noise from s.sup.0 and p.sup.0, respectively; 2. {circumflex over (m)} monochromatic image is unchanged; and 3. image frequencies outside band frequencies are unchanged, where R(p) and R(s) are roughness penalties or regularization terms for p and s, respectively, {circumflex over (m)} is an energy level parameter in keV unit, and α is an algorithm control parameter;
{circumflex over (p)}=p.sup.0−Â and ŝ=s.sup.0+Â, where $\hat{A}={\hat{A}}_{L_1}+\mathrm{ScaleUp}\left({\hat{A}}_{L_2}^d-\mathrm{ScaleDown}\left({\hat{A}}_{L_1},d\right),d\right),.Math.{\hat{A}}_{L_1}=\underset{A}{\mathrm{argmin}}.Math.R\left(s^0+A\right)+R\left(p^0-A\right)+\int {\lambda }_1.Math.h_{\delta }\left(A\right),.Math.{\hat{A}}_{L_2}^d=\underset{A}{\mathrm{argmin}}.Math.R\left(s_d^0+A\right)+R\left(p_d^0-A\right)+\int {\lambda }_2.Math.A^2,$ where s.sub.d.sup.0=ScaleDown (s.sup.0, d), p.sub.d.sup.0=ScaleDown (p.sup.0, d), d is a scale parameter, R(.Math.) is a roughness penalty or regularization term, λ.sub.1 and λ.sub.2 are weights, Â is the estimated anti-correlated noise image, A is a prior estimate of the anti-correlated noise image, h.sub.δ(A)=δ.sup.2(√{square root over (1+(A/δ).sup.2)}−1) is the pseudo-Huber penalty function and δ is the pseudo-Huber parameter; or
{circumflex over (p)}=p.sup.0−Â and ŝ=s.sup.0+Â, where $\hat{A}={\hat{A}}_{L_1}+\mathrm{ScaleUp}\left({\hat{A}}_{L_2}^d,d\right),.Math.{\hat{A}}_{L_1}=\underset{A}{\mathrm{argmin}}.Math.R\left(s^0+A\right)+R\left(p^0-A\right)+\int \frac{{\lambda }_1}{n}.Math.h_{\delta }\left(A\right),\mathrm{and}$ ${\hat{A}}_{L_2}^d=\underset{A}{\mathrm{argmin}}.Math.R\left(s_d^0+A\right)+R\left(p_d^0-A\right)+\int \left(\frac{{\lambda }_2}{n^2}+{\lambda }_3.Math.F\right).Math.A^2,$ where s.sub.d.sup.0=ScaleDown (s.sup.0+Â.sub.L.sub.1, d), p.sub.d.sup.0=ScaleDown (p.sup.0−Â.sub.L.sub.1, d), d is a scale parameter, R(.Math.) is a roughness penalty or regularization term, λ.sub.1, λ.sub.2 and λ.sub.3 are weights, Â is the estimated anti-correlated noise image, A is a prior estimate of the anti-correlated noise image, h.sub.δ(A)=δ.sup.2(√{square root over (1+(A/δ).sup.2)}−1) is the pseudo-Huber penalty function, δ is the pseudo-Huber parameter, n is an estimated noise map, and $F=\mathrm{ScaleDown}\left(\frac{\sigma \left(s_d^0+{\hat{A}}_{L_1}\right).Math.\sigma \left(p_d^0-{\hat{A}}_{L_1}\right)}{1+{\sigma \left({\hat{A}}_{L_1}\right)}^2}\right),$ where σ(x) is the local standard deviation of the image .sub.x.

Description

BRIEF DESCRIPTION OF THE DRAWINGS

[0014] The invention may take form in various components and arrangements of components, and in various steps and arrangements of steps. The drawings are only for purposes of illustrating the preferred embodiments and are not to be construed as limiting the invention.

[0015] FIG. 1 schematically illustrates an example computing system with an anti-correlated noise filter.

[0016] FIG. 2 schematically illustrates an example of the spectral data manipulation and filtering with the anti-correlated noise filter and a first algorithm.

[0017] FIG. 3 schematically illustrates an example of the spectral data manipulation and filtering with the anti-correlated noise filter and a second algorithm.

[0018] FIG. 4 flowcharts a method of filtering anti-correlated noise in image pairs.

[0019] FIG. 5 flowcharts a method of filtering anti-correlated noise in projection data set pairs.

DETAILED DESCRIPTION OF EMBODIMENTS

[0020] Initially referring to FIG. 1, an imaging system 100 includes an imaging scanner, such as a computed tomography (CT) scanner 102. The scanner is configured to generate projection data and/or images which decompose into anti-correlated portions. The scanner 102 includes one or more radiation sources 104, such as an x-ray tube, which emits radiation that traverses an examination region 106. In one instance, a mean or peak emission voltage of the radiation source 104 is switched between an emission voltage of two or more emission voltages (e.g., 80 and 140 kVp, 100 and 120 kVp, etc.). In another variation, the radiation source 104 includes a single broad spectrum x-ray tube.

[0021] A detector array 108 opposite the radiation source 104 detects the emitted radiation that traverses the examination region 106 and generates projection data 110 indicative of the object or subject in the examination region 106. Where the radiation source voltage is switched between at least two emission voltages and/or includes two or more x-ray tubes emit radiation at two different emission voltages, the detector array 108 generates projection data 110 for each of the radiation source voltages. For a single broad-spectrum x-ray tube, the detector array 108 includes an energy-resolving detector (e.g., multi-layered, photon counting, etc.) that produces the spectral projection data 110.

[0022] The projection data 110 can be represented as projection images, sinograms, and the like. The projection data 110 can include data organization, such as file and/or data set organization, database organizations, object and/or element definition, and the like. The projection data 110 can be stored in random access memory, such as computer memory, local memory, server memory, cloud storage, server storage, local storage, solid state storage, flash storage and the like. In one instance, the projection data 110 is stored in a Picture Archiving and Communication System (PACS) 112, Radiology Information System (RIS), Hospital Information System (HIS), Electronic Medical Record (EMR) or other system or device communicatively connected to the system 100.

[0023] A decomposition unit 114 decomposes the projection data into at least a first portion 116 and a second portion 118, such as a pair of data sets or data structures. Each portion includes noise, which is anti-correlated with noise in the other portion. For example, the decomposition unit 114 decomposes spectral projection data into photoelectric and Compton-scatter components, water and Iodine components, water and Calcium components, or acetal homopolymer resin, e.g., Delrin® and tin components, and/or other basis material sets.

[0024] An anti-correlated noise filter 120 jointly processes the first portion 116 and second portion 118. The anti-correlated noise filter 120 iteratively filters the first portion 116 and second portion 118. The anti-correlated noise filter 120 receives two projection data sets or image data sets as input and outputs two anti-correlated noise filtered projection data sets or image data sets. The iteration can be performed until a stopping criteria is reached, such as a predetermined number of iterations, a predetermined elapsed time, a difference between the input pairs and the output pairs satisfying a predetermined difference, combinations thereof, and the like. As described in greater detail below, the anti-correlated noise filter 120 jointly filters according to a function of a first algorithm or a second algorithm, which minimizes the anti-correlated noise. The functions are implemented by one or more methods, which filter out the anti-correlated noise in the output projection data sets or image data sets.

[0025] A reconstructor 122 reconstructs the filtered first portion 116 and second portion 118 into each into an image and/or a combined image. In one instance, the reconstructor 122 reconstructs each of the filtered first portion 116 and the filtered second portion 118, each into separate images. The separate images can be combined to form a combined image which is displayed on a display device 124. In another instance, the filtered first portion 116 and filtered second portion 118 are combined in projection space and then reconstructed as one image, which is displayed on the display device 124 and/or stored, such as in the PACS 112.

[0026] In one instance, the reconstructor 122 reconstructs the first portion 116 and second portion 118, e.g., from projection space, into images, e.g., into image space, and the anti-correlated noise filter 120 jointly filters the reconstructed images of the reconstructed first portion 116 and the reconstructed second portion 118. The filtered reconstructed first portion 116 and the filtered reconstructed second portion 118 can be combined into an image or otherwise manipulated. The combined image is displayed on the display device 124 and/or stored in the PACS 112.

[0027] The decomposition unit 114, the anti-correlated noise filter 120 and the reconstructor 122 are suitably embodied by a data processor 126, such as a electronic data processor, microprocessor, digital processor, optical processor, and the like, configured to execute computer readable instructions stored in a non-transitory computer readable storage medium or computer readable memory, e.g., software.

[0028] The processor 126 can also execute computer readable instructions carried by a carrier wave, a signal or other transitory medium to perform the disclosed techniques. The processor 126 can receive parameters through one or more input devices 128, such as a keyboard, mouse, touch screen, microphone, and the like. The processor 126, display device 124 and input device 128 can comprise a computing device 130, such as a desktop computer, laptop, smartphone, body worn device, distributedly connected computing devices, such as a server and a communicatively connected peer or client computer, and the like.

[0029] With reference to FIG. 2, an example of spectral data 200 manipulation and filtering with the anti-correlated noise filter and a first algorithm is schematically illustrated. The spectral data 200 can include either the projection data 110 or reconstructed image data. The decomposition unit 114 decomposes the spectral data 200 into the first portion 116, such as Compton scatter or scatter, and into the second portion 118, such as photoelectric or photo.

[0030] A Structure Propagation (SP) filter 202 can be used to initially filter the decomposed portions individually. The SP filter 202 can be a bi-lateral filter that uses information on edges of the combined image 200 or a weighted combination. The output is an initially filtered first portion 204 and an initially filtered second portion 206. In a variation, the SP filter 202 is omitted.

[0031] The anti-correlated noise filter 120 iteratively filters the initially filtered first portion 204 and an initially filtered second portion 206. In one instance, the initially filtered first portion 204 and an initially filtered second portion 206 are the first portion 116 and the second portion 118, e.g., without the SP filter 202. The anti-correlated noise filter 120 uses a minimization function in a first algorithm to iteratively filter the decomposed image volume pair s.sup.0 204 and p.sup.0 206 from the combined projection data or image volume u.sup.0. An example minimization function is shown in EQUATION 1:

[00006]$\begin{array}{cc}\left(\hat{p},\hat{s}\right)=\arg .Math..Math.{\min }_{\left(p,s\right)}.Math.R\left(p\right)+R\left(s\right)+\frac{1}{2}.Math.\int {{\lambda }^u\left(p+s-u^0\right)}^2+\frac{1}{2}.Math.\int {{\lambda }^p\left(p-p^0\right)}^2+\frac{1}{2}.Math.\int {{\lambda }^s\left(s-s^0\right)}^2,& \mathrm{EQUATION}.Math..Math.1\end{array}$

where R(p) and R(s) are roughness penalties or regularization terms for p and s, respectively, u.sup.0 is an image volume where the correlated noise maximally cancels out with the initially decomposed portions, p.sup.0 and s.sup.0, e.g., u.sup.0=p.sup.0+s.sup.0, p and s are the filtered image volumes, and λ.sup.u, λ.sup.p and λ.sup.s are weights. In one instance, the roughness penalty includes a Huber penalty, such as

[00007]$R\left(p\right)=\int h_{\delta }\left(\nabla p\right),\mathrm{where}.Math..Math.h_{\delta }\left(\nabla p\right)=\{\begin{array}{cc}\frac{1}{2}.Math.{\left(\nabla p\right)}^2& \mathrm{if}.Math..Math..Math.\nabla p.Math.\le \delta \\ \delta \left(.Math.\nabla p.Math.-\frac{1}{2}.Math.\delta \right)& \mathrm{otherwise}\end{array}$

and δ is a Huber parameter. R(s) is similarly constructed. In another instance, R(p)is a total variation penalty, such as R(p)=∫|Vp|. The constraint λ.sup.u (p+s−u.sup.0) ensures that the edges are maintained in the original minimal noise (sum) image volume u.sup.0, while the constraints λ.sup.p (p−p.sup.0) and λ.sup.s (s−s.sup.0) reduce the amount of crosstalk between the image volumes p and s. The weights provided by the lambda's can be adjusted to provide different trade-offs between smoothing and edge preservation and can vary across each image volume. If noise in the initial image volumes p.sup.0 and s.sup.0 is not perfectly anti-correlated,the image volumes can be balanced by scaling both portions, such as b.sub.1p.sup.0 and b.sub.2s.sup.0, where b.sub.1 and b.sub.2 are scaling factors are picked so that the anti-correlated noise completely cancels out in the sum u.sub.0=b.sub.1p.sup.0+b.sub.2s.sup.0.

[0032] The minimization function of EQUATION 1 can be implemented by iterative update functions shown in EQUATIONS 2 and 3:

[00008]$\begin{array}{cc}{s}_{i,j,k}^{n+1}=\frac{{\lambda }_{i,j,k}^u\left(u_{i,j,k}^0-p_{i,j,k}^n\right)+{\lambda }_{i,j,k}^s.Math.s_{i,j,k}^0+{\mathrm{\delta \Sigma }}_D.Math.{\sigma }_D^n.Math.s_D^n}{{\lambda }_{i,j,k}^u+{\lambda }_{i,j,k}^s+{\mathrm{\delta \Sigma }}_D.Math.{\sigma }_D^n}& \mathrm{EQUATION}.Math..Math.2\\ p_{i,j,k}^{n+1}=\frac{{\lambda }_{i,j,k}^u\left(u_{i,j,k}^0-s_{i,j,k}^n\right)+{\lambda }_{i,j,k}^p.Math.s_{i,j,k}^0+{\mathrm{\delta \Sigma }}_D.Math.{\varphi }_D^n.Math.s_D^n}{{\lambda }_{i,j,k}^u+{\lambda }_{i,j,k}^p+{\mathrm{\delta \Sigma }}_D.Math.{\varphi }_D^n},& \mathrm{EQUATION}.Math..Math.3\end{array}$

where λ.sup.u, λ.sup.p and λ.sup.2 are weights, is the initial combined image or projection data, p.sup.0 and s° are the initial input images or projection data, p.sup.n and s.sup.n are current values of the n.sup.th iteration, p.sup.n+1 and s.sup.n+1 are the next iteration filtered projection data or images, D includes a set of each orthogonal direction {E(ast), W(est), S(outh), N(orth), U(p), and (d)O(wn))}, and i,j,k represent a current voxel in the image or position in projection space volume, and δ is a Huber parameter. The weights, σ of s and φ of p, for each direction, such as E(east), are given by EQUATION 4 and EQUATION 5:

[00009]$\begin{array}{cc}{\sigma }_{E,1}^n=\sqrt{\begin{array}{c}\frac{{\left(s_E-s\right)}^2}{d_y^2}+\frac{{\left(s_N+s_{\mathrm{NE}}-s_S-s_{\mathrm{SE}}\right)}^2}{4.Math.d_x^2}+\\ \frac{{\left(s_U+s_{\mathrm{EU}}-s_O-s_{\mathrm{EO}}\right)}^2}{4.Math.d_z^2}\end{array}};& \mathrm{EQUATION}.Math..Math.4\\ {\varphi }_{E,1}^n=\sqrt{\begin{array}{c}\frac{{\left(p_E-p\right)}^2}{d_y^2}+\frac{{\left(p_N+p_{\mathrm{NE}}-p_S-p_{\mathrm{SE}}\right)}^2}{4.Math.d_x^2}+\\ \frac{{\left(p_U+p_{\mathrm{EU}}-p_O-p_{\mathrm{EO}}\right)}^2}{4.Math.d_z^2}\end{array}},& \mathrm{EQUATION}.Math..Math.5\end{array}$

where s and p are the corresponding current voxels, e.g., s.sub.i,j,k.sup.n and p.sub.i,j,k.sup.n, and subscripts E, N, NE, S, SE, U, EU, O, EO refer to the voxels or pixels east, north, northeast, south, southeast, up, east up, down, and east down, respectively of the current voxel or pixel, e.g., the voxels in orthogonally adjacent voxels, e.g., E, N, S, U, and O, which are not opposite or west, and two additional voxels east of the current voxel, which are not orthogonal and in the same direction, such as NE and SE, d.sub.xis the voxel distance in the North/South direction in the image volume, d.sub.y is the voxel distance in the East/West direction in the image volume, d.sub.z is the voxel distance in the Up/Down direction in the image volume. The weighing is derived analogously for each remaining direction. In each direction D ∈ (E, W, N, S, U, O), the weights σ.sub.D,1 and φ.sub.D,1 are adjusted to

[00010]${\sigma }_D=\frac{1}{d_D.Math.\max \left({\sigma }_{D,1},\delta \right)}.Math..Math.\mathrm{and}.Math..Math.{\varphi }_D=\frac{1}{d_D.Math.\max \left({\varphi }_{D,1},\delta \right)}$

and used in update equations 2 and 3, where d.sub.s, d.sub.N:=d.sub.x; d.sub.E, d.sub.W=d.sub.y; and d.sub.U, d.sub.O=d.sub.z.

[0033] The output of each iteration is an anti-correlated noise filtered first portion 208 and an anti-correlated second portion 210, such as the s.sup.n+1 and the p.sup.n+1 projection data or image volume. The output is used as the input in a next iteration. The anti-correlated noise filtered first portion 208 and the anti-correlated second portion 210 as projection data can be reconstructed into separate images and then combined as an anti-correlated filtered image 212. The anti-correlated noise filtered first portion 208 and the anti-correlated second portion 210 as reconstructed images can be combined as the anti-correlated filtered image 212.

[0034] With reference to FIG. 3, an example of the spectral data manipulation and filtering with the anti-correlated noise filter 120 and a second algorithm is schematically illustrated. The spectral data 200 can include either the projection data 110 or image data. The decomposition unit 114 decomposes the spectral data 200 into the first portion 116, such as the Compton scatter or scatter, and into the second portion 118, such as the photoelectric or photo.

[0035] The denoising Structure Propagation (SP) filter 202 can be used to initially filter the decomposed portions individually. In one instance, the SP filter is omitted. The output is an initially filtered first portion 204 and an initially filtered second portion 206. The anti-correlated noise filter 120 iteratively filters the initially filtered first portion 204 and an initially filtered second portion 206. In one instance, the initially filtered first portion 204 and an initially filtered second portion 206 are the first portion 116 and the second portion 118, e.g., without the SP filter 202. The anti-correlated noise filter 120 uses a minimization function of a second algorithm to iteratively filter the decomposed image pair s.sup.0 204 and p.sup.0 206 from the combined projection data or image u.sup.0. An example of the constrained minimization function of the second algorithm is shown in EQUATION 6:

[00011]$\begin{array}{cc}\left(\hat{s},\hat{p}\right)=\arg .Math..Math.\underset{\left(s,p\right)}{\min }.Math.\alpha .Math..Math.R\left(s\right)+\left(1-\alpha \right).Math.R\left(p\right)& \mathrm{EQUATION}.Math..Math.6\end{array}$ [0036] subject to the constraints that (s.t.) [0037] 1. s and p are obtained by removing negatively correlated estimated noise from s.sup.0 and p.sup.0, respectively; [0038] 2. {circumflex over (m)} monochromatic image is unchanged; and [0039] 3. image frequencies outside band frequencies are unchanged,
where R(p) and R(s) are roughness penalties or regularization terms for p and s, respectively, such as R(x)=∫|Vx|, {circumflex over (m)} is an energy level parameter in keV unit, and α is an algorithm control parameter with a default value, such as equal to 0.5.

[0040] The minimization constrained by the condition that the virtual monochromatic image or projection data with energy {circumflex over (m)} does not change. The energy {circumflex over (m)} 300 is selected such that its anti-correlated noise is minimized. In one instance, {circumflex over (m)} is defined by:

[00012]$\begin{array}{cc}\hat{m}=\arg .Math..Math.\underset{m}{\min }.Math.\frac{R\left(c_s\left(m\right).Math.s+c_P\left(m\right).Math.p\right)}{c_s\left(m\right).Math.R\left(s\right)+c_p\left(m\right).Math.R\left(p\right)},& \mathrm{EQUATION}.Math..Math.7\end{array}$

where R is a regularization function, c.sub.s(m) and c.sub.p(m) are the coefficients of the s and p, respectively, to obtain the monochromatic image for energy m in keV.

[0041] In another instance, the spectral virtual monochromatic image, {circumflex over (m)}, is detected by defining a selection region of the combined spectral data 200 using a predetermined threshold value, such as −200 HU. A local standard deviation is calculated for a neighborhood of size ne for the combined spectral data 200. A set, q of locations is created of the r smallest local standard deviations located in the selection region and an example of {circumflex over (m)} is defined in EQUATION 8:

[00013]$\begin{array}{cc}\hat{m}=\arg .Math..Math.\underset{m}{\min }.Math.\Sigma .Math..Math.\mathrm{localstddev}\left(c_s\left(m\right).Math.s+c_p\left(m\right).Math.p,q,\mathrm{ne}\right),& \mathrm{EQUATION}.Math..Math.8\end{array}$

where the local standard deviation is calculated only over the set q and ne specifies the neighbourhood of the local standard deviation.

[0042] The minimization function of EQUATION 6 is constrained to process image frequencies only within a band of frequencies, this is obtained as follows: the input p and s images are down scaled by a factor of d, which includes the range (0,1] and 1 is no scaling, 0.5 is scaling by a factor of 2, 0.25 is scaling down by a factor of 4, etc. Scaling down the images combines pixels through averaging and/or interpolation. The following optimization, EQUATION 9, is performed, using for example the lagged diffusivity fixed point iteration algorithm:

[00014]$\begin{array}{cc}\hat{A}=\arg .Math..Math.\underset{A}{\min }.Math.\alpha .Math..Math.R\left(s_d+c_p\left(\hat{m}\right).Math.A\right)+\left(1-\alpha \right).Math.R\left(p_d-c_s\left(\hat{m}\right).Math.A\right)& \mathrm{EQUATION}.Math..Math.9\end{array}$

where α is an algorithm control parameter, R is the roughness penalties or regularization terms, c.sub.s({circumflex over (m)}) and c.sub.p({circumflex over (m)}) are the coefficients of the downscaled images s.sub.d and p.sub.d, respectively, that enable the algorithm to obtain the monochromic image for energy {circumflex over (m)} keV, d is the downscaling factor, and  is the estimate of the anti-correlated noise image. The denoised down scaled s and p image are defined in EQUATION 10 and 11:

s.sub.d.sup.denoised=s.sub.d+c.sub.p({circumflex over (m)})Â and p.sub.d.sup.denoised=p.sub.d−c.sub.s({circumflex over (m)})Â.   EQUATION 10 and 11

The denoised down scaled monochromatic image at energy {circumflex over (m)} in keV doesn't change, as shown in an example EQUATION 12:

[00015]$\begin{array}{cc}{I}_{d,\hat{m}.Math.\mathrm{ke}.Math..Math.V}^{\mathrm{denoised}}=c_p\left(\hat{m}\right).Math.p_d^{\mathrm{denoised}}+c_p\left(\hat{m}\right).Math.s_d^{\mathrm{denoised}}=c_p\left(\hat{m}\right).Math.\left(p_d-c_s\left(\hat{m}\right).Math.\hat{A}\right)+c_s\left(\hat{m}\right).Math.\left(s_d+c_p\left(\hat{m}\right).Math.\hat{A}\right)=c_p\left(\hat{m}\right).Math.p_d+c_s\left(\hat{m}\right).Math.s_d=I_{d,\hat{m}.Math.\mathrm{ke}.Math..Math.V}& \mathrm{EQUATION}.Math..Math.12\end{array}$

The denoised down scaled monochromatic images are scaled up to generate the {circumflex over (p)} and ŝ output images, in which only image frequencies within a band of frequencies were processed. Examples of scaling up the denoised down scaled monochromatic images are shown in EQUATIONS 13-16:

s.sub.t=s+ScaleUp(s.sub.d+c.sub.p ({circumflex over (m)})Â−ScaleDown (s, d), d)=s+ScaleUp(c.sub.p ({circumflex over (m)})Â, d) and ŝ=s.sub.t+ScaleUp (ScaleDown(s−s.sub.t, u), u)   EQUATION 13 and 14

p.sub.t=p+ScaleUp(p.sub.d−c.sub.s({circumflex over (m)})Â−ScaleDown (p,d),d)=p−ScaleUp(c.sub.s({circumflex over (m)})Â, d) and {circumflex over (p)}=p.sub.t+ScaleUp (ScaleDown(p−p.sub.t, u), u),   EQUATION 15 and 16

where u is a scale factor parameter, such as u=1/16.

[0043] The operator s.sub.u=ScaleUp(s,u) returns the image s.sub.u that is

[00016]$\frac{1}{u}$

times the size of the image s, and the operator s.sub.d =ScaleDown(s,d) returns the image s.sub.d that is d times the size of the image s.sub.d=ScaleDown(s,d) returns the image s.sub.d that is d times the size of the image s. Note that s.sub.t and p.sub.t represent intermediate results in which the image high frequency are preserved, i.e., the high frequencies of the images s.sub.t and p.sub.t are very similar to the high frequencies of the input images s and p, respectively.

[0044] The output of each iteration is an anti-correlated noise filtered first portion 208 and an anti-correlated second portion 210, such as the Ŝ and the {circumflex over (P)} projection data or image. The output is used as the input in a next iteration. The anti-correlated noise filtered first portion 208 and the anti-correlated second portion 210 as projection data can be reconstructed into separate images and then combined as an anti-correlated filtered image 212. The anti-correlated noise filtered first portion 208 and the anti-correlated second portion 210 as reconstructed images can be combined as the anti-correlated filtered image 212. In another instance, the filter includes minimization of the following two functionals:

[00017]$\begin{array}{cc}\hat{A}=\underset{A}{\arg .Math..Math.\min }.Math..Math.R\left(s^0+A\right)+R\left(p^0-A\right)+\int {\lambda }_1.Math.h_{\delta }\left(A\right).Math..Math.\hat{A}=\underset{A}{\arg .Math..Math.\min .Math.}.Math.R\left(s^0+A\right)+R\left(p^0-A\right)+\int {\lambda }_2.Math.A^2,& \mathrm{EQUATION}.Math..Math.17.Math..Math.\mathrm{and}.Math..Math.18\end{array}$

where R(.Math.) is a roughness penalty or regularization term, λ.sub.1 and λ.sub.2 are weights, Â is the estimated anti-correlated noise image, h.sub.67 (A)=δ.sup.2(√{square root over (1+(A/δ).sup.2)}−1) is the pseudo-Huber penalty function and δ is the pseudo-Huber parameter. EQUATION 17 is applied to the p.sup.0 and s.sup.0 images using, for example, a lagged diffusivity fixed point iteration algorithm according to EQUATION 19.

[00018]$\begin{array}{cc}{\hat{A}}_{L_1}=\underset{A}{\arg .Math..Math.\min }.Math..Math.R\left(s^0+A\right)+R\left(p^0-A\right)+\int {\lambda }_1.Math.h_{\delta }\left(A\right).& \mathrm{EQUATION}.Math..Math.19\end{array}$

EQUATION 18 is applied to on the scaled down p.sup.0 and s.sup.0 images, with an optimization performed using, for example, the lagged diffusivity fixed point iteration algorithm according to EQUATION 20.

[00019]$\begin{array}{cc}{\hat{A}}_{L_2}^d=\underset{A}{\arg .Math..Math.\min }.Math..Math.R\left(s_d^0+A\right)+R\left(p_d^0-A\right)+\int {\lambda }_2.Math.A^2,& \mathrm{EQUATION}.Math..Math.20\end{array}$

where s.sub.d.sup.0=ScaleDown (s.sup.0, d) and p.sub.d.sup.0=ScaleDown (p.sup.0, d). The filtered images are obtained by combining the estimates from EQUATIONS 19 and 20 according to EQUATION 21.

Â=Â.sub.L.sub.1+ScaleUp (Â.sub.L.sub.2.sup.d−ScaleDown (Â.sub.L.sub.1, d),d), {circumflex over (p)}=p.sup.0−Â and ŝ=s.sup.0+Â,   EQUATION 21

where  is the final estimate for the anti-correlated noise image, {circumflex over (p)} and ŝ are the denoised photo and scatter images, respectively.

[0045] A third algorithm uses an estimated noise map and filters the anti-correlated noise with an minimization function as defined in EQUATION 22. The estimate noice map, n, is estimated using noise modeling techniques, such as with a Monte Carlo estimate, by analytical methods such as by Wunderlich and Noo in “Image Covariance and Lesion Detectability in Direct Fan-Beam X-Ray Computed Tomorgraph”, or direct extraction techniques, such as described in U.S. Pat. No. 8,938,110.

Â=Â.sub.L.sub.1+ScaleUp(Â.sub.L.sub.1.sup.d, d) {circumflex over (p)}=p.sup.0−Â and ŝ=s.sup.0+Â,   EQUATION 22

where  is the estimated anti-correlated noise image, {circumflex over (p)} and ŝ are the denoised photo and scatter images, respectively, and p.sup.0 and s.sup.0 are the initial or input photo and scatter images, respectively, d is a scale parameter, and Â.sub.L.sub.1 and Â.sub.L.sub.2.sup.d are defined by EQUATION 23 and 24 respectively.

[00020]$\begin{array}{cc}{\hat{A}}_{L_1}=\underset{A}{\arg .Math..Math.\min }.Math..Math.R\left(s^0+A\right)+R\left(p^0-A\right)+\int \frac{{\lambda }_1}{n}.Math.h_{\delta }\left(A\right).& \mathrm{EQUATION}.Math..Math.23\end{array}$

where R is a roughness penalty or regularization term, such as ∫ h.sub.δ(|Vp|), A is the image, n is the estimated noise map, λ.sub.1 is a weight, and h.sub.67 (A)=δ.sup.2(√{square root over (1+(A/δ).sup.2)}−1) is the pseudo-Huber penalty function and δ is the pseudo-Huber parameter. EQUATION 23 is applied to on the scaled down p.sup.0 and s.sup.0 images, with an optimization performed using, for example, the lagged diffusivity fixed point iteration algorithm according to EQUATION 24.

[00021]$\begin{array}{cc}{\hat{A}}_{L_2}^d=\underset{A}{\arg .Math..Math.\min }.Math..Math.R\left(s_d^0+A\right)+R\left(p_d^0-A\right)+\int \left(\frac{{\lambda }_2}{n^2}+{\lambda }_3.Math.F\right).Math.A^2,& \mathrm{EQUATION}.Math..Math.24\end{array}$

where s.sub.d.sup.0=ScaleDown (s.sup.0+Â.sub.L.sub.1, d) and p.sub.d.sup.0=ScaleDown (p.sup.0−Â.sub.L.sub.1, d), λ.sub.2 and λ.sub.3 are weights and F is defined by EQUATION 25.

[00022]$\begin{array}{cc}F=\mathrm{ScaleDown}\left(\frac{\sigma \left(s_d^0+{\hat{A}}_{L_1}\right).Math.\sigma \left(p_d^0-{\hat{A}}_{L_1}\right)}{1+{\sigma \left({\hat{A}}_{L_1}\right)}^2}\right),& \mathrm{EQUATION}.Math..Math.25\end{array}$

where σ(x) is the local standard deviation of the image .sup.x.

[0046] With reference to FIG. 4, a method of filtering anti-correlated noise in image pairs is flowcharted. It is to be appreciated that the ordering of the acts is not limiting. As such, other orderings are contemplated herein. In addition, one or more acts may be omitted and/or one or more additional acts may be included.

[0047] At 400, spectral image data is received, such as the spectral image data 200 described in reference to FIGS. 2 and 3. The spectral image data can be received from memory, such as the PACS 112, or from the reconstructor 122.

[0048] The spectral image data is decomposed at 402 into the first decomposed portion, such as Compton scatter, and the second decomposed portion, such as the photoelectric absorption. The decomposed first and second portions include anti-correlated noise.

[0049] At 404, each of the decomposed first and second portions can be individually filtered. For example, the first and second portions are each filtered separately with the SP filter 202 using information from the combined image 200 or spectral image data.

[0050] At 406, the decomposed pair of images from 402 or the individual filtered pair from 404 are iteratively filtered as a pair using the anti-correlated noise filter 120 as described in reference to FIG. 1 and as examples in reference to FIGS. 2 and 3. The iterative algorithms can be repeated until a stopping criteria is reached, such as a number of iterations, an elapse time, or a threshold change between characteristics of the input images and the output images, a combination, and the like.

[0051] At 408, the output images are combined in memory. The output images can be further manipulated, such as other filtering, computations, segmentation, and the like before combining.

[0052] The combined image is displayed on a display device 124 at 410. The display can include displaying the output images, i.e., anti-correlated noise filtered images. The combined image can also be output on film by a filmer or the like.

[0053] The above may be implemented by way of computer readable instructions, encoded or embedded on non-transitory computer readable storage medium, which, when executed by a computer processor(s), cause the processor(s) to carry out the described acts. Additionally or alternatively, at least one of the computer readable instructions is carried by a signal, carrier wave or other transitory medium.

[0054] With reference to FIG. 5, a method of filtering anti-correlated noise in projection data set pairs is flowcharted. At 500, spectral image data is received, such as the spectral image data 200 described in reference to FIGS. 2 and 3. The spectral image data can be received from memory, such as the PACS 112, or from the imaging scanner, such as the CT scanner 102.

[0055] The spectral image data is decomposed at 502 into the first decomposed portion, such as Compton scatter, and the second decomposed portion, such as the photoelectric absorption. The decomposed first and second portions include anti-correlated noise.

[0056] At 504, each of the decomposed first and second portions can be individually filtered. For example, the first and second portions are each filtered separately with the SP filter 202 using information from the combined image 200 or spectral image data.

[0057] At 506, the decomposed pair of projection data sets from 502 or the individual filtered pair from 504 are iteratively filtered as a pair using the anti-correlated noise filter 120 as described in reference to FIG. 1 and as examples in reference to FIGS. 2 and 3. The iterative algorithms can be repeated until a stopping criteria is reached, such as a number of iterations, an elapse time, or a threshold change between characteristics of the input images and the output images, a combination, and the like.

[0058] At 508, the output data sets are combined in memory. The output data sets can be further manipulated, such as other filtering, computations, segmentation, and the like before combining. At 510, the combined output data set is reconstructed by the reconstructor 122 into an image. Alternatively, 510 and 508 can be reversed, and the output data sets are reconstructed each into an image, and the two reconstructed images are combined.

[0059] The combined image is displayed on a display device 124 at 512. The display can include displaying the output images, i.e., anti-correlated noise filtered images. The combined image can also be output on film by a filmer or the like.

[0060] The above may be implemented by way of computer readable instructions, encoded or embedded on non-transitory computer readable storage medium, which, when executed by a computer processor(s), cause the processor(s) to carry out the described acts.

[0061] Additionally or alternatively, at least one of the computer readable instructions is carried by a signal, carrier wave or other transitory medium.

[0062] The invention has been described with reference to the preferred embodiments. Modifications and alterations may occur to others upon reading and understanding the preceding detailed description. It is intended that the invention be constructed as including all such modifications and alterations insofar as they come within the scope of the appended claims or the equivalents thereof.