Login| Sign Up| Help| Contact|

Patent Searching and Data


Title:
METHODS OF DETERMINING LOCAL SPECTRUM AT A PIXEL USING A ROTATIONALLY INVARIANT S-TRANSFORM (RIST)
Document Type and Number:
WIPO Patent Application WO/2013/078451
Kind Code:
A1
Abstract:
An image processing device and methods for performing Rotationally Invariant Stransform (RIST) for an image are provided herein. An example method of determining the RIST magnitude at a pixel is provided herein. Further, an example method of determining RIST magnitudes and statistics in a region of interest is provided herein.

Inventors:
CHENG CHUN HING (CA)
MITCHELL JOSEPH ROSS (US)
Application Number:
PCT/US2012/066450
Publication Date:
May 30, 2013
Filing Date:
November 23, 2012
Export Citation:
Click for automatic bibliography generation   Help
Assignee:
MAYO FOUNDATION (US)
International Classes:
G06K9/56
Foreign References:
US20090161945A12009-06-25
US20060083407A12006-04-20
US6803919B12004-10-12
US20100008589A12010-01-14
Attorney, Agent or Firm:
AARONSON, Lawrence et al. (Meunier Carlin & Curfman, LLC,Suite 500,817 W. Peachtree Street N, Atlanta Georgia, US)
Download PDF:
Claims:
WHAT IS CLAIMED:

1. A method of determining rotational invariant local spectrum at a pixel in an image processing device, comprising:

receiving an input image;

receiving an input coordinate of the pixel; and

determining the values of a rotational invariant form of two-dimensional S-Transform (RIST) at the input coordinate.

2. The method of claim 1, determining RIST further comprising:

determining the S-Transform (ST) magnitudes (Al) using positive discretization at the input coordinate of the pixel;

flipping the input image along x direction;

determining the ST magnitudes (A2) using positive discretization at the coordinate of the corresponding pixel in the flipped image; and

determining the average of the above two sets Al and A2 of magnitudes.

3. The method of claim 1, determining RIST further comprising:

determining RIST at a pixel using a modified form of the method in FTFT-2D.

4. The method of claim 1, determining RIST further comprising:

determining RIST values and statistics in a region of interest (ROI) using a modified form of the method in FTFT-2D.

5. The method of claim 3, further comprising:

setting parameters;

preparing basis values;

receiving an input image;

determining a two-dimensional Fourier Transform (FT) of the image as a matrix H; receiving an input coordinate of the pixel;

determining the ST magnitudes (Bl) using positive discretization at the input coordinate of the pixel using the matrix H and the parameters,

flipping the input image along x direction; determining the ST magnitudes (B2) using positive discretization at the coordinate of the corresponding pixel in the flipped image; and

determining the average of the above two sets Bl and B2 of magnitudes.

6. The method of claim 4, further comprising:

setting parameters;

preparing basis values;

receiving an input image;

determining a two-dimensional Fourier Transform (FT) of the image as a matrix H; receiving an indication of the region on interest (ROI);

determining the ST magnitudes (CI) using positive discretization in the ROI using the matrix H and the parameters;

flipping the input image along x direction;

determining the ST magnitudes (C2) using positive discretization in the corresponding ROI in the flipped image; and

determining the average of the above two sets CI and C2 of magnitudes.

7. The method of claim 5, further comprising:

if the width Nx and height Ny of the input image are not both equal to N, wherein N is a power of 2, then:

determine a smallest integer M such that Nx≤ 2M and Ny≤ 2M;

set N = 2M; and

adjust a size of the input image by expanding the input image into an NxN image by optimized Hanning window.

8. The method of claim 7, preparing basis values for each of the low band, the medium band and the high band further comprising:

determining support intervals for each pure complex sinusoid;

determining a range of PCS, the range being for ST values for values of frequency index k = 0 through N/2 - 1;

identifying a low set of PCS with a relatively small frequency index q, wherein the ST are copied into the basis; identifying a medium set of PCS with a frequency index between the relatively small frequency index q of the low set of PCS and a relatively large frequency index q, wherein the Offset TT-Transform (OTT) are used in the basis;

determining crop limits for each pure complex sinusoid in the medium set;

identifying basis nodes for each pure complex sinusoid in the medium set;

identifying a high set of PCS with the relatively large frequency index q, wherein the Offset TT-Transform (OTT) are used in the basis;

determining crop limits for each pure complex sinusoid in the high set;

identifying basis nodes for each pure complex sinusoid in the high set;

subsampling along a time axis; and

determining basis values for each pure complex sinusoid in the high set, the medium set and the low set.

9. The method of claim 5, determining the ST magnitudes further comprising:

multiplying a matrix of basis values for N to the matrix H on the left to form an intermediate matrix product; and

multiplying a transpose of matrix of basis values for N to the intermediate matrix on the right to form a matrix product of compressed ST magnitudes for the pixel.

10. The method of claim 5, further comprising:

interpolating the matrix of compressed ST values along an x direction; and

interpolating a result along a y direction to obtain a matrix of semi-compressed ST values for the pixel.

11. The method of claim 10, further comprising:

decompressing the matrix of semi-compressed ST values for the pixel along the x direction; and

decompressing a result along the y direction to obtain a matrix of the ST values at the input coordinate.

12. The method of claim 6, preparing basis further comprising: determining the basis values for the image width Nx using the primary parameters along an x direction; and

determining the basis values for the image height Ny using the primary parameters along a y direction.

13. The method of claim 6, determining the ST values further comprising determining a bounding rectangle of the ROI.

14. The method of claim 13, wherein if an x-length of the ROI is greater than a y-length, then the method further comprises:

forming an intermediate matrix product for all ix in an x-projection of the ROI;

traversing a pixel tree; and

for each node P(ix, iy), if it is in the ROI and not computed before, then multiplying a matrix of basis values for iy to the intermediate matrix product on the right to form a matrix of compressed ST values for the pixel.

15. The method of claim 13, wherein if an x-length of the ROI is not greater than a y- length, then the method further comprising:

forming an intermediate matrix product for all iy in a y-projection of the ROI;

traversing a pixel tree; and

for each node P(ix, iy), if it is in the ROI and not computed before, then multiplying a matrix basis values for iy to the intermediate matrix product on the left to form a matrix of compressed ST values for the pixel.

16. The method of claim 6, determining ST in the ROI further comprising determining a local spectrum at each pixel (ix, iy) in the ROI.

17. The method of claim 6, determining ST in an ROI further comprising augmenting weights and updating statistics.

18. The method of claim 6, further comprising: determining a low band, a medium band and a high band of frequency components; and

selecting a skipping strategy to skip computing predetermined ones of the ST values.

19. The method of claim 18, further comprising:

building a forest of quad-trees with two levels;

selecting pixels at every other x position and every other y position;

for a first two leaves of each tree, corresponding to a pair of diagonally opposite pixels, computing ST values for the low band, the medium band and the high band;

determining an upper-difference between ST values of these two pixels at each (kx, ky) in an upper quadrant of a 2D frequency index space; and

if the upper-difference is less than a predetermined threshold, skipping computing ST values in the low band, the medium band and the high band for other two leaves in that tree.

20. The method of claim 18, further comprising:

determining low band ST values for each 2 x 2 square of the ROI; and

skipping determining the ST values for the medium band and the high band if a predetermined selection of high band ST magnitude is less than a threshold.

21. The method of claim 18, further comprising:

determining low band ST values for each 4 x 4 square of the ROI;

determining medium band ST values for each 2 x 2 square of the ROI;

building a forest of quad-trees having three levels, wherein at a top level, every fourth x position and every fourth y position is selected;

traversing children from a selected x position and y position; and

determining a ST value of a pixel in accordance with:

if that node is the top level of the tree, then determine its ST values for the low band, the medium band and the high band;

if that node is in a middle level, then determine the ST values for the medium band and the high band; and

if that node is in a lower level, then determine ST values for the high band.

22. The method of claim 18, further comprising performing an automatic selection of a skipping strategy.

23. The method of claim 6, further comprising applying a weight to the ST values.

24. The method of claim 1, further comprising determining the RIST value as a complex number at a point (nx, ny) wherein the input image is an N x N square image.

25. The method of claim 24, further comprising:

determining the complex number in accordance with the relationship: if *χ > 0, * > 0

-kx, ky] if kx < 0, ky > 0

kx, -k ] if k > 0, k < 0

- l - ny,-kx,-ky] if kx < 0, ky < 0

26. The method of claim 25, wherein kx and ky ca n take positive and negative values within N/2 - 1, -1, 0, 1, N/2 - 1.

27. The method of claim 25, further comprising expressing the relationship in a simplifie

28. The method of claim 24, further comprising:

displaying a semicircle; and

averaging over the semicircle of radius r to determine a texture curve.

29. A tangible computer readable medium having computer executable instructions thereon that when executed by a processor of a computing device perform a method of determining rotational invariant local spectrum at a pixel in an image processing device, comprising:

receiving an input image; receiving an input coordinate of the pixel; and

determining the values of a rotational invariant form of two-dimensional S-Transform (RIST) at the input coordinate.

30. The tangible computer readable medium of claim 24, determining RIST further comprising instructions for:

determining the S-Transform (ST) magnitudes (Al) using positive discretization at the input coordinate of the pixel;

flipping the input image along x direction;

determining the ST magnitudes (A2) using positive discretization at the coordinate of the corresponding pixel in the flipped image; and

determining the average of the above two sets Al and A2 of magnitudes.

31. The tangible computer readable medium of claim 24, determining RIST further comprising instructions for:

determining RIST at a pixel using a modified form of the method in FTFT-2D.

32. The tangible computer readable medium of claim 24, determining RIST further comprising instructions for:

determining RIST values and statistics in a region of interest (ROI) using a modified form of the method in FTFT-2D.

33. The tangible computer readable medium of claim 26, further comprising instructions for:

setting parameters;

preparing basis values;

receiving an input image;

determining a two-dimensional Fourier Transform (FT) of the image as a matrix H;

receiving an input coordinate of the pixel;

determining the ST magnitudes (Bl) using positive discretization at the input coordinate of the pixel using the matrix H and the parameters,

flipping the input image along x direction; determining the ST magnitudes (B2) using positive discretization at the coordinate of the corresponding pixel in the flipped image; and

determining the average of the above two sets Bl and B2 of magnitudes.

34. The tangible computer readable medium of claim 27, further comprising instructions for:

setting parameters;

preparing basis values;

receiving an input image;

determining a two-dimensional Fourier Transform (FT) of the image as a matrix H;

receiving an indication of the region on interest (ROI);

determining the ST magnitudes (CI) using positive discretization in the ROI using the matrix H and the parameters;

flipping the input image along x direction;

determining the ST magnitudes (C2) using positive discretization in the corresponding ROI in the flipped image; and

determining the average of the above two sets CI and C2 of magnitudes.

35. The tangible computer readable medium of claim 28, further comprising

instructions for:

if the width Nx and height Ny of the input image are not both equal to N, wherein N is a power of 2, then:

determine a smallest integer M such that Nx≤ 2M and Ny≤ 2M;

set N = 2M; and

adjust a size of the input image by expanding the input image into an NxN image by optimized Hanning window.

36. The tangible computer readable medium of claim 30, preparing basis values for each of the low band, the medium band and the high band further comprising instructions for: determining support intervals for each pure complex sinusoid;

determining a range of PCS, the range being for ST values for values of frequency index k = 0 through N/2 - 1; identifying a low set of PCS with a relatively small frequency index q, wherein the ST are copied into the basis;

identifying a medium set of PCS with a frequency index between the relatively small frequency index q of the low set of PCS and a relatively large frequency index q, wherein the Offset TT-Transform (OTT) are used in the basis;

determining crop limits for each pure complex sinusoid in the medium set;

identifying basis nodes for each pure complex sinusoid in the medium set;

identifying a high set of PCS with the relatively large frequency index q, wherein the Offset TT-Transform (OTT) are used in the basis;

determining crop limits for each pure complex sinusoid in the high set;

identifying basis nodes for each pure complex sinusoid in the high set;

subsampling along a time axis; and

determining basis values for each pure complex sinusoid in the high set, the medium set and the low set.

37. The tangible computer readable medium of claim 28, determining the ST

magnitudes further comprising instructions for:

multiplying a matrix of basis values for N to the matrix H on the left to form an intermediate matrix product; and

multiplying a transpose of matrix of basis values for N to the intermediate matrix on the right to form a matrix product of compressed ST magnitudes for the pixel.

38. The tangible computer readable medium of claim 28, further comprising

instructions for:

interpolating the matrix of compressed ST values along an x direction; and

interpolating a result along a y direction to obtain a matrix of semi-compressed ST values for the pixel.

39. The tangible computer readable medium of claim 33, further comprising

instructions for:

decompressing the matrix of semi-compressed ST values for the pixel along the x direction; and decompressing a result along the y direction to obtain a matrix of the ST values at the input coordinate.

40. The tangible computer readable medium of claim 29, preparing basis further comprising instructions for:

determining the basis values for the image width Nx using the primary parameters along an x direction; and

determining the basis values for the image height Ny using the primary parameters along a y direction.

41. The tangible computer readable medium of claim 29, determining the ST values further comprising instructions for determining a bounding rectangle of the ROI.

42. The tangible computer readable medium of claim 36, wherein if an x-length of the ROI is greater than a y-length, then the method further comprises instructions for:

forming an intermediate matrix product for all ix in an x-projection of the ROI;

traversing a pixel tree; and

for each node P(ix, iy), if it is in the ROI and not computed before, then multiplying a matrix of basis values for iy to the intermediate matrix product on the right to form a matrix of compressed ST values for the pixel.

38. The tangible computer readable medium of claim 36, wherein if an x-length of the ROI is not greater than a y-length, then the method further comprising instructions for:

forming an intermediate matrix product for all iy in a y-projection of the ROI;

traversing a pixel tree; and

for each node P(ix, iy), if it is in the ROI and not computed before, then multiplying a matrix basis values for iy to the intermediate matrix product on the left to form a matrix of compressed ST values for the pixel.

39. The tangible computer readable medium of claim 29, determining ST in the ROI further comprising instructions for determining a local spectrum at each pixel (ix, iy) in the ROI.

40. The tangible computer readable medium of claim 29, determining ST in an ROI further comprising instructions for augmenting weights and updating statistics.

41. The tangible computer readable medium of claim 29, further comprising

instructions for:

determining a low band, a medium band and a high band of frequency components; and

selecting a skipping strategy to skip computing predetermined ones of the ST values.

42. The tangible computer readable medium of claim 41, further comprising

instructions for:

building a forest of quad-trees with two levels;

selecting pixels at every other x position and every other y position;

for a first two leaves of each tree, corresponding to a pair of diagonally opposite pixels, computing ST values for the low band, the medium band and the high band;

determining an upper-difference between ST values of these two pixels at each (kx, ky) in an upper quadrant of a 2D frequency index space; and

if the upper-difference is less than a predetermined threshold, skipping computing ST values in the low band, the medium band and the high band for other two leaves in that tree.

43. The tangible computer readable medium of claim 41, further comprising instructions for:

determining low band ST values for each 2 x 2 square of the ROI; and

skipping determining the ST values for the medium band and the high band if a predetermined selection of high band ST magnitude is less than a threshold.

44. The tangible computer readable medium of claim 41, further comprising

instructions for:

determining low band ST values for each 4 x 4 square of the ROI;

determining medium band ST values for each 2 x 2 square of the ROI;

building a forest of quad-trees having three levels, wherein at a top level, every fourth x position and every fourth y position is selected; traversing children from a selected x position and y position; and

determining a ST value of a pixel in accordance with:

if that node is the top level of the tree, then determine its ST values for the low band, the medium band and the high band;

if that node is in a middle level, then determine the ST values for the medium band and the high band; and

if that node is in a lower level, then determine ST values for the high band.

45. The tangible computer readable medium of claim 41, further comprising instructions for performing an automatic selection of a skipping strategy.

46. The tangible computer readable medium of claim 29, further comprising instructions for applying a weight to the ST values.

47. The tangible computer readable of claim 29, further comprising instructions for determining the RIST value as a complex number at a point (nx, ny) wherein the input image is an N x N square image.

48. The tangible computer readable medium of claim 47, further comprising instructions for:

determining the complex number in accordance with the relationship:

49. The tangible computer readable medium of claim 48, wherein kx and ky can take positive and negative values within N/2 - 1, -1, 0, 1, N/2 - 1.

50. The tangible computer readable medium of claim 48, further comprising expressing the relationship in a simplified as:

51. The tangible computer readable medium of claim 47, further comprising instructions for:

displaying a semicircle; and

averaging over the semicircle of radius r to determine a texture curve.

Description:
METHODS OF DETERMINING LOCAL SPECTRUM AT A PIXEL USING A ROTATIONALLY

INVARIANT S-TRANSFORM (RIST)

CROSS-REFERENCE TO RELATED APPLICATION

[0001] This application claims the benefit of U.S. Provisional Patent Application No. 61/562,504, filed on November 22, 2011, entitled "RIST Patent Detailed Description," the disclosure of which is expressly incorporated herein by reference in its entirety.

BACKGROUND

[0002] Continuous S-transform (ST) ca n be regarded as a hybrid of Gabor and continuous wavelet transforms, providing a "time frequency representation" (TFR) of a signal by localizing with a Gaussian window that depends on the frequency. Its discrete 1- dimensional form (ID ST) is finding many applications in processing signals and time series, while its discrete 2-dimensional form (2D ST) is used for processing 2-dimensional data and images, where it should be more correctly called a "space frequency representation" (SFR), as it represents the localized frequency spectrum at each point in the 2-dimensional data set or at each pixel in the image.

[0003] Fast Time Frequency Transform tools have been developed, such as a FTFT-1D and FTFT-2D (Fast Time Frequency Transform), that generate discrete ID ST values and 2D ST magnitudes fast and accurately. The FTFT-2D can produce local ST magnitudes at each pixel in a medical image, as well as ST statistics over a region of interest (ROI) in the image. However, the discretization of 2D ST renderings are not rotationally invariant. By rotational invariance of an SFR, it is meant that when the image is rotated by any angle, the radial component of the SFR is unchanged. This is desirable as the pathology inferred from this radial component should not be affected when the patient is positioned at a different orientation on the imaging couch. SUMMARY

[0004] A method of determining rotational invariant local spectrum at a pixel in an image processing device. The method may include receiving an input image; receiving an input coordinate of the pixel; and determining the values of a rotational invariant form of two- dimensional S-Transform (RIST) at the input coordinate.

[0005] In some implementations, the method further includes determining the S- Transform (ST) magnitudes (Al) using positive discretization at the input coordinate of the pixel; flipping the input image along x direction; determining the ST magnitudes (A2) using positive discretization at the coordinate of the corresponding pixel in the flipped image; and determining the average of the above two sets Al and A2 of magnitudes.

[0006] The RIST algorithm may be implemented using a modified form of a FTFT-2D method.

[0007] In some implementations, the method may be implemented by a computing device executing the method as computer-executable instructions read from a tangible computer-readable medium.

[0008] It should be understood that the above-described subject matter may also be implemented as a computer-controlled apparatus, a computer process, a computing system, or an article of manufacture, such as a computer-readable storage medium.

[0009] Other systems, methods, features and/or advantages will be or may become apparent to one with skill in the art upon examination of the following drawings and detailed description. It is intended that all such additional systems, methods, features and/or advantages be included within this description and be protected by the accompanying claims.

BRIEF DESCRIPTION OF THE DRAWINGS

[0010] The components in the drawings are not necessarily to scale relative to each other. Like reference numerals designate corresponding parts throughout the several views.

[0011] Figures 1A-1C show the magnitudes of the ST values for a real chirp signal and the values found by a "positive discretization" and a "symmetric discretization" definition;

[0012] Figure 2A illustrates a 256 x 256 MRI image of a diseased brain with a white cross at pixel P(174, 176);

[0013] Figure 2B illustrates the 256 x 256 MRI image of the diseased brain of Figure 2A rotated by about -42°, with a white cross at pixel P'(134, 68); [0014] Figures 2C and 2D illustrate ST magnitudes at P and P' obtained by positive discretization, with x and y frequency indices on the axes;

[0015] Figures 2E and 2F illustrate radial components at P and P', with a radius on horizontal axis and radial ST magnitude on vertical axis;

[0016] Figure 2G illustrates radial components at P and P' put together for comparison;

[0017] Figures 2H and 2L illustrate the same features a Figures 2C-2G, however using symmetric discretization;

[0018] Figure 3A illustrates 256 x 256 MRI image I with a white cross at pixel P(174,

176);

[0019] Figure 3B illustrates the image of Figure 3A flipped with a white cross at corresponding (81, 176);

[0020] Figure 4A illustrates a 256 x 256 MRI image of a diseased brain, with a white cross at pixel P(174, 176);

[0021] Figure 4B illustrates a 256 x 256 MRI image of the diseased brain rotated by about -42° with a white cross at pixel P'(134, 68);

[0022] Figures 4C and 4D illustrate RIST magnitudes at P and P', with x and y frequency indices on the axes;

[0023] Figures 4E and 4F illustrate radial components at P and P', with a radius on a horizontal axis and radial ST magnitude on vertical axis;

[0024] Figure 4G illustrates radial components at P and P' put together for comparison;

[0025] Figures 4H-4L illustrates the same images as Figures 4C-4G, but using a FTFT- 2D modified for RIST.

[0026] Figure 5A illustrates a 256 x 256 MRI image of a uniform line pattern with sinusoidal intensities inclined at 45° with a white cross at any point P;

[0027] Figure 5B illustrates a 256 x 256 MRI image of the uniform line pattern with same separation inclined at 22.5° with a white cross at any point P';

[0028] Figures 5C and 5D illustrates RIST magnitudes of the two patterns, with x and y frequency indices on the axes;

[0029] Figures 5E and 5F illustrates radial components of the two patterns, with a radius on a horizontal axis and radial ST magnitude on vertical axis; [0030] Figure 5G illustrates radial components of the two patterns put together for comparison;

[0031] Figure 6A illustrates 256 x 256 MRI image of a diseased brain with a white cross at pixel P(174, 176);

[0032] Figure 6B illustrates a 256 x 256 MRI image of the diseased brain rotated by about -42° with a white cross at pixel P'(134, 68);

[0033] Figures 6C and 6D illustrate ST magnitudes at P and P' obtained by RIST* with x and y frequency indices on the axes;

[0034] Figures 6E and 6F illustrate radial components at P and P' from RIST* with a radius on a horizontal axis and radial ST magnitude on vertical axis. ;

[0035] Figures 7A-7D illustrate screenshots of a FTFT-RIST tool for an artificial image of concentric circles. It is rotated by four angles such that a specific pixel (shown as a cross) makes angles 15, 61, 127 and 143 degrees with the x-axis; and

[0036] FIG. 8 is a block diagram of an example computing device.

DETAILED DESCRIPTION

[0037] Unless defined otherwise, all technical and scientific terms used herein have the same meaning as commonly understood by one of ordinary skill in the art. Methods and materials similar or equivalent to those described herein can be used in the practice or testing of the present disclosure. As used in the specification, and in the appended claims, the singular forms "a," "an," "the" include plural referents unless the context clearly dictates otherwise. The term "comprising" and variations thereof as used herein is used synonymously with the term "including" and variations thereof and are open, non-limiting terms. While implementations will be described for performing an S-transform in the context of performing image processing techniques, it will become evident to those skilled in the art that the implementations are not limited thereto.

[0038] 1. Introduction

[0039] Below, the present disclosure describes a variant of a 2D S-transform (ST), called a "Rotationally Invariant S-Transform" (RIST), that is substantially rotationally invariant. Regarding the usage of RIST, while the 2D ST is a complex value, the formula of RIST provides a magnitude (modulus) of the complex number, but not the phases. RIST may be used for square images; as such because most medical images are square or can be made so by cropping and padding the image, RIST has applicability to such images. Moreover, the RIST values obtained by the original formulae are inherently not smooth.

[0040] As determining the RIST value directly may take a long period of time and/or utilize large amounts of memory a FTFT-2D may be used to generate RIST magnitudes for pixels and RIST statistics for regions of interest quickly and accurately. The FTFT-2D algorithm and tools are disclosed in U.S. Provisional Patent Application No. 61/562,486, filed on November 22, 2011, entitled "FTFT-2D Patent Detailed Description," and U.S. Provisional Patent Application No. 61/562,498, filed on November 22, 2011, entitled "FTFT-2D Patent Detailed Description," the disclosures of which are expressly incorporated herein by reference in their entireties.

[0041] RIST magnitudes produced by the FTFT-2D tool may be used for SRF in many medical applications, such as virtual biopsy. Also described herein is another rotationally invariant ST, called RIST*. RIST* may be used in both SFT visualization and spectral analysis. In an implementation, a FTFT-RIST tool displays the values and graphs of RIST* for each pixel or a region of interest (ROI). It also outputs a vector of texture and spectral features based on RIST*.

[0042] Below is a discussion of the algorithms from which the RIST and RIST* are derived.

[0043] 2. Discretization of 1-Dimensional S-Transform

[0044] The 1-dimensional Continuous ST of a complex function of time h(t) is a joint complex function of time t and fre uency/:

[0045] The discrete ST for a signal or time series can be found using the frequency domain, derived by the Convolution Theorem. There are two ways to perform the above, which differ in the summation endpoints.

[0046] A first is as follows:

[0047] A second is as follows:

Here, h[n] = h{n) is the discrete time series and H[k] = H{k/N) is its Fourier Transform, assuming that the sampling interval is 1. The values n and k are the time and frequency indices respectively. The value k is equal to /V/where / is the frequency. Herein, the usage of " [ ]" is for discrete functions of integers, while "( )" is for continuous functions of rea l or complex numbers. In practice, by Nyquist Theorem, the present disclosure seeks to find the ST for/from

0 to 1/2, i.e. for k = 0, 1, N/2 - 1, as there may not be information to find ST for higher frequencies. The following terms, "positive discretization" and "symmetric discretization" are used respectively to signify that the values taken by the summation index are mostly positive in the former and are almost symmetric in the latter.

[0048] Figures 1A-1C show the magnitudes of the ST values for a real chirp signal (Figure 1A) and the values found by the "positive discretization" definition (Figure IB) and the "symmetric discretization" definition (Figure 1C) (shown are ST magnitudes for k = N/2, N -

1 for completeness). I n some instances, Symmetric discretization provides better results than positive discretization, as such estimates of symmetric-discretization ST values are used.

[0049] 3. Discretization of 2-Dimensional S-Transform

[0050] The 2-dimensional Discrete S-Transform (2D ST) of a complex 2-dimensional N x x N y data set or image is a simple extension of ID ST. It is assumed that the intensity

h\n n 1

function x ' y in the image is real. 2D ST is a means of performing SFR. Like ID ST, its frequency-domain formula has two forms: With positive discretization, the following relationship applies:

S P [n x ,n y , k x , k y ]

(4)

whereas with s mmetric discretization, it becomes: ■ (5) [0051] Here, n x , k x , n y , k y are the time and frequency indices respectively in each H\k k 1

direction, and x ' y is the 2-dimensional Fourier Transform. In practice, by Nyquist Theorem, the present disclosure seeks to find the 2D ST for frequency f x and f y from 0 to 1/2, i.e. for k x = 0, 1, N x /2 - 1, and k y = 0, 1, N y /2 - 1. The above equations are applicable when k x , k y are positive.

[0052] Figure 2A illustrates a 256 x 256 MRI image of a diseased brain with a white cross at pixel P(174, 176). Figure 2B illustrates the 256 x 256 MRI image of the diseased brain of Figure 2A rotated by about -42°, with a white cross at pixel P'(134, 68). Figures 2C and 2D illustrate ST magnitudes at P and P' obtained by positive discretization, with x and y frequency indices on the axes. Figures 2E and 2F illustrate radial components at P and P', with a radius on horizontal axis and radial ST magnitude on vertical axis. Figure 2G illustrates radial components at P and P' put together for comparison. Figures 2H and 2L illustrate the same features a Figures 2C-2G, however using symmetric discretization.

[0053] Figures 2C and 2H show the 2D ST magnitudes of a pixel in a 256 x 256 MRI of Figure 2A found using relationship (4) and relationship (5) respectively. Figures 2E and 2J depict the radial component of 2D ST in the k-space. It is only for square images (N = N x = N y ), defined as

In relationship (6), r is the radius in the k-space. | .... | stands for the magnitude of the complex ST value. round() means the nearest integer of a real number. In implementations, the ST magnitudes for those points in the k-space whose magnitudes do not exceed N/2 are considered. Thus, in Figures 2A-2L and all other figures, only a circular sector is displayed. The points outside the sector do not contribute to the texture curve.

[0054] 3. Transformational Invariances of Space Frequency Representation

[0055] For an SFR of a square image several types of transformational invariance may be defined. They are imposed on magnitudes only, as 2-dimensional phases are usually not useful. In practice, it is difficult for any transformational invariance to be satisfied by any SFR exactly (except for reflectional and right-angle rotational invariance of RIST as described in Section 6, below), because of the following concerns: The image is finite with edge effects; the image may not be square (as noted above, it is assumed that the image is square, as in most medical images); the pixel on a rotated image cannot be found that correspond exactly to a given pixel on the original; and a rotated image is a little blurred compared to the original one, due to the interpolation of pixel gray levels during the rotation operation.

[0056] 3.1 Translational I nvariance

[0057] An SFR possesses a "translational invariance" property if the following is true: For any image I and its translation Γ by any vector (u, v), and for any pixel P(n x , n y ) on I and the corresponding pixel P'(n x + u, n y + v) on Γ, the SFR magnitude at every (k x , k y ) in the k-space for P on I is equal to the SFR magnitude at (k x , k y ) for P' on . Tra nslational invariance is well satisfied by most SFR. It is easy to show that ST magnitude is translationally invariant (except for the edge effects), and so for RIST, which is formed in terms of ST.

[0058] 3.2 Rotational Invariance

[0059] An SFR possesses "rotational invariance" property if the following is true: For any image / and its rotation /' by any angle ^ about any point (a, b), and for any pixel P on I and the corresponding pixel P' on , the radial component of SFR magnitudes at any radius r in the k-space for P on I is identical to that for P' on . Thus, given translational invariance, an SFR that is rotationally invariant about a point (a, b) is also rotationally invariant about any other point (α', b').

[0060] I n accordance with the present disclosure, the image in Figures 2A and 2B is rotated and the ST magnitudes are computed at corresponding points using positive and symmetric discretization. The radial component graphs show that neither discretization of ST can satisfy this rotational invariance. Unlike the other invariances defined in Sections 3.1, 3.3 and 3.4, only a weak rotational invariance can be defined in terms of radial components. For instance a strong rotational invariance is as follows: For any image I and its rotation Γ by any angle ^ about any point (a, b), and for any pixel P(n x , n y ) on I and the corresponding pixel P'(n' x , n' y ) on , the SFR magnitude at every (k x , k y ) in the k-space for P on I is equal to the SFR magnitude at (k' x , k' y ) for P' on Γ, where (k' x , k' y ) is the point obtained by rotating (k x , k y ) in the k-space by the same angle Accordingly, below a RIST* is introduced which is an alternative form of RIST, which roughly satisfies this strong requirement.

[0061] 3.3 Reflectional I nvariances

[0062] An SFR possesses a "reflectional invariance about x" property if the following is true: For any image I and its x-reflection I x about any line x = c (with intensity function h i n x -> n y ] h[2c n x ,n y ] ^ ^ P(n x ,n y ) Qn ^ ^ corresponding pixel

P x (2c— n n ) T x

x ' y J on i , the SFR magnitude at every (k x , k y ) in the k-space for P on I is identical to that at same point (k x , k y ) for P x on I x .

[0063] An SFR possesses a "reflectional invariance about y" property if the following is true: For any image I and its y-reflection I Y about any line y = d (with intensity function h ¥ [n x , n v ] = h[n x , 2d - n v ] . . . F(n x ,n v ) . , ,. . .

x x y ), and for any pixel x y ' on I and the corresponding pixel

P Y ? 2d— n ) Y

x ' y ' on I , the SFR magnitude at every (k x , k y ) in the k-space for P on I is identical to that at the same point (k x , k y ) for P Y on I Y . Thus, given translational invariance, an SFR that is reflectionally invariant about a line x = c is also reflectionally invariant about any other line x = c'. Similarly for reflectional invariance about y.

[0064] An SFR possesses a "diagonal reflectional invariance" property if the following is true: For any image I and its reflection I D about the diagonal x = y (with intensity function h D [n x , n v ] = h[n v ,n x ] < . . F(n x ,n v ) . , ,. . . F D (n v , n x ) x y y x ), and for any pixel x y ' on I and the corresponding pixel y x ' on I , the SFR magnitude at every [k x , k y ) in the k-space for P on I is identical to that at the diagonally flipped point (k y , k x ) for P D on I D . Thus, reflectional invariance about x (respectively y) and diagonal reflectional together imply reflectional invariance about y (respectively x).

[0065] 3.4 Right-angle Rotational Invariance

[0066] An SFR possesses a "right-angle rotational invariance" property if the following is true: For any image I and its rotation by ± 90° any point (a, b) and for any pixel P on I and the corresponding pixel P' on , the SFR magnitude at every (k x , k y ) in the k-space for P on I is equal to the SFR magnitude at the diagonally flipped point (k y , k x ) in the k-space for P' on . Thus, given translational invariance, an SFR that is right-angle rotationally invariant about a point (a, b) is also right-angle rotationally invariant about any other point (α', b').

[0067] It is implied by the conjunction of reflectional invariance about x or y, and the diagonal reflectional invariance, because a rotation by +90° is equivalent to diagonal reflection followed by x-reflection, or to y-reflection followed by the diagonal reflection, and similarly for -90°.

[0068] 5. Rotationally Invariant S-Transform (RIST) [0069] It is only defined for square images {N = N x = N y ). For an N x N image I with h\n n 1

intensity x ' y , the 2-dimensional Discrete Rotationally Invariant S-Transform (RIST) magnitude is defined b :

S x \n n k k 1 T x where p x ' y ' x ' y stands for the ST value in positive discretization for the image I obtained by flipping the given image along x, i.e. the intensity in I x is given by

h X [n x ,n y ] = h[N - l - n x ,n y ]

[0070] In the present disclosure, the magnitude of RIST has been defined in terms of the magnitudes of ST, without first defining the complex value of RIST, , itself.

(k k )

As such, relationship (7) can be expressed in words: For each x ' y ' in the k-space, the RIST

V(n n )

magnitude of an image I at pixel x ' y/ is equal to the arithmetic mean of the positive- discretization ST magnitude of the given image at that pixel and that of the flipped image I x at

P x j Y _ l _ n

the corresponding pixel x ' y/ .

[0071] Figures 3A and 3B show the two images with Λ/=256, n * =174 and ^ =176. In particular, Figure 3A illustrates 256 x 256 MRI image I with a white cross at pixel P(174, 176). Figure 3B illustrates the image of Figure 3A flipped with a white cross at corresponding (81, 176). Below, it is shown that relationship (7) satisfies the reflectional (and hence right-angle rotational) invariances defined in Section 3. From examples, relationship (7) roughly satisfies rotational invariance in general.

[0072] Figure 4A illustrates a 256 x 256 MRI image of a diseased brain, with a white cross at pixel P(174, 176). Figure 4B illustrates a 256 x 256 MRI image of the diseased brain rotated by about -42° with a white cross at pixel P'(134, 68). Figures 4C and 4D illustrate RIST magnitudes at P and P', with x and y frequency indices on the axes. Figures 4E and 4F illustrate radial components at P and P', with a radius on a horizontal axis and radial ST magnitude on vertical axis. Figure 4G illustrates radial components at P and P' put together for comparison. Figures 4H-4L illustrates the same images as Figures 4C-4G, but using a FTFT-2D modified for RIST. [0073] Figure 5A illustrates a 256 x 256 MRI image of a uniform line pattern with sinusoidal intensities inclined at 45° with a white cross at any point P. Figure 5B illustrates a 256 x 256 MRI image of the uniform line pattern with same separation inclined at 22.5° with a white cross at any point P'. Figures 5C and 5D illustrates RIST magnitudes of the two patterns, with x and y frequency indices on the axes. Figures 5E and 5F illustrates radial components of the two patterns, with a radius on a horizontal axis and radial ST magnitude on vertical axis. Figure 5G illustrates radial components of the two patterns put together for comparison.

[0074] As can be shown, the same results can be achieved if the image is flipped along y instead of alon x, i.e.,

S p [n x , n v , k x , k v ]

[0075] where F x ' x ' y stands for the positive-discretization ST value for the image I Y obtained by flipping the given image along y, i.e. the intensity in I Y is given by % \n n 1 = h\n N— 1— n 1

x' y x ' y . To prove that relationship (8) is also true, the 2-dimensional

Fourier Transforms of an image I and its flipped counterparts I x , I Y are related by:

H [k x , k y ] = H[-k x , k y ] e

(9) and

ϋπ , ΙΝ

H Y [k x , k y ] = H[k x , -k y ] e

(10) Hence, from relationship (4),

ι2π(- ilnk IN

∑ ∑ H[-K x - k x ,K y + k y e N + N

i2 k x IN

[0076] The last equality comes from the fact that e can be taken outside the double summation and its ma nitude is unity. Similarly:

S P Y [n x ,N -

0077] The right-hand sides of relationships (11) and (12) are equal because their

and therefore relationship (8) is equivalent to relationship (7).

[0078] While there is no rigorous mathematical proof why relationship (7) attains some degree of rotational invariance, experimental results conclude that it is true. However, RIST satisfies right-angle rotational invariance. If the smoothness and small variation of the error function is assumed, then the error should vary from 0 at rotation angle 0° to 0 at angle 90°, through small values at intermediate rotation angles between 0 and 90°. A demonstration of rotational invariance of RIST will be given in Section 9. As RIST magnitude is based on positive discretization, the result is not smooth, as in Figures 4C-4F. The FTFT-2D algorithm has a smoothing effect especially for high frequencies, as seen in Figures 4H-4K.

[0079] 6. Reflectional and Right-angle Rotational I nvariances of RIST

[0080] By relationship (7), RIST satisfies reflectional invariance exactly about the middle line x = (N - l)/2. By the alternative relationship (8), it also satisfies reflectional invariance about the middle line y = (N - l)/2. The diagonal reflectional invariance holds for RIST as well. To this end, x and y can be interchanged for each term on the right-hand side of relationship (7) without changing their values, because the order of double summations in relationship (4) and in the formula for 2-dimensional Fourier Transform inside (4) can be swapped. So relationship (7) becomes: When x, and y are interchanged, the images are reflected about the diagonal, so in relationship Z) oDY

(14), p may be used instead of p . Also, the second term is for reflection along y, so p „_ w

Γ Tid y be used- 0081] From relationship (8), the right-hand side of relationship (14) is exactly , so the proof is complete. As explained in Section 3.4, above, these reflectional invariances imply right-angle rotational invariance about the centre about the centre ((/V - l)/2, (N - l)/2) for RIST.

[0082] By virtue of translational invariance, it can be deduced that RIST has reflectional invariance about any line x = c and about any line y = d, as well as right-angle rotational invariance about any point. But all these are not exact since translational invariance is not.

[0083] 7. FTFT-2D Modified for RIST

[0084] The FTFT-2D may be modified so that they compute RIST values fast and accurately. By "accurately", it is meant that the results obtained are a reasonable approximation of relationship (7). In particular, given a square image \, the flipped image I x is created and both images pre-processed. Then, to find the RIST value for a pixel P in I, the FTFT-

2D algorithm is applied twice, to find the ST magnitude at P in I and that at the corresponding x x (k k )

pixel ? in I , for each x ' y } in the k-space. Finally, the magnitudes are averaged. Like RIST, this modified form of FTFT-2D for RIST satisfies the reflectional (and hence right-angle rotational) invariances. From Figure 4(h)-(l), it roughly satisfies rotational invariance in general.

The time taken to find an RIST value by the new FTFT-2D is slightly more than double that to find ST, but is still very short.

[0085] Figures 4H-4K show the RIST magnitudes found by this new FTFT-2D, and they appear smoother than that by the exact relationship (7), thanks to the cropping operation inside FTFT-2D that helps remove the noise. For example, the process time is 0.039 seconds, while that using the FTFT-2D to find ST magnitude is 0.018 seconds.

[0086] 8. Complex RIST* with Signed Frequencies

[0087] I n accordance with some implementations, an improved form of RIST, called RIST* will now be described. It differs from RIST in several ways. First, it is defined as a complex number, whereas with RIST only a magnitude is defined by relationship (7). Second, it allows the frequency indexes k x and k y to be signed, thus enabling a more comprehensive visualization and analysis of the spectral characteristics of the image. Third, it provides a more convincing demonstration of the rotational invariance of RIST. The RIST* value at a point (n x , n y ) in an N x N square image may be defined as a complex number: \

S P [n x ,n y ,k x ,k y ] if k x >0, k y >0

S p [N-l-n x ,n y ,-k x ,k y ] ifk x <0, k y ≥0

SRisT[ n x> n v> >K] - )

S P Y [n x ,N -l-n y ,k x ,-k y ] ifk x ≥0, k <0

S?[N-l-n x ,N-l-n ,-k x -k ] if k x <0, k <0

(15) where k x and k y can take positive and negative values within N/2 ., N/2 - 1, and p means x-reflection followed by y-reflection.

[0088] Usually the magnitudes of RIST* are sufficient. Thus, only the cases with non- negative k y are needed:

cases with negative k y are redundant because by relationship (13)

and, by replacing P by P in relationship (17):

Hence, only the upper half of RIST* need be drawn.

[0089] Similarly to 2D and RIST, the RIST* magnitudes are of interest for those points in the k-space whose magnitudes do not exceed N/2. So the following figures, only a semicircle is displayed. The points outside the semicircle do not contribute to the texture curve. The texture curve for RIST* is formed in the same way as for RIST, by relationship (6), except that RIST* only averages over the semicircle of radius r, not over the quadrant there. From relationship (7) and relationship (16) that the texture curves of RIST and RIST*, shown in Figures 4E and 6E, are identical, as they are averaging exactly the same quantities. Like RIST, FTFT-2D can be applied to RIST* as well, producing accurate and smoother results in very short time. [0090] Figure 6A illustrates 256 x 256 MRI image of a diseased brain with a white cross at pixel P(174, 176). Figure 6B illustrates a 256 x 256 MRI image of the diseased brain rotated by about -42° with a white cross at pixel P'(134, 68). Figures 6C and 6D illustrate ST magnitudes at P and P' obtained by RIST* with x and y frequency indices on the axes. Figures 6E and 6F illustrate radial components at P and P' from RIST* with a radius on a horizontal axis and radial ST magnitude on vertical axis. Figures 6A-6F show the RIST* magnitudes and texture curves for the example of Figures 4A-4F. As can be seen, the semicircular of RIST* provides more spectral information than the quadrant of RIST. The above may be provided by a FTFT- RIST tool that that displays the magnitude of RIST* in the semicircle, as well as the texture curve of RIST*. It also divides the semicircle into s equal sectors, where s is specified by the user. Then it finds the sector with the largest average RIST* magnitude and draws the texture curve for this major sector.

[0091] Figures 7A-7D illustrate screenshots of a FTFT-RIST tool for an artificial image of concentric circles. It is rotated by four angles such that a specific pixel (shown as a cross) makes angles 15, 61, 127 and 143 degrees with the x-axis. Using s = 90, it determines the major sector at each rotation. As shown, the major sector is at the same four angles. The texture curve in Figure 7(c) represents the radial variation of RIST* in all directions, while the major texture curve in Figure 7(d) is only for that major sector.

[0092] Figure 7 also provides evidence of the strong rotational invariance mentioned in Section 3.2: the semicircular RIST* diagram rotates with the image. As the strong condition implies the weak one, is demonstrates that RIST and RIST* are fairly rotationally invariant.

[0093] 9. Conclusion

[0094] Thus, described above are two methods of formulating a substantially rotationally invariant 2-dimensional discrete SFR based on 2D ST. These new representations, called RIST and RIST*, are different from the discrete ST, but are better in quantifying, visualizing and analyzing localized frequency content in the image. Moreover, the

representations provide a very fast way to compute them for a pixel or for an ROI, using modified forms of the FTFT-2D tool. They are useful for spectral analysis of medical images.

[0095] 10. Exemplary Computing Environment

[0096] FIG. 8 shows an exemplary computing environment in which example embodiments and aspects may be implemented. The computing system environment is only one example of a suitable computing environment and is not intended to suggest any limitation as to the scope of use or functionality.

[0097] N umerous other general purpose or special purpose computing system environments or configurations may be used. Examples of well known computing systems, environments, and/or configurations that may be suitable for use include, but are not limited to, personal computers, server computers, handheld or laptop devices, multiprocessor systems, microprocessor-based systems, network personal computers (PCs), minicomputers, mainframe computers, embedded systems, distributed computing environments that include any of the above systems or devices, and the like.

[0098] Computer-executable instructions, such as program modules, being executed by a computer may be used. Generally, program modules include routines, programs, objects, components, data structures, etc. that perform particular tasks or implement particular abstract data types. Distributed computing environments may be used where tasks are performed by remote processing devices that are linked through a com munications network or other data transmission medium. I n a distributed computing environment, program modules and other data may be located in both local and remote computer storage media including memory storage devices.

[0099] With reference to FIG. 8, an exemplary system for implementing aspects described herein includes an image processing device, such as computing device 800. I n its most basic configuration, computing device 800 typically includes at least one processing unit 802 and memory 804. Depending on the exact configuration and type of computing device, memory 804 may be volatile (such as random access memory (RAM)), non-volatile (such as read-only memory (ROM), flash memory, etc.), or some combination of the two. This most basic configuration is illustrated in FIG. 8 by dashed line 806.

[00100] Computing device 800 may have additional features/functionality. For example, computing device 800 may include additional storage (removable and/or nonremovable) including, but not limited to, magnetic or optical disks or tape. Such additional storage is illustrated in FIG. 8 by removable storage 808 and non-removable storage 810.

[00101] Computing device 800 typically includes a variety of computer readable media. Computer readable media can be any available media that ca n be accessed by device 800 and includes both volatile and non-volatile media, removable and non-removable media. [00102] Computer storage media include volatile and non-volatile, and removable and non-removable media implemented in any method or technology for storage of information such as computer readable instructions, data structures, program modules or other data.

Memory 804, removable storage 808, and non-removable storage 810 are all examples of computer storage media. Computer storage media include, but are not limited to, RAM, ROM, electrically erasable program read-only memory (EEPROM), flash memory or other memory technology, CD-ROM, digital versatile disks (DVD) or other optical storage, magnetic cassettes, magnetic tape, magnetic disk storage or other magnetic storage devices, or any other medium which can be used to store the desired information and which can be accessed by computing device 800. Any such computer storage media may be part of computing device 800.

[00103] Computing device 800 may contain communications connection(s) 812 that allow the device to communicate with other devices. Computing device 800 may also have input device(s) 814 such as a keyboard, mouse, pen, voice input device, touch input device, etc. Output device(s) 816 such as a display, speakers, printer, etc. may also be included. All these devices are well known in the art and need not be discussed at length here.

[00104] It should be understood that the various techniques described herein may be implemented in connection with hardware or software or, where appropriate, with a combination of both. Thus, the methods and apparatus of the presently disclosed subject matter, or certain aspects or portions thereof, may take the form of program code (i.e., instructions) embodied in tangible media, such as floppy diskettes, CD-ROMs, hard drives, or any other machine-readable storage medium wherein, when the program code is loaded into and executed by a machine, such as a computer, the machine becomes an apparatus for practicing the presently disclosed subject matter. In the case of program code execution on programmable computers, the computing device generally includes a processor, a storage medium readable by the processor (including volatile and non-volatile memory and/or storage elements), at least one input device, and at least one output device. One or more programs may implement or utilize the processes described in connection with the presently disclosed subject matter, e.g., through the use of an application programming interface (API), reusable controls, or the like. Such programs may be implemented in a high level procedural or object- oriented programming language to communicate with a computer system. However, the program(s) can be implemented in assembly or machine language, if desired. In any case, the language may be a compiled or interpreted language and it may be combined with hardware implementations.

[00105] Although the subject matter has been described in language specific to structural features and/or methodological acts, it is to be understood that the subject matter defined in the appended claims is not necessarily limited to the specific features or acts described above. Rather, the specific features and acts described above are disclosed as example forms of implementing the claims.