WO2004104611A2

WO2004104611A2 - Mri probe design and tracking, and efficient mri reconstruction and deblurring

Info

Publication number: WO2004104611A2
Application number: PCT/US2004/014128
Authority: WO
Inventors: Eddy Y. Wong; Claudia M. Hillenbrand; Jonathan S. Lewin; Jeffrey L. Duerk; Daniel R. Elgort; Chris Flask; Frank K. Wacker; Brian M. Dale; Hisamoto Moriguchi
Original assignee: Case Western Reserve University
Priority date: 2003-05-05
Filing date: 2004-05-05
Publication date: 2004-12-02
Also published as: EP1623245A2; WO2004104611A3

Abstract

A probe suitable for attachment to, or incorporation in, a medical interventional device, such as a catheter, and which may be employed for tracking, imaging, or both, includes a first material having an MR resonance frequency distinct from a resonance frequency of a second material adjacent to the first material.

Description

TITLE: MRI PROBE DESIGN AND TRACKING, AND EFFICIENT MRI RECONSTRUCTION AND DEBLURRING

CROSS-REFERENCE TO RELATED APPLICATIONS

[0001] This application incorporates by reference in entirety, and claims priority to and benefit of, the U.S. Provisional Patent Applications having the following numbers and respective filing dates: 60/468,172 (filed on 05 May 2003); 60/468,173 (filed on 05 May 2003); 60/468,177 (filed on 05 May 2003); 60/483,219 (filed on 27 Jun 2003); and 60/485,823 (filed on 08 Jul 2003).

BACKGROUND

[0002] A goal of endovascular MRI-guided interventions is the combination of MRI's diagnostic capabilities (e.g., angiography, morphology, plaque analysis, perfusion imaging and others) with therapeutic interventions such as angioplasty, catherectomy and stent placement. A successful MRI-guided endovascular therapeutic procedure incorporates a subset of the following steps: MR guidance of the interventional device to the target region, high-resolution imaging at the target location in order to diagnose disease within the vessel wall, performance of a therapeutic intervention, and evaluation of the efficacy of therapy. The requirements for successful guidance and high-resolution imaging are generally quite different. To date, intravascular MRI devices have been designed primarily for either active and passive tracking or for high-resolution intravascular imaging. Minimally invasive image-guided therapy requires a quick and robust method for device localization, as well as the ability for tracking and guidance of interventional devices inserted into the body. [0001] In a typical magnetic resonance imaging (MRI) system, a subject such as a human body is placed in a static magnetic field that orients the proton magnetic dipoles. A field gradient is imposed along the z-axis in the direction of the main , magnetic field such that a narrow plane of protons resonate within a band of frequencies. A phase encoding gradient along the x-axis is activated for a short time during which the dipoles acquire a different phase. A frequency encoding gradient along the y-axis is then activated to frequency encode the positions of the dipoles while a receiver coil is activated to record the signal. Typically, 128, 256, or 512 data points are recorded along the frequency encoded axis. As each recorded data point corresponds to a respective pixel in the image to be generated, the number of recorded data points determines the resolution of that image.

[0002] A Fourier Transform (FT) algorithm is used to decode frequency information contained in a proton signal at each location in the imaged plane to corresponding intensity levels. The resulting image is then displayed as shades of gray in a matrix arrangement of pixels. MRI may be used to reconstruct images for any standard orientation, such as transverse, coronal, or sagittal slices for example, or for any oblique orientation. MRI can also be used for interventional procedures, such as guiding medical devices through vessels and placing such devices inside vessels.

[0003] Interventionalists can perform procedures efficiently and safely if provided access to relevant image information. For example, an interventionalist may require information about a current position/orientation of a device, such as a catheter, in the body, which could be overlaid with a three-dimensional (3D) reference map showing the catheter's current position in real time. Real-time imaging can also overcome problems associated with patient movement.

[0004] In imaging most tissues with MRI, the hydrogen protons from water are preferably detected as most soft tissues are composed of greater than approximately eighty percent water. A device with a different resonance frequency can be attached to or incorporated in the catheter, with software in the scanner alternating between localizing the catheter and collecting image data.

[0005] Adjusting the MR image parameters (e.g., slice position, orientation, resolution, TE, TR, etc.) during an interventional catheter-based procedure can be a cumbersome process that may require the interventionalist or technologist to leave the magnet room and use a keyboard and a mouse in combination with a graphical user interface. Other conventional adaptive image parameter systems verify the feasibility of using feedback based on a catheter's insertion speed to adjust, in realtime, the value of specific adaptive image parameters. These systems, however, provide little clinical utility because they are neither sufficiently flexible nor robust for use in intravascular MR guided procedures. [0003] It would therefore be desirable to provide a system and process that automatically and continuously adjusts specific image parameters in real time, based on a catheter's speed of insertion.

[0004] In one aspect, the present invention relates to magnetic resonance imaging ("MRI"). It finds particular application in conjunction with reconstruction and deblurring of MRI images.

[0005] Magnetic resonance imaging is a diagnostic imaging modality that does not rely on ionizing radiation. Instead, it uses strong (ideally) static magnetic fields, radio-frequency ("RF") pulses of energy and magnetic field gradient waveforms. More specifically, MR imaging is a non-invasive procedure that uses nuclear magnetization and radio waves for producing internal pictures of a subject. Three- dimensional diagnostic image data is acquired for respective "slices" of an area of the subject under investigation. These slices of data typically provide structural detail having a resolution of one (1 ) millimeter or better.

[0006] Programmed steps for collecting data, which is used to generate the slices of the diagnostic image, are known as an MR image pulse sequence. The MR image pulse sequence includes magnetic field gradient waveforms, applied along three axes, and one or more RF pulses of energy. The set of gradient waveforms and RF pulses are repeated a number of times to collect sufficient data to reconstruct the slices of the image.

[0007] For image reconstruction, the collected k-space data are typically reconstructed by performing an inverse Fourier transform (IFT). However, in certain experimental settings, such as spiral acquisition techniques and non-rectilinearly sampled data, image reconstruction is not simple and artifacts, such as blurring due to off-resonance effects have to be corrected. In addition, a large number of 2D- FFTs have to be performed if the data set is large, which may cause impractical and unacceptable delays in image processing.

[0008] It would therefore be desirable to provide more efficient methods for image reconstruction that render MR images with an image quality that is practically indistinguishable from that obtained with conventional image reconstruction methods, such as gridding algorithms and frequency-segmented off-resonance correction method used in spiral imaging.

SUMMARY OF THE INVENTION

[0009] There is a continuing need for an improved dedicated tracking device, in particular in conjunction with active tracking, and of a dedicated imaging device. There also is a need for a high-performance device suitable for both tracking and imaging, and for an improved tracking device that reduces or eliminates tissue damage caused by local heating of wired probes. Accordingly, the MRI tracking and imaging systems described herein are directed at addressing these needs. The systems and methods of the invention are directed, in at least one aspect, at diagnostic and therapeutic MRI imaging applications. More particularly, the systems and methods disclosed herein are directed at designs of coils and wireless probes that can be used for tracking and imaging applications.

[0010] According to one aspect of the invention, an MR imaging coil for intravascular imaging includes a coil assembly comprising coaxially arranged first and second coils, wherein the first coil is wound in a first direction and the second coil is at a distance from the first coil and wound in a direction opposite of the first direction. In one configuration, the connecting wires connecting the first coil with the ^• second coil are disposed either longitudinally near or at the center between the first and second coils, or substantially collinearly longitudinally along an imaginary line connecting a periphery of the first and second coils. These coil arrangements provide high image quality with a low number of artifacts introduced by the connection wires.

[0011] In an alternate embodiment, each coil can be connected to a separate receive channel that can be accessed independently during an MR imaging experiment. According to one practice, each coil can be individually frequency- tuned, and the signals received from each coil can be processed separately and used to localize the imaging coil and/or to image the vessel walls.

[0012] According to another aspect of the invention, a wireless (or coil-less) MRI tracking probe includes a lumen containing a first material having a magnetic resonance (MR) frequency distinct from a resonance frequency of a second material adjacent to the lumen. The first material may include a liquid, such as acetic acid and a contrast material which may include, for example, a rare-earth compound, such as gadolinium. The second material can include anatomic tissue, producing, for example, a fat signal, a water signal, or a combination thereof.

[0013] The probe described above can be incorporated in a catheter. The probe and/or the catheter can be tracked by applying a conventional FLASH (Fast Low Angle SHot) sequence to excite the magnetic resonance (MR) frequency for obtaining an anatomic image of the second material, and by applying a CHESS (Chemical shift selective excitation) pulse to excite the MR resonance frequency of the first material for obtaining a catheter-selective image. A suitable image of the probe in the anatomic tissue can be obtained by comparing the anatomic image and the catheter-selective image. Repeating the FLASH and CHESS-based image acquisition at substantially discrete temporal instances produces a sequence of catheter-selective and anatomic images that facilitate tracking of the catheter over a length of time. [0006] In one embodiment, the invention is directed to a system and method that, inter alia, incorporates real-time imaging sequences, flexible catheter tracking methods, adaptive parameter modes, and a user-friendly interface for the interventional physician.

[0007] This adaptive tracking system uses, in various embodiments, real-time tracking techniques to continually monitor a catheter tip's 3D position (including, for example, position relative to a target anatomy or position relative to the MR imager's receive coils), orientation, insertion speed, and a combination of physiological parameters, such as, without limitation, breathing rate, heart rate, etc. The device position and orientation information is used to automatically adjust the scan plane for real-time imaging. Insertion speed may be used to automatically adjust pre-specified adaptive image parameters in real-time. Image resolution, FOV, Bandwidth, TE, TR, and temporal resolution are employed as adaptive parameters; many other image parameters can also be used.

[0008] The systems and methods described herein represent a "hands free" interventional MRI systems and methods that automatically and adaptively adjust imaging parameters and react to changing clinical requirements in real-time. The systems and methods disclosed herein provide the interventional radiologist with a means to effect a change in specific image parameters, in real-time, using the catheter itself.

[0009] According to one embodiment, the invention is directed at a method of adaptively adjusting at least one MR imaging parameter for an interventional procedure. The method includes (a) adaptively tracking an MR-guided probe inserted into an object by (i) locating the probe by acquiring first probe coordinates with respect to a reference coordinate system, and (ii) calculating a velocity of the probe relative to the object by acquiring second probe coordinates with respect to the reference coordinate system (typically the MR scanner/imager's coordinate system, which may or may not be an orthogonal coordinate system), and (b) based at least partially on the calculated velocity, adjusting a subset of the at least one MR imaging parameter to adaptively track the probe, image a target region of the object, or both.

[0010] According to one practice, acquiring the first probe coordinates includes acquiring a plurality of one-dimensional frequency-encoded projections to determine at least one of a three-dimensional position of the probe and an orientation of the probe. Various imaging parameters that may be adjusted by the systems and methods disclosed herein include a subset of: field of view, image spatial resolution, image scan plane position, scan plane orientation, temporal resolution, bandwidth, slice thickness, imaging pulse sequence, image contrast, TR, TE, active receiver channels, k-space trajectory, excitation flip angle, MR scanner table position (e.g., to keep the probe proximal to an isocenter of the table).

[0011] According to one aspect, adjusting the one or more imaging parameters is at least partially based on an auxiliary parameter. The parameter includes an element belonging to a subset of: position of the probe relative to a target region of the object, position of the probe relative to the MR imager's receive coils, probe orientation, a physiological parameter associated with the object, and a combination these.

[0014] According to another embodiment, the invention is directed at A method of adaptively adjusting at least one MR imaging parameter for an interventional procedure. The method includes (a) locating an MR-guided probe inserted into an object by acquiring first probe coordinates with respect to a reference coordinate system, and (b) based at least partially on a parameter associated with the located probe, adjusting a subset of the at least one MR imaging parameter to adaptively to track the probe, to image a target region of the object, or both.

[0002] In one aspect, the invention is directed to methods for efficient reconstruction and deblurring of MRI data obtained with non-rectilinearly acquired k-space data using, for example, spiral imaging.

[0003] According to one aspect, the invention is directed at a method of reconstructing a magnetic resonance image from non-rectilinearly-sampled k-space data. The method includes, (a) distributing the sampled k-space data on a rectilinear k-space grid, (b) inverse Fourier transforming the distributed data, (c) setting to zero a selected portion of the inverse-transformed data, (d) Fourier transforming the zeroed and remaining portions of the inverse-transformed data, at grid points associated with the selected portion, (e) replacing the Fourier-transformed data with the distributed k-space data at corresponding points of the rectilinear k-space grid, thereby producing a grid of updated data, (f) inverse Fourier transforming the updated data, and (g) applying an iteration of steps b through f to the inverse Fourier-transformed updated data until a difference between the inverse Fourier- transformed updated data and the inverse Fourier-transformed distributed data is sufficiently small.

[0004] Accord ing to one practice, a 2D-IFT is performed on a large rescaled matrix after k-space data are distributed without density compensation. After the initial image is reconstructed, all the matrix elements except the central Λ/ x N region are replaced by zeros. This matrix is 2D-Fourier transformed, wherein the obtained matrix is the result of a convolution of the rescaled matrix with a 2D sine function. The matrix coordinates where the original data exist, i.e., the non-zero element locations in the rescaled matrix, the data are replaced by the original data values, and a 2D-IFT is subsequently performed on this matrix, leading to the updated reconstructed image. These procedures are repeated until the difference between the updated reconstructed image and the image at the previous iteration becomes sufficiently small.

[0005] In one practice, the aforedescribed process can be improved by applying consecutively increasing scaling factors. In an alternative practice, the acquired k- space region can be partitioned into several blocks, and the aforementioned algorithms can be applied to each block.

[0006] According to another aspect of the invention, a Fourier transform technique is used for variable density spiral (VDS) sampling. VDS trajectories acquire the central region of k-space data so that the Nyquist criterion is satisfied. However, outer k- space regions are undersampled to reduce the acquisition time. The original target grid is an N x N matrix. K-space data are distributed into a larger matrix rescaled by a factor of s. The location of each datum in the large rescaled matrix is determined by multiplying the original k-space coordinate by s. An Inverse Fourier Transform (I FT) is performed on the rescaled matrix, leading to an image matrix, with the intermediate reconstructed image appearing in the center of the N x N matrix, on which region a phase constraint is imposed. The region outside of the central N x N matrix is set to zeros, and a FT is performed on the matrix which includes the zeros, leading to a transformed matrix which is an estimate of the phase constrained raw data. Then, a data-consistency constraint is imposed on this transformed data matrix, i.e., the data where the original data exist are replaced by the original data values, and an IFT is subsequently performed, with an updated reconstructed image appearing in the central N x N matrix. The procedures are repeated until the difference between the updated image and the image at the previous iteration becomes sufficiently small.

[0007] According to one aspect, the invention is directed at a method of reconstructing a magnetic resonance image from non-rectilinearly-sampled k-space data. The method includes (a) distributing the sampled k-space data on a rectilinear k-space grid, (b) convolving the distributed data with a sine function, (c) at least partially based on a characteristic of the sine function, replacing a portion of the convolved data with a corresponding portion of the k-space data distributed on the rectilinear k-space grid, thereby producing a grid of updated data, and (d) applying an iteration of steps b through c to the updated data until a difference between the updated data and the distributed data is sufficiently small.

[0008] According to one practice, the invention is directed at a method of reconstructing a magnetic resonance image from non-rectilinearly-sampled k-space data. The method includes (a) distributing the sampled k-space data on a rectilinear k-space grid, (b) partitioning the k-space grid into blocks, (c) inverse Fourier transforming distributed data of at least one of the blocks, (d) setting to zero a selected portion of the inverse-transformed data in the at least one of the blocks, (e) Fourier transforming the zeroed and remaining portions of the inverse Fourier- transformed data, (f) at grid points associated with the selected portion, replacing the Fourier-transformed data with the distributed k-space data at corresponding points of the rectilinear k-space grid, thereby producing a grid of updated block data, (g) inverse Fourier transforming the updated block data; and (h) applying an iteration of steps b through g to the inverse Fourier-transformed updated block data until a difference between the inverse Fourier-transformed updated block data and corresponding inverse Fourier-transformed distributed block data is sufficiently small.

[0009] According to one practice, in a method for reconstructing deblurred spiral images, off-resonance correction proceeds block-by-block through a reconstructed image, and FFTs are performed on matrices (M M) that are smaller than the full image matrix (N x N, with M<N). A small block region Mx M is extracted and a 2D- FFT is performed on the M M image matrix. After frequency demodulation, with the mean off-resonance frequency of the central rM x rM pixels (0 < r ≤ 1 , and r is typically 0.5.) in the Mx M phase image matrix used as demodulation frequency, the M x M k-space data is subsequently 2D-inverse Fourier transformed. Only the central rM x rM pixels of the Mx M deblurred image matrix are kept for the final reconstructed image to remove artifacts. This procedure is repeated until the entire scanned object is deblurred.

[0010] Accord ing to another aspect, the invention is directed at a method of reconstructing a magnetic resonance image from non-rectilinearly-sampled k-space data. The method includes (a) distributing the sampled k-space data on a rectilinear k-space grid, (b) partitioning the k-space grid into blocks, (c) convolving distributed data of at least one of the blocks with a sine function, (d) at least partially based on a characteristic of the sine function, replacing a portion of the convolved block data with a corresponding portion of the k-space data distributed on the rectilinear k- space grid, thereby producing a grid of updated block data, and (e) applying an iteration of steps b through d to the updated block data until a difference between the updated block data and corresponding distributed block data is sufficiently small.

[0011] Further features and advantages of the present invention will be apparent from the following description of preferred embodiments and from the claims.

BRIEF DESCRIPTION OF THE DRAWINGS

[0001] The following figures depict certain illustrative embodiments of the invention in which like reference numerals refer to like elements. These depicted embodiments are to be understood as illustrative of the invention and not as limiting in any way.

[0002] Fig. 101 shows a single-element (A) and double-element (B) tracking antenna with and without internal signal source and (C) a single-element tracking antenna mounted on a balloon catheter;

[0003] Fig. 102 depicts a bi-plane radial localization algorithm;

[0004] Fig. 103 depicts vessel phantom imaging with a two-element tracking antenna having an internal signal source;

[0005] Fig. 104 shows a temporal sequence of in vivo porcine imaging experiments with a two-element tracking antenna without a signal source;

[0006] Fig. 105 shows an approximation of a solenoid using a single square and octagon winding;

[0007] Fig. 106 (A-E) show five different embodiments of opposed-solenoid imaging coils;

[0008] Fig. 107 depicts uniform phantom imaging studies conducted with the five opposed-solenoid imaging coils of FIGs. 6A-6E;

[0009] Fig. 108 shows Biot-Savart simulation results obtained for different placements of wires connecting opposed solenoid windings;

[0010] Fig. 109 (A-D) show Biot-Savart simulation results for different coil construction parameters of the coils depicted in FIGS. 6A-6E;

[0011] Fig. 110 (A-D) show in vivo and in situ images obtained from a porcine imaging experiment;

[0012] Fig. 111 shows the opposed solenoid array (a) and a single-channel opposed solenoid antenna (b);

[0013] Fig. 112 shows schematically the implemented device localization method based on two-channel acquisition;

[0014] Fig. 113 shows axial and sagittal coil sensitivity profiles in a vessel phantom, wherein

(a) and (c) show sagittal slices acquired simultanously with both elements;

(b) and (d) show the corresponding signal intensity profiles along the marked lines;

(e) shows the combined image calculated from (a) and (c),

(g) shows the corresponding slice acquired with the conventional opposed solenoid coil of the same size;

(f) and (h) provide axial slices through the center of the device for an array coil and a conventional coil;

[0015] Fig. 114 shows porcine imaging tracking experiments in vivo using the opposed solenoid phased array coil;

[0016] Fig. 115 shows In vivo TrueFISP imaging from inside the abdominal artery of a pig;

[0017] Fig. 116 shows In situ TrueFISP imaging of the abdominal aorta using the catheter array coil;

[0018] Fig. 117 shows vascular phantom images, (a, b) Axial and sagittal slice- selective FLASH images with acetic acid catheter (arrows), (c, d) Axial and sagittal CS-FLASH images with selectively excited catheter clearly visible; and [0019] Fig. 118 shows human volunteer images, (a, b) Axial and sagittal slice- selective FLASH images with acetic acid catheter (arrows), (c, d) Axial and sagittal CS-FLASH images with catheter selectively excited in projection images.

[0012] FIG. 201 A is a functional diagram of an exemplary adaptive tracking system architecture;

[0013] FIG. 201 B is a functional diagram of an exemplary adaptive tracking system architecture, wherein a 3D MIP roadmap is acquired in addition to real-time 2D images of the catheter and the surrounding tissue;

[0014] FIG. 202 shows one embodiment of a tuned resonance circuit mounted on a catheter;

[0015] FIG. 203 shows a second embodiment with two tuned resonance circuits mounted on a catheter;

[0016] ^■ FIG. 204 shows an exemplary binary function for adjusting image parameters; [0017] FIG. 205 shows an exemplary continuous function for adjusting image parameters;

[0018] FIG. 206 shows an exemplary user interface for real-time updating of parameters;

[0019] FIG. 207 shows a temporal sequence of image frames using a phantom experiment; and

[0020] FIG. 208 shows a temporal sequence of image frames using a porcine experiment.

[0007] Fig. 301 is a flow chart of a basic INNG algorithm;

[0008] Fig. 302 is a flow chart of a facilitated INNG algorithm;

[0009] Fig. 303 is an exemplary partition scheme of a BINNG algorithm;

[0010] Fig. 304 is a numerical phantom used in computer simulation experiments; [0011 ] Fig. 305 shows profiles at the 34th row and 102nd column of the reconstructed images, wherein

[0012] Fig. 305(a) is the basic INNG algorithm with s = 8 after 30 iterations;

[0013] Fig. 305(b) is the basic INNG algorithm with s = 8 after 101 iterations;

[0014] Fig. 305(c) is the facilitated INNG algorithm after a total of 30 iterations;

[0015] Fig. 305(d) is the BINNG algorithm with (d₂, of₄, d₈) = (0.01 , 0.01 , 0.001 );

[001.6] Fig. 305(e) is the conventional gridding algorithm with Voronoi DCF; and

[0017] Fig. 305(f) is the BURS algorithm with p = 0.8;

[0018] Fig. 306 shows the RMS error of the image with ideal data at each iteration;

[0019] Fig. 307 shows the measured SNR of the image with noisy data at each iteration;

[0020] Fig. 308 shows reconstructed in-vivo brain images;

[0021] Fig. 308(a) shows the facilitated INNG algorithm after 25 total iterations (6 and 19 iterations for s = 2 and 4, respectively).

[0022] Fig. 308(b) shows the BINNG algorithm with (d₂, d , d₈) = (0.01 , 0.01 , 0.003);

[0023] Fig. 308(c) shows the conventional gridding algorithm with Voronoi DCF;

I [0024] Fig. 308(d) shows the BURS algorithm with p = 0.8;

[0025] Fig. 309 shows reconstructed images: (a) standard reconstruction from reduced data, (b) the PFSR algorithm from reduced data, and (c) the standard reconstruction from full data;

[0026] Fig. 310 shows an exemplary block diagram of Block Regional Off- Resonance Correction (BRORC); [0027] Fig. 311 shows axial brain images;

[0028] Fig. 311(a) is an image before off-resonance correction;

[0029] Fig. 311 (b) is the image using FSORC with L = 12;

[0030] Fig. 311(c) is the image using BRORC with (M, r) = (8, 0.25);

[0031] Fig. 311 (d) is the image using BRORC with (M, r) = (32, 0.5);

[0032] Fig. 312 shows heart images, with a 160 x 160 matrix centered on the heart being cropped;

[0033] Fig. 312(a) is the image before off-resonance correction;

[0034] Fig. 312(b) is the image using FSORC with L = 15; and

[0021] Fig. 312(c) is the image using BRORC with (M, r) = (32, 0.5).

DETAILED DESCRIPTION OF CERTAIN ILLUSTRATED EMBODIMENTS

[0022] The invention generally is directed at systems and methods for MRI tracking and imaging applications. In particular, the systems and methods described herein can be both actively and passively tracked and may also be suitable for providing high-resolution images.

[0023] Various implementations of tracking applications will now be described. Several different catheter designs with tracking coils and/or antennas are considered: each varies by the number and arrangement of active loop elements and whether or not an internal signal source is incorporated into the device. A coil-less tracking device is also disclosed.

[0024] Referring to Fig. 101 , a single-element (A) and double-element (B) tracking antenna with and without internal signal source and (C) a single-element tracking antenna was mounted on a balloon catheter. In one embodiment, active loop elements are wound from a copper magnet wire of about 30 AWG. Dimensions of the loop elements are approximately 4 mm along the long axis and 2.5 mm along the short axis. The loops of the double-loop tracking coils were wound with a center-to- center distance of about 23mm. Tuning and matching of the resonant circuit was accomplished using surface mount capacitors. Capacitive coupling to the MR receiver system was made utilizing a micro-coaxial cable. A plastic tube was then affixed over the active antenna elements and secured into place with epoxy. The tube was then filled with an internal signal source and sealed. Devices utilizing a signal source have a maximum dimension of about 11 F in diameter while those that do not have a signal source have a maximum diameter of about 8F. A single- element device is shown in Fig. 101 A and a two-element device in Fig. 101B. A single-element tracking antenna with no signal source was also mounted on a balloon catheter. This device is shown in FIG. 101C.

[0025] The use of a loop coil design provides for sensitivity in most orientations during typical device use. An internal signal source provides several advantages in active tracking experiments. The signal source allows for tracking to be performed with very low flip angles (approx. 1° - 2°). Tip angle amplification results in the signal source seeing an effective tip angle greater than the rest of the surrounding tissue, making identification from the surrounding tissue easier. The use of a low tip angle also eliminates the need for dephaser gradients to eliminate tissue signal during tracking, allowing for increased temporal resolution in tracking experiments.

[0026] Phantom and in vivo porcine tracking experiments were conducted using a Siemens 1.5T Sonata clinical scanner. Tracking experiments examined the reliability and robustness of each antenna design for detection and tracking. Evaluation criteria included the ability to reliably follow catheter insertion and retraction in both a vessel phantom and in an animal model (e.g., anatomic vessel, tissue, etc.). Accuracy measurements were also performed to evaluate the ability to properly position the scan plane. Software implementing a subset of the aspects of the systems and methods described herein was employed to cause the scanner to automatically track the catheter in real-time using a limited number of projections. The system also allows an imaging slice location and orientation to follow the catheter by alternating between localization and imaging modes. For use with the single-element catheter, the software collects three projections (orthogonal or non- orthogonal) and updates the scan plane position. For use with the double-element catheter, the tracking software implements a bi-plane radial localization algorithm and updates the scan plane position and orientation. This process is depicted pictorially in Fig. 102. The tracking software is combined with several fast imaging sequences to collect image data between localizations. The tracking experiments performed using the systems and methods described herein may utilize both fast, low long-angle (FLASH) sequences and steady-state (True-FISP) sequences.

[0027] Catheter advancement and retraction were tracked in all imaging experiments. The contrast-to-noise ratio (CNR) of the internal signal source was approximately 10 during the tracking phase of the experiments; the system accuracy was better than 3mm in displacement error and 2° orientation error.

[0028] The use of a loop antenna design is advantageous for several reasons. In catheter-based interventions on a horizontal B0 field MR system, the loop antenna design provides the greatest B1 sensitivity as compared with other designs such as the solenoid. Additionally, the loop design conforms well to the shape of the catheter with minimal increase in dimensions; moreover, a sheath for the internal signal source can be easily applied over the antenna. An internal signal source can provide several advantages in active tracking experiments. When used for tracking applications inside blood vessels, however, the internal source does not offer a distinct advantage as the blood itself serves as a signal source close to the antenna elements. However, in interventions that occur in locations where there may be a lack of excitable spins near the antenna element, such as in the digestive tract or lung, an internal signal source is advantageous.

[0029] As mentioned above, in various embodiments, configurations with single- and double-active loop elements are used. The advantage of a single-loop design is that the localization system is more robust, generally not having to extract more than one signal peak from the tracking data. With respect to speed, tracking can be performed faster on a single loop system than tracking with a two-antenna system. However, there exists a distinct advantage in having a two-loop design, as a two- loop device allows the tracking system to provide both device orientation and device position.

[0030] Fig. 103 depicts vessel phantom imaging with a two-element tracking antenna having an internal signal source. Fig. 104 shows a temporal sequence of in vivo porcine imaging experiments performed with a two-element tracking antenna without a signal source.

[0031] Whereas the aforedescribed coil design is well-suited for tracking, a coil design with good radial homogeneity is advantageous for applications aimed at characterization of the vessel wall or locating vulnerable plaque. Due to the complexity and the generally lengthy process required for construction of a catheter- based imaging coil, simulations are performed to determine the coil parameters that provide the optimal SNR and field homogeneity.

[0032] A common way to quantitatively characterize the magnetic field produced by a distribution of steady-state currents is by employing the Biot-Savart law, a derivation from Maxwell's equations. For a segment of wire of length dl carrying a current I, the Biot-Savart law states that the magnetic field dB produced by that wire segment at a distance r from the segment is given by the following equation: μ_QI disυxQ dB = 4π r [0033] Solving this for an arbitrary configuration of current elements, such as with a coil, the magnetic field vector at any point in space can be calculated. To calculate B, the equation can be integrated as follows:

[0034] For a wire element, carrying current I, the integration leads to the following formula for B:

B =

4^r where r is the radial distance from the wire element at which B is to be calculated.

[0035] From these equations, it is possible to calculate the magnetic field for a given configuration of wires.

[0036] Biot-Savart simulations were performed using MATLAB (The Mathworks, Natick, Massachusetts). B-field components were calculated for prescribed imaging planes with respect to the location of the coil. Separate files were created for each opposed solenoid coil configuration to be examined. Solenoid coils were approximated using octagons. This was done as a compromise between an accurate approximation of a solenoid and the ease by which the variations in each parameter to be optimized could be incorporated. Octagons were deemed sufficient for representing a solenoid through inspection of a field plot generated with both octagon and square winding, as shown in Fig. 105. Field lines were obtained for a distance from the windings equivalent to the location of the enhancement region of an opposed solenoid coil. A circle was inscribed in the field plots and distortions were noted.

[0037] Two areas of coil design and their effect on SNR and field homogeneity were examined. The placement of the electrical wires to connect the two opposed solenoids and their effect on the field homogeneity was first investigated. Here, five different configurations for wire location were examined. The first two configurations involved placing the wires on opposite sides of the interventional device. In the first configuration, an additional half turn was added to one of the solenoids so that the connecting wires would be located on opposite sides of the interventional device, while in the second configuration, the additional half-turn is eliminated by feeding the wire straight through the interventional device. The third configuration involved running the wires along the inside wall of the interventional device. The fourth configuration involved running the wires inside the device but substantially along the center instead of the inside walls and the fifth configuration involved running the wires along the same side of the device. For reference, a simulation was generated for an opposed solenoid coil with no connecting wires. [0038] The second area of simulations examined different solenoid construction parameters for their effect on the B-field in the imaging region of an opposed solenoid coil. Four parameters were examined; the number of windings in each solenoid, diameter of the solenoid, the pitch of the solenoid windings and the separation distance between the opposed solenoid coils. The number of windings was varied between one and nine, the diameter of the solenoid was varied between about 3F and about 18F, the pitch of the solenoid windings was varied between about zero and about 2.5 times the wire thickness and the separation distance was varied between about zero and about five times the diameter using a solenoid with a diameter of about 5F.

[0039] The Biot-Savart simulation results obtained for the opposed solenoid coil should follow closely those trends prescribed by the two equations defining the B- field of a conventional solenoid. Inside a solenoid, the equation for calculating the transverse magnetization, Bxy, is given as:

B xy -ώ _LH

[0040] From this equation, a direct relationship between the number of turns (N) and Bxy is observed while an inverse relationship between the pitch of the windings, which relates to overall length (L), and Bxy is also observed. For calculating the field outside of the solenoid, which examines the radial "reach" of an imaging coil, the following equation applies:

2

B_xv(z) ^ A,IR

2z³

[0041] From this equation, a direct relationship between the diameter of the coil and B_Xy at a given distance away from the coil can be observed. Simulation results obtained follow the relationships as described by the equations.

[0042] Five opposed solenoid coils, each with different variants in the four simulation parameters, were constructed for validation of the simulation results. The five coils constructed included one with an optimized coil design (FIG. 106A), one with a large separation distance between the opposed solenoid windings (FIG. 106B), one with zero pitch spacing (FIG. 106C), one with twice the number of windings as the optimized coil design (FIG. 106D) and one with a larger diameter than that of the optimized coil design (FIG. 106E). The coils were constructed on plastic formers that measured about 9 French or about 14 French in diameter. Copper wire of about 30 AWG was used to wind the solenoid elements. The opposed solenoids were tuned to 63.6 MHz and matched to 50Ω using variable capacitors. Active detuning was provided via a choke and PIN diode.

[0043] Loaded Q measurements were obtained for each coil using a reflection-type measurement from Sn polar plots. Sn measurements were obtained from an HP 3577A Network Analyzer with S-parameter test set (Hewlett Packard, Palo Alto, CA). Measurements were performed with the coils placed in the lumen of a uniform saline phantom.

[0044] Imaging experiments were conducted using the same saline phantom on the Siemens Magnetom Sonata 1.5T whole body clinical imager (Siemens Medical Solutions, Erlangen, Germany) described above. Small diameter transmission lines and connectors were used to link the opposed solenoid imaging coil to a Siemens flex-loop interface which, in turn, interfaced with two of the eight RF receivers. The RF receiver channels on the MR system have a ±250kHz maximum bandwidth. The MR system has 40 mT/m maximum amplitude gradients and a minimum gradient slew rate of 200mT/m/ms. RF excitation was performed using the standard body coil. The system was configured via the standard user interface. The saline phantom was placed substantially in the center of the magnet and the imaging coils were placed substantially horizontally, along Bo, substantially in the center of the saline phantom. A Segmented EPI sequence was used for imaging (TR/TE = 1520/38 ms, α= 90°, SL = 3 mm, Matrix = 128x128, FOV = 40mm2). Line intensity plots were taken from the image and grayscale intensities were compared for each of the different coils (FIG. 107).

[0045] An optimized, opposed-solenoid imaging coil was created and mounted on a 5-French catheter. Each solenoid included 5 windings of 30 AWG copper wire with a pitch spacing of approximately one wire diameter. The individual counter-wound solenoid coils were placed about 10 mm apart. The antenna was tuned and matched using surface mount capacitors, and a micro-coaxial cable was utilized to provide capacitive coupling to the MR receiver. Porcine imaging experiments were also conducted on the 1.5 T Siemens Sonata imager. The catheter was placed in both the vena cava and the iliac artery, and imaging was performed both in vivo and in situ. True FISP imaging techniques that have been previously defined as suitable for micro-imaging applications were utilized for high-resolution vessel wall imaging.

[0046] Simulation results were obtained for each axis and transverse sensitivity plots (B_Xy) were calculated using the root sum square of the transverse components (B_x and B_y). B_xy plots obtained from simulation results for each of the different wire configurations are shown in FIG. 108. Distortions from a true circular pattern in the radial homogeneity as well as inhomogeneities in the sensitivity of the opposed solenoid coil are observed when placing the connecting wires on opposite sides of the interventional device. Field distortions are also observable when the wires are placed inside the interventional device. When the wire was placed in the center of the device and along the same side of the device, only very slight distortions are observed. Examining the opposed solenoid without connecting wires, no field distortions were observed.

[0047] Bio-Savart simulation results based on the different coil design parameters show the following trends: (1) a direct relationship between coil radius and B_xy (FIG. 109A) (2) a direct relationship between number of windings and B_xy (FIG. 109B) (3) a parabolic relationship between pitch and B_xy with a maximum occurring at a spacing equal to about one diameter of the wire (FIG. 109C) (4) and a parabolic relationship between B_xy at a given radial imaging depth and coil separation with a maximum of about 2 to about 3 times the diameter of the coil (FIG. 109D).

[0048] When examining the test coils for simulation validation, QLOADED measurements obtained for each of the different coils was measured to be approximately 30 for each configuration. Values for each coil are shown in Table 101. Line plots taken from phantom imaging studies reveal the maximum signal intensity and signal roll-off characteristics for each of the coils examined (FIG. 107).

Table 101

[0049] In vivo and in situ porcine imaging experiments show vessel wall structures and surrounding vasculature being resolved (FIG. 110). In vivo images (10A-10B) were obtained with a temporal resolution of about 15 seconds per slice and an achieved in-plane resolution of about 240μm. In situ imaging results (10C-10D) show an in-plane resolution of about 160μm, and images were acquired with a temporal resolution of about 9 seconds per slice. In both imaging experiments, image resolution and SNR were sufficient to depict a thickening of the adventitia and the media of the vessel wall. [0050] Biot-Savart simulations model the performance of receiver coils used in intravascular imaging applications. The ease of which parameters can be adjusted make this method of investigation preferred over the construction of several coils to examine how individual parameters can affect coil SNR and homogeneity.

[0051] The results of Biot-Savart simulations for connecting wire placement show that placement of the wires either in the center of the interventional device or along the same edge of the device is preferable for maintaining field homogeneity. There is a potential with these two methods of a loss of a usable lumen in a standard catheter if the connecting wires are to be placed in the center, and the presence of small field inhomogeneities near the wires if they are placed externally along the same side of the interventional device. Even with this potential limitation, these two methods still provide superior field homogeneity over configurations wherein the wires are placed on opposite sides of the device or along the inside wall of the interventional device. The field inhomogeneities observed in these configurations may be due to the "loop antenna"-like configuration created when the wires are placed in such a fashion.

[0052] Two of the simulation's variables show monotonic increases in B_xy without theoretical bounds, but instead are limited by physical constraints. The dimensions of the vasculature to be investigated ultimately limit the diameter of the solenoid. While maximizing device size, it is also important to maintain sufficient space for avoiding issues such as vessel occlusions and/or loosening plaque components. The number of windings of the solenoid is also subject to physical limitations. These include the increase in resistance that comes from increased use of wire necessary to increase loops, which lowers the overall Q of the circuit. The inductance of the solenoid also increases as the number of windings increases. In order to create a resonant circuit, the inductor is matched with a capacitor to define a resonant frequency. A balance may be struck with the capacitor; extreme values generally are not to be utilized in a circuit. Therefore, about 5 to about 10 windings are used to maintain this balance. A parabolic relationship is observed when altering the pitch of the windings, with a maximum occurring at pitch spacing of one wire diameter. Less spacing results in field line cancellation between winding elements, while increasing the spacing beyond optimum results in the loss of flux linkage between the winding elements of the solenoid. With the spacing of the two coils in an opposed solenoid configuration, a maximum is seen when the coils are placed at a distance approximately two to three times that of the diameter of the solenoid coil. The distance between the opposed solenoid coils may be maximized to increase the length of vessel that can be imaged without necessitating the repositioning of the imaging coil. However, too great of a separation distance may result in signal losses due to lack of field-line coupling between the two opposed solenoid coils.

[0053] Variable capacitors, for circuit tuning and matching, were utilized in the construction of coils for simulation validation. They avoided inexact tuning and matching that would result from the use of surface mount capacitors with discrete capacitance values. Substantially accurate tuning and matching through variable capacitors was achieved at the cost of slightly lower overall Q values. However, similar Q values were obtained from measurements taken for each of the opposed solenoid coils. Uniform phantom imaging experiments reveal that the coil design provides high signal intensity and less signal roll-off as compared with other coil designs of the same diameter. When compared with the larger diameter coil, the performance of the optimized solenoid is nearly equivalent in a device that is 35% smaller than the non-optimized 14F opposed-solenoid imaging coil. Obtaining equivalent performance in a smaller-diameter coil is advantageous, as it would allow for the interrogation of smaller blood vessels. [0054] Porcine imaging experiments reveal the ability of the coil to resolve vessel wall structures. In vivo imaging reveals the vessel wall and surrounding vasculature with minimal motion artifact and good radial homogeneity. In situ imaging allows for even smaller FOV imaging at the expense of longer acquisition times. Structures of the vessel wall are depicted at a resolution of approximately 200 μm. [0055] While the aforedescribed coils are well suited for either tracking or imaging, a phased-array coil that advantageously combines both features will now be described. Referring now to FIG. 111 , two independent solenoid coils were wound in opposite direction on a cylindrical form. When considered together, the design and geometry may appear to be similar to that of a conventional opposed solenoid antenna. However, each coil may be individually tuned, for example to 63.6 MHz, matched to 50 Ω, and connected to a separate receive channel of the MR system. Each coil could therefore be turned on or off independently during the MR experiments, and signals from the single coils could be processed together or independently depending on whether the tracking or imaging mode was being chosen. Each coil may therefore have a unique and well localized sensitivity for tracking. The combined coils have an extended length of high radial homogeneity between the individual elements for high-resolution imaging.

[0056] A copper wire of about 30 AWG was used for both coils. Each solenoid probe had a 5F diameter, included 5 windings, and had a length of 4.5 mm. The gap between the coils was chosen to be about 1 cm based on Biot Savart simulations using the boundary conditions of a 5 F coil diameter. Tuning, matching, and passive decoupling with crossed diodes was performed on the tip of the catheter device in order to reduce electrical losses. The coil was encased in a biocompatible polymer shrink tubing over a length of about 40 cm. The maximum outer diameter of the prototype imaging/tracking catheter was about 12 F at the location of the tune, match and passive decoupling circuitry. [0057] In addition, a single channel opposed solenoid coil (i.e., not phased array) of identical geometric size and primary winding distributions was built for performance comparison (FIG. 111 b). Both the conventional opposed solenoid and the two elements of the phased array were connected using micro-coaxial cable.

[0058] All experiments were performed on a 1.5 T whole body scanner (Magnetom Sonata, Siemens Medical Solutions, Erlangen, Germany), equipped with eight RF receivers each with ±250 kHz maximum bandwidth, as described above. The MR system had 40 mT/m maximum amplitude gradients and a minimum gradient slew rate of 200mT/m/ms. RF excitation was performed using the standard body coil. The integrated panoramic array in combination with the real-time scanner interface allowed variable coil selection and combination on separate receiver channels including individual or combined image data processing, respectively. The system was configurable via the standard user interface. Imaging and tracking software implementing a subset of aspects of the systems and methods described herein was employed (e.g., on a Siemens Integrated Development Environment for Applications (IDEA) platform and Image Calculation Environment (ICE) platform) for pulse sequence design and image reconstruction. IDEA and ICE are built upon the C++ programming language, which provides a high level of software flexibility.

[0059] Measurements of the quality factor, Q, were performed for the single channel opposed solenoid coil and for each individual element of the micro-coil array, and compared to each other. Unloaded and loaded coil Q measurements were made. Unloaded Q's were measured by suspending the coil in air away from conductive material; loaded Q was measured by placing the coil in a saline filled phantom. Both unloaded and loaded Q was calculated from the polar plot of an S11 reflection type measurement on a network analyzer.

[0060] To compare the general signal characteristics of the coil array with the conventional single-channel opposed solenoid coil, multi-slice proton density weighted FLASH (TE 5ms, TR 20 ms, flip angle 40°, matrix 5122, FOV 120 mm, SL 3 mm) images were acquired with the catheter probes placed along B0 within a vessel phantom. The phantom consisted of two concentric NMR sample tubes, that were sealed at one end, and inserted into a 200 ml plastic bottle. The cylindrical cavities were filled with different solutions of saline and copper sulfate to establish different contrasts between the inner lumen,^' the simulated vessel wall and the simulated non-vascular tissue.

[0061] The data received simultaneously from both elements of the coil array were reconstructed individually to assess the spatial separation of signal peaks, and also used to calculate combined images using Roemer's sum of square algorithm. These combined images were compared to those acquired with the single-channel opposed solenoid coil using the same imaging and slice parameters. The spatial separation of signals from both elements was measured from slices strictly coronal to the device. Profiles were plotted parallel to the catheter's longitudinal axis and the location of the highest peak intensity was identified for each element. If the distance between these peaks was in the order of the geometrical distance between the elements, the separation was considered to be sufficient.

[0062] Tracking and imaging capabilities of the coil configuration were tested in phantoms, and in vivo in 5 domestic farm pigs (25-42 kg bodyweight). All in vivo endovascular procedures were performed under general anesthesia in accordance with protocols approved by the institutional animal care and use committee (IACUC) at our institution. The induction of anesthesia was achieved via intramuscular injection with 4-6 mg/kg of a mixture of Tiletamine HCL and Zolazepam HCL

(Telazol, 100mg/ml, Lederle Parenterals, Inc., Carolina, Puerto Rico). Anesthesia was maintained through continuous intravenous infusion of a mixture of 2 mg/kg Xylazine (Xyla-Ject, 20 mg/ml, Phoenix Pharmaceutical, Inc. St. Joseph, MO) and 20 mg/kg Ketamine hydrochloride (Ketaject, 100 mg/ml, Phoenix Pharmaceutical, Inc, St Joseph, MO). After local anesthesia, a 14F sheath was inserted into the proximal femoral artery under ultrasound guidance. The coil was then advanced through this sheath and placed in the abdominal aorta at different levels close to the origin of the renal arteries.

[0063] Real-time device tracking and automated imaging slice positioning was performed using, for example, ICE/IDEA tracking software previously described. The tracking software incorporates three main components: (a) a fast, active device localization module, (b) data processing software that calculates a position and direction vector that define the 3D position and orientation of the catheter; these values were input into (c) a FISP or True FISP real-time imaging module, that used the updated scan plane parameters to acquire and reconstruct new image data online, that depicted the catheter and gross anatomical features. Due to the flexibility of IDEA and ICE no additional hardware or software was needed to perform all tracking, localization and slice position updating, other than the catheter coils.

[0064] For active device localization, nonselective RF pulses were employed to excite all spins within the FOV. In subsequent scans, 1 D projection signals (FOV 400mm, matrix 1*256, TE 2.5ms, TR 5ms) from all 3 axes of the scanner were acquired. No coil other than the two-channel catheter coil array was used for signal reception during this first module (total acquisition time = 15ms). Due to the limited sensitivity of the solenoid coils mounted on the catheter, only spins close to the respective coil contributed to the signal acquired from the respective, individual channel (FIG. 112). Three projections (e.g. in the physical x-₍ y-, and z-directions of the magnet coordinate system or another possibly non-orthogonal reference coordinate system), and simultaneous signal reception via the two independent receiver coils/channels led to unambiguous determination of the position of each receiver coil in three dimensions. These positions were used to calculate the location and orientation of the subsequent imaging plane for device guidance.

[0065] A real-time FISP sequence with TE 3ms, TR 6ms, flip angle 15°, FOV 350x350 mm², matrix 128*128, slice thickness 5mm, and a bandwidth of 250Hz Pixel was used for device guidance. The acquisition time per image was 768ms, resulting in a frame rate of about 1.5 images per second. Spine- and body phased array coils were selected for signal reception. In addition, the elements of the catheter coil were also enabled during device guidance for a better visualization of the catheter's tip. They appeared as areas of high signal amplitude in the images.

[0066] In vivo tracking experiments in 5 domestic pigs were conducted using the previously described localization and guidance protocols. The two-channel micro- coil was inserted into the proximal femoral artery and advanced from the iliac artery to the abdominal aorta. The following experiments were conducted in a phantom and in the 5 pigs to test the tracking software performance with the newly proposed device, both in vitro and in vivo. First, the orientation and position of the tip of the catheter coil were determined via accurate manual localization using standard scout imaging. Then, 100 near-real-time image frames were collected with the device stationary and the system recalculating and refreshing the slice parameters for every frame from only the current localization information. Ideally, slice position and orientation should remain constant in its original state (i.e. centered on the tip of the stationary catheter coil) and should be identical to the values obtained manually with the scout imaging. Variations from the manually determined location and orientation were used to determine the standard deviation or accuracy of the tracking system. Finally, these tests were repeated with the catheter being moved within the water phantom and along the abdominal aorta, respectively. An imaging protocol option was chosen that was designed to position the tip of the catheter in the center of the updated frame. Deviations of the actually-depicted device position and the center of FOV were used as a measure of the system's accuracy during catheter advancement.

[0067] After successful guidance of the device to a remote site in the abdominal aorta, several high-resolution images of the vessel wall were acquired. Most imaging protocols used were based on the known steady state free precession technique TrueFISP. The following imaging parameters were selected:

(a) In vivo imaging: TR 13ms, TE 6.5ms, SL 2.5mm, 240 μm in plane resolution, FOV 23*30mm, matrix 96*128, TA 15sec/image.

(b) In situ imaging: TR 14.6ms, TE 7.3ms, SL 2.5mm, BW 130

Hz/Pixel, FOV 13*20mm, matrix 80*128, in plane resolution ~160μm, TA 9sec/image.

[0068] Unloaded Q was 40 for the proximal element and 34 for the distal element; loaded Q's of 39, and 34, were measured with the coil placed in a saline filled phantom. The conventional opposed solenoid coil had an unloaded Q of 45 and a loaded Q of 40.

[0069] FIG. 113 shows representative slices from FLASH experiments employed to explore the spatial sensitivity pattern of the phased array device in comparison to a conventional opposed solenoid coil. The images in FIGs. 113a and 113c show the spatial sensitivity of the individual coil elements along the longitudinal axis of the device - here, in a coronal plane through the center of the device. Corresponding profiles parallel to the long axis are depicted in FIGs. 113b and 113c, respectively. These profiles show a high-amplitude peak and a separation of the area of highest sensitivity for each coil element which is greater than the gap between the elements. This feature supports the advanced device localization method based on simultaneous projection data sets from both RF channels. Combined images, formed from the square root of the sum of squares of the images in FIG. 113a and 113c, are depicted in FIG. 113e and compared to the respective coronal slice from the conventional opposed solenoid coil (FIG. 113g). The sensitivity patterns are found to be similar for both coil designs. Both share the comparable spatial profiles, sensitivities, and longitudinal B1 inhomogeneity including the typical regions of zero sensitivity that are inherent to opposed solenoid designs. A comparison of the axial sensitivity of both devices is performed in FIGs. 113f and 113h. The sensitivity of the micro-coil array and the regular coil is highest adjacent to their outer wall, both have B1 sensitivities that drop off rapidly, and both have profiles that are nearly identical in radial symmetry. There is substantially no observable difference in axial and longitudinal sensitivity. Thus, the array coil is well suited for high-resolution imaging, as was the single channel coil.

[0070] Catheter advancement and retraction were successfully tracked in both vessel phantom and porcine imaging experiments. FIG. 114 shows representative slices from a real-time movie acquired in the abdominal aorta of a pig. Device tracking and automated slice positioning (location and orientation) was judged to be reliable, robust, and accurate in over 1000 tracking/imaging frames. The software automatically updates the scan plane and orientation and the elements can be individually switched on or off during image acquisition. The measured success rate was 100% for the motionless catheter; the error rate for the moving catheter in the aorta was less than 3%. In phantom and in vivo experiments, the inaccuracy was found to be less than 2mm of displacement error and of 2° orientation error.

[0071] FIG. 115 shows an in vivo example where the device location information was successfully used for automated slice positioning and then for vessel wall imaging. An in-plane resolution of 240 μm was achieved within 15 sec per slice in these images of the vessel wall. Substantially, no artifacts from arterial flow or device motion compromise the image quality, which confirms the array coil's capability for high-resolution endovascular imaging.

[0072] FIG. 116 depicts representative slices out of a dataset acquired in situ using the phased array micro-coil. These images have an in plane resolution of 160 μm and were acquired in 9 sec per slice. Resolution and signal-to-noise ratio were sufficient to depict in the deceased animal (a) the collapsed arterial vessel wall (compare to the substantially round shape of the aorta in vivo (FIG. 115)), and (b) a pronounced thickening of two layers of the vessel wall, the adventitia (outermost black ring structure (arrow)) and media (adjacent grey rim (arrowhead)).

Furthermore, a non-uniform signal-to-noise level was observed in the individual images of this series due to longitudinal B1 inhomogeneity of opposed solenoid devices.

[0073] The disclosed dual-purpose device for simultaneous MR tracking through the vascular system and high-resolution imaging of the vessel wall incorporates concepts from array coil technology, since two solenoid coils, that are wound counter to each other, are mounted on a conventional catheter, and connected to individual receiver channels of the MR system. In combination with 3 gradient projections, each coil provides the unique opportunity for fast device localization. The complete localization process is accomplished within about 15ms.

[0074] Device position and orientation can be determined easily and unambiguously if tracking is done with both coils. Peak ambiguity is a known problem and reported in earlier work on active device tracking, especially when single-channel devices are used. Previously, more than 6 projections were acquired and additional post-processing (e.g. cluster analysis) was performed to correctly identify peaks. This is more than double the acquisition time compared to the localization method disclosed herein. These additional acquisitions and computations may be avoided by using a two-channel approach. Furthermore, the tracking accuracy will not be compromised when using the dual-channel device, as demonstrated when comparing with accuracy measurements results obtained from of similar tracking software and a dedicated 2-marker single channel tracking antenna.

[0075] Phantom and in vivo experiments as well as direct comparison with a conventional opposed solenoid coil of equal geometry demonstrate that the proposed phased array design provides opposed solenoidal imaging performance if the signals from the two independent receiver channels are combined appropriately. [0076] According to another practice, a coil-less and/or wireless tracking device can be advantageously employed for tracking and visualization of a catheter in interventional MRI. Passive tracking methods using susceptibility artifacts can provide adequate catheter visualization. However, these methods are dependent on the orientation of the device (coil) in the magnetic field and on the slice thickness of the acquired image acquisition. Active tracking methods utilizing tuned micro-coils provide accurate localization for scan plane determination, but are not as useful for steering catheters through the complex vasculature because of their point-like nature. Guide wire antennas provide a larger field-of-view for catheter visualization/profiling, but suffer from local heating which may result in significant tissue damage.

[0077] The disclosed device is based on the selective excitation of a chemically- shifted NMR signal source within a catheter. This tracking/profiling method can provide the necessary selectivity and large field of view for catheter visualization necessary to allow catheter steering when overlaid onto a previously acquired roadmap image. This new device provides the same capabilities as guide wire antennas without the risks of localized tissue heating. [0078] A prototype catheter was created by infusing 1 ml of concentrated acetic acid (σ ~ 7ppm) doped with 1 mM Gd contrast (Magnevist™, Schering AG) into a plastic 1 ml syringe (ID = 4.7mm). This solution provides a signal source with a proton chemical shift frequency distinct from typical tissue protons. The syringe was placed into a vascular phantom and near a volunteer's head to develop the tracking/profiling sequence. A FLASH (Fast Low Angle SHot) sequence was developed with a 10ms chemical shift selective excitation (CHESS) pulse to excite the off-resonance spins (TR/TE/FA = 20ms/10ms/30°). No slice-select gradients were applied in this sequence. Gradient shimming was applied prior to the image acquisition to limit the effects of field inhomogeneities. Imaging acquisitions toggled between CS-FLASH and conventional FLASH to provide alternate catheter and anatomic images. The images from the chemical shift-selective FLASH sequence (CS-FLASH) were compared with standard slice-selective FLASH (TH=5mm) to demonstrate the capability for acquiring catheter-only images to be overlaid onto previously acquired anatomic images needed for catheter steering. [0079] Axial and sagittal views of the vascular phantom and volunteer head images are shown in FIGs. 117 and 118, respectively. FIGs. 117a,b and 118a,b are FLASH images with a slice-selective excitation pulse. The acetic acid syringe is identifiable within the vascular portion of the phantom and near the volunteer's left ear in the axial image (FIG. 118a). FIGS. 117c,d and 118c,d are from the CS-FLASH sequence with the CHESS pulse resulting in the "catheter-selective" images. In an in-vivo intravascular procedure the sagittal catheter-selective image (FIG. 118d) would be overlaid onto the corresponding anatomic image (FIG. 118b).

[0080] A contrast ratio (or selectivity) of approximately 2:1 was observed for the aforedescribed chemically-shifted probe using acetic acid. This contrast ratio may be insufficient to track the probe in vivo, in particular when using active tracking. The contrast ratio is hereby a measure of the signal from the probe relative to the signal from the surrounding tissue. A low-contrast ratio may lengthen the data acquisition time which would be unacceptable for vascular imaging. Interfering noise from other sources may add to the difficulty of extracting a useful signal for tracking purposes. In addition, probes used in vivo have to meet certain biocompatibility standards.

[0081] Adverse health effects can result from exposure to the materials from which a device is made. The biocompatibility of a device depends on several factors, especially the type of patient tissue that will be exposed to device materials and the duration of the exposure. The tracking device should be biocompatible at least in the following aspects: (1) no introduction of sublethal or lethal effects as observed at the cellular level (Cytotoxicity); (2) no localized reaction of tissue to leachable substances (Intracutaneous Reactivity); (3) no adverse effect occurring within a short time after administration of a single dose of a substance (Acute Systemic Toxicity); and (4) no undesirable changes in the blood caused directly by a medical device or by chemicals leaching from a device (Hemocompatibility). Undesirable effects of device materials on the blood may include hemolysis, thrombus formation, alterations in coagulation parameters, and immunological changes.

[0082] Suitable biocompatible materials for vascular tracking and imaging are, for example, fluorinated compounds, such as fluorinated ethylene, polyether urethanes, and more particularly propylene polytetrafluoroethylene (PTFE). These compounds can provide an enhanced signal over the signal derived from the chemically-shifted probe. The contrast ratio can be further enhanced by labeling the compounds with stable isotopes. For example, compounds labeled with stable isotopes such as Glucose-1-¹³C and Glutamic-¹³C acid are used in Magnetic Resonance Imaging (MRI) techniques to render visible metabolic changes. MRI tracking probes made of, for example, fluorinated hydrocarbon compounds can be labeled with ¹³C, whereas biocompatible azo-compounds can be labeled with ¹⁵N. In this way, the MRI system can discriminate between the MR signal from the carbon atoms in the catheter itself and the MR tracking signal from the labeled PTFE of the probe. A catheter probe could then be fabricated entirely of a solid material that can be attached to or integrated with the catheter itself. When using fluorocarbons or boric acid instead of acetic acid, a contrast ratio or selectivity of greater than 5:1 , and even 10:1 or 20:1 can be achieved.

[0083] In another embodiment, the acetic acid in the liquid-filled lumen described above can be replaced with boric acid. Boric acid produces a signal with a greater bandwidth, for example, 3 MHz at 1.5 Tesla, and can in addition be tagged with ¹¹B. The boric-acid-filled lumen can also contain a contrast material. The improved contrast ratio or selectivity are in particular beneficial for spiral imaging and other MRI signal acquisition techniques where k-space data have to be acquired within 100 msec or less.

[0084] Two FLASH sequences are used to acquire both the anatomical and catheter-selective images. However, other combinations of sequences and materials can be used to optimize the anatomic and catheter imaging independently for the particular interventional procedure (i.e., optimized for speed, catheter/background contrast, tissue/vessel contrast, resolution variations between images, etc.). The bandwidth, center frequency, and the magnitude and phase of the CHESS pulse can also be modified to generate images with selected suppression bands allowing for better visualization of the catheter within the vasculature. The method is easy to implement on conventional scanners and requires only a single receiver channel. The sequence and catheter design can be further optimized, including construction with fully biocompatible materials and a reduction in catheter size, facilitating real-time in-vivo active tracking for intravascular catheters. [0020] The invention is directed to a system and method for real-time catheter tracking and adaptive imaging using MRI. In particular, the system and method described herein can be used to track and position catheters and stents in a human body.

[0021] Referring now to the figures, and in particular to FIG. 201 , an adaptive parameter system software architecture is shown that was interfaced with a 1.5 Tesla Siemens Sonata scanner (Siemens Medical Solutions, Erlangen Germany) using the Siemens Integrated Development Environment for Applications (IDEA) and Image Calculation Environment (ICE) for pulse sequence design and image reconstruction. IDEA and ICE are built around the C++ programming language to afford the developer a sufficient degree of software flexibility. The software architecture has three main components: (1 ) a fast device localization method, (2) data processing software that performs velocity calculations and updates image parameter values, and (3) a real-time imaging technique that automatically incorporates the updated scan plane parameters, performs appropriate pulse sequence revision, acquires new image data, and reconstructs the image online. These components form a closed feedback loop system that, when continuously repeated, provides an interventional environment with real-time imaging and image parameters that adapt to the changing clinical circumstances.

[0022] As seen in FIG. 201A, the system alternates between acquiring (substantially in real-time) a two-dimensional (2D) catheter-selective image and an image of the surrounding tissue (e.g., a tissue map image). The system then localizes the catheter-based markers to localize the device, and updates the scan plane position based on the location, trajectory, and/or orientation of the tracking coil. The image parameter values are adjusted based on the velocity of the tracking coil calculated from two or more successive position measurements, using Eqs. [201] or [202] described below. The image resolution and/or the field of view (FOV) are adjusted using the adjusted image values and a new image is acquired reflecting the new image parameters.

[0023] FIG. 201 B depicts an alternative embodiment of the systems and methods described herein. The process pictorially depicted by FIG. 201 B includes acquiring a 3D roadmap based on maximum intensity projections (MIP), such as, without limitation, an angiogram roadmap, in addition to the 2D real-time images referred to in FIG. 201A.

[0024] Various device localization techniques can be used, with one exemplary technique incorporating a single coil tip tracking method depicted in FIG. 202, which can provide information about the catheter's three-dimensional position within the magnet. FIG. 202 depicts a tuned resonant circuit that is capacitively coupled to the MR system and is mounted on the tip of a catheter to provide information about the three-dimensional position of the device.

[0025] The second tracking method, depicted in FIG. 203, includes an analytic radial tracking method that uses two active regions or markers on a catheter and can provide information about the three-dimensional position and orientation of the catheter. In particular, FIG. 203 depicts a tuned resonant circuit with two active regions. Additional detail of the coil designs will be described later. [0026] Both localization methods require the collection of a limited number of 1 D projections (the first tracking method requires 3 projections for catheter position, whereas the second tracking method requires 8 projections for catheter position and orientation). These localization projections are collected prior to the acquisition of each set of image data. The raw k-space projection data is sent to the image reconstruction computer where the data processing software (written in C++ and developed in ICE) determines the location of the tracking markers. The software then computes six values: three positional values that define the 3D position of the catheter and three angles that define the 3D orientation of the catheter. When using the first localization method, the three angles are set to fixed values that define a standard transverse, sagittal or coronal plane. Values are defined within the scanner's X, Y, and Z coordinate space and the positions have units of pixels.

[0027] Using standard methods provided in ICE, a dedicated real-time link is established between the image reconstruction computer and the hardware control computer, which executes the pulse sequence software on the scanner. The six position and orientation values are sent via this real-time link, where the pulse sequence software (written in C++ and developed in IDEA) accepts and stores them for use. The position and orientation values are then converted from magnet XYZ coordinate space, to a coordinate system defined by the patient using Read, Phase, and Slice-Shift axis. The 3D positions are also converted from units of pixels, to millimeters. A real-time software kernel ensures that these computations are performed within a 20 ms pause located after the collection of the localization projections and before acquisition of the next set of rapid image data. The new catheter position and orientation information is then used to automatically define the new scan plane just prior to image data acquisition.

[0028] The pulse sequence software uses localization data from multiple time points to calculate the speed of the device. A variable-point finite difference digital filter is used for this calculation. The number of time points used in the digital filter is adjustable via the user interface, allowing the clinician to control the system's sensitivity to sudden changes in catheter speed. The catheter speed is then used to adjust the value of selected image acquisition parameters (e.g., image resolution, temporal resolution, bandwidth, field of view, slice thickness). A variable image parameter, P(V), is expressed as a function of the device speed, V, with limits for the device speed set, for example, before each procedure. The value of this function determines how each adaptive parameter is set, relative to its full range of acceptable values. P(V) need not be calculated more than once before all of the selected image parameters can be updated. [0029] Referring now to FIGs. 204 and 205, the system has two types of functions describing the relationship of the catheter speed to the variable image parameter: a step function (e.g., a binary step function is shown in FIG. 204, but a multi-step function is also allowable), and a continuous sigmoidal function (FIG. 205). The binary set function uses a velocity threshold, which is adjustable via the user interface, to determine if the device is moving or stationary (Eq. [201]). When using this function, the adaptive image parameters will be set to one of two values:

-P_min ϊf V < V_Thresh

P(V) =

P max Otherwise

[0030] If the current speed of the catheter is found to be less than the designated user-defined threshold Vr_hres_h, then the selected image parameters will be set to one predetermined value P_m\n (which may be ideal for imaging with a stationary or slowly- moving catheter). Similarly, if the calculated catheter speed is larger than the designated threshold Vr_hresh, then the selected image parameters will be set to a different predetermined value P_max (which may be better suited for a faster moving catheter). [0031] A binary function may be inadequate to adjust to a changing insertion speed. Accordingly, a continuous mode can be employed that uses a smoothly- varying function of catheter speed, such as a sigmoidal function, to adjust the adaptive image parameters (Eq. [202]):

[0032] This function is depicted in FIG. 205. The upper and lower asymptotes prevent image acquisition parameters from being set in ranges that would cause the pulse sequence to exceed hardware limitations. This function's center velocity, V₀, and static sensitivity, S, are also adjustable via the user interface. Any image parameter can be automatically and separately adjusted by this system, and the system can incorporate virtually any MR imaging technique.

[0033] A user-defined range for each specific image parameter is set prior to the experiment. Subsequently, during the intervention, the binary and/or continuous set functions are used to set the values of the parameters relative to their specified respective ranges. The clinician is able to designate which image parameters to treat as adaptive parameters via the user interface. The parameter set function calculations and the image parameter update are performed by the real-time kernel immediately prior to the acquisition image data. [0034] In the exemplary system and method, a variety of pulse sequences have been integrated into the adaptive tracking software. One pulse sequence is a TrueFISP pulse sequence (TR= 5ms, TE= 2.5ms, FA= 70°, matrix= [128x128 to 1024x1024], FOV= [150mm to 300mm]). Another pulse sequence is a FISP sequence (TR= 5ms, TE= 2.5ms, FA= 70°, matrix= [128x128 to 1024x1024], FOV= [150mm to 300mm]). Another pulse sequence is a FLASH sequence (TR= 8ms, TE= 4ms, FA= 15°, matrix= [128x128 to 1024x1024], FOV= [150mm to 300mm]). Yet another pulse sequence is a Radial True-FISP sequence (TR = 3.6 ms, TE = 1.8 ms, FA = 70 degrees, Matrix = 128 x 128, FOV = [50 mm to 400 mm], Radial Lines = [64 to 512]). [0035] FIG. 206 shows an exemplary user interface for the system with an online display that updates images in real-time as soon as new image data is reconstructed, or shortly thereafter. The interface also allows the clinician to interactively toggle and configure the device tracking and adaptive parameters; this can be done, for example, using a dialogue card. [0036] To evaluate the system's feasibility and performance in vitro and in vivo, trials were conducted in two vessel phantoms and eight porcine imaging experiments, using methods approved by our Institutional Animal Care and Use Committee. To quantify the system's reliability and precision, 100 image frames were collected (in vivo) in which the system updated the slice position and orientation so that it was centered on the tip of a stationary catheter. This procedure was then repeated on a catheter that was moving at clinically relevant speeds within a vessel phantom. The experiments measured image data with two active receive channels: the imaging coil and the tracking coil. The data from these two channels was combined before displaying the image so that the tracking markers appeared in the image as areas of high signal amplitude. The distance between the position of the tracking markers within the image and the center of the image was used as a measure of the system's accuracy.

[0037] FIG. 207 shows results from in vivo porcine experiments in which the catheter was inserted throughout the length of the abdominal aorta. FIG. 208 shows image data collected in a vessel phantom experiment. In both the phantom and porcine trials, the resolution and FOV were automatically varied. Both sets of images illustrate a temporal sequence in which the catheter is slowed to a stop.

[0038] The system was able to accurately localize a motionless catheter 100% of the time and a moving catheter 98% of the time. The small error rate is substantially due to the system misidentifying flowing spins outside the imaging plane as the signal from the tracking markers. When using the single coil tracking method, the system collected all of the necessary tracking data within 15ms, whereas the figure for the two-marker method was about 25ms. An additional 20ms was then required to perform the localization, velocity calculations, and updating of the image parameter values. The system ran continuously and responded in real-time to calculated changes in all eight in-vivo and both phantom trials. Following parameter determination, the system successfully responded to changes in device speed by dynamically adjusting specified image parameters. In all cases, the slice plane location and orientation was automatically placed at the catheter tip using information extracted during the localization phase.

[0039] The adaptive interface was further tested by performing MR imaging- guided renal artery stent placement, which was performed in two pigs using a catheter-based system that interactively adjusts the scan plane and automatically adjusts a number of imaging parameters in a manner described above.

[0040] Although renal MR angiograms have become state-of-the-art for the diagnosis of renovascular disease, conventional fluoroscopy remains the guidance method of choice for renovascular intervention. However, MRI-guided stent placement has recently become recognized as a feasible alternative. [0041] Miniature radiofrequency coils were implemented into a catheter to determine the exact 3D position of the device in real time. As mentioned above, multiple coils can be implemented to derive information about the device orientation. [0042] In the experiment, multiple MR image-guided procedures were performed in two fully anesthetized pigs using a 1.5 T scanner (Magnetom Sonata, Siemens Medical Solutions, Erlangen, Germany). Gadomer-17 (Schering AG, Berlin , Germany), an intravascular contrast agent for prolonged vascular enhancement, was used to acquire a gradient echo baseline 3D-MRA and to improve artery opacification during the intervention. A 5F C1 -catheter equipped with two short single loop coils at the tip was then used to catheterize the renal arteries using a transfemoral approach.

[0043] The three-dimensional position and orientation of the catheter micro-coils were determined every 300ms by means of an analytic radial tracking method using the two active catheter-based coils. The position data was used to define the MR scan plane position and orientation of a steady-state free-precession (SSFP) sequence acquired with a frame rate of three images per second. The field of view (FOV) of the images was adjusted according to the speed of the catheter movement: an FOV was enabled for a slowly-moving catheter and a larger FOV was used at higher catheter speeds. Using a balloon catheter equipped with one short single loop coil proximal to the balloon (e.g., immediately proximal thereto), MR guidance was then used to place a stent over a wire into the ostium of the renal arteries of both pigs.

[0044] The procedure time was measured, and the stent position was verified using conventional angiography. The high-amplitude signal from the coil of the instrumented balloon catheter was used to exactly position the stent at the level of the renal artery ostium in both pigs. The stent deviation as measured with conventional angiography was less than 3mm. The procedure times were 14 and 19 minutes, including the MRI acquisition.

[0045] The systems and methods described herein demonstrate a method for intravascular imaging that significantly improves the performance of image-guided intravascular procedures. With this, the physician has means with which to control the scanner and dynamically alter specific image parameters during an MR procedure. Allowing the MR scanner to respond to a moving catheter by adjusting the value of imaging parameters creates a more natural interface with the MR scanner for the clinician during intravascular procedures by eliminating the need for manual adjustment of the scan plane position or specific image parameters during the intervention. These scan plane adjustments are applied to imaging slices that are also following the location and orientation of the catheter. Hence, advancing the catheter more slowly will automatically improve the resolution or SNR properties of the images, or can even effect a total change in tissue contrast to allow more accurate characterization of vessel wall pathology if the clinician wishes to see more detail in a certain region of interest.

[0046] The interventional systems and methods described herein can be incorporated into virtually any imaging protocol. The tuned coils operating as capacitively-coupled localization markers can be manufactured to be small enough to fit a variety of clinical catheters. Other fiducial marker designs may also be incorporated with ease into the systems and methods disclosed herein. For example, and without limitation, multiple inductor types can be used (e.g. single-loop, solenoid, and opposed solenoid inductors); capacitively-coupled tuned resonant circuit markers may be connected to a single receiver channel or multiple channels as phased-arrays; these phased-array markers may also be used for catheter-based vessel wall imaging). According to one practice, inductively-coupled tuned resonant circuit markers may be used. According to another practice, markers filled with a distinct signal source, such as a fluid with a large chemical shift relative to water, can also be used. An underlying theme is that a localized signal is available that may be readily segmented from a non-marker background signal.

[0047] The systems and methods according to the invention have been able to reliably (with accuracy of at least 98%) localize a catheter to within 2mm and 1 degree of rotational error; these values are comparable to existing commercially- available MR-tracking technology which is unable to provide real-time tracking with adaptive imaging. The observed small localization error rate of at most about 2% occurs at least in part because the imaging plane contains the tracking markers and tends to saturate the spins surrounding the markers, making it possible for signals outside the imaging plane (and farther away from the tracking markers) to be mistakenly identified as the tracking markers. [0048] The problem can be addressed by allowing more time for signal recovery, introducing a non-selective saturation between imaging and localization, or constraining the extent of rotational or positional shifts between sequential images. During normal use, the system minimizes localization errors by monitoring the distance that the catheter has moved between each image frame. The slice position and adaptive parameters are not updated if the detected change in catheter position distance is greater than a pre-determined value corresponding to the preset maximum clinically permissible insertion speed. In practice, this ensures that a catheter localization failure during a given frame will generally not cause the slice to be placed at an incorrect position and the adaptive image parameters to be set to incorrect values.

[0049] The systems and methods described herein use standard clinical hardware, with the exception of the small markers that were affixed to the catheter. The systems and methods also use standard clinical gradient and RF hardware, control computer hardware, and reconstruction computer.

[0050] The software interface was also merged into the standard clinical interface provided by the vendor. These features make the systems and methods according to the invention easy and inexpensive to implement and intuitive for most who are experienced with the operation of the MR imager.

[0051] The adaptive tracking system requires, at most, an additional 60ms per image. This allows real-time imaging sequences to continue to operate in real time, while providing a great deal of increased functionality and flexibility. The ability to respond, in real time, to changes in device velocity allows the scanner to automatically adjust a number of image parameters. Without requiring the clinician to intervene, the scanner can automatically increase resolution and decrease frame rate as the catheter slows, or increase the field of view (FOV) as the catheter's speed increases.

[0085] Moreover, according to one practice an adaptive tracking system uses markers with a resonant frequency that is distinct from a resonant frequency of surrounding tissue; in this practice, the MR scanner employed has RF hardware configured to transmit and/or receive signals at this distinct frequency and a marker has an internal signal source characterized by a resonant frequency substantially equal to the distinct frequency at the scanner's field strength.

[0035] In one aspect, the methods described herein are directed, inter alia, to efficient reconstruction of high-quality MR images. In particular, the methods described herein can be applied to non-rectilinearly sampled data and spiral MRI sampling schemes.

[0036] In one exemplary practice, an approach for optimal reconstruction using rescaled matrices from non-uniformly sampled k-space data is described. Non- rectilinear data acquisition methods have advantages over rectilinear data sampling schemes and hence are often performed in magnetic resonance imaging (MRI). For example, projection reconstruction, i.e., radial trajectories, shows reduced motion artifacts, and spiral trajectories are insensitive to flow artifacts. Image reconstruction from non-rectilinearly sampled data is not simple, because 2D-lnverse Fourier Transform (IFT) cannot be directly performed on the acquired k-space data set. K- space gridding is commonly used as it is an efficient reconstruction method. Gridding is the procedure by which non-rectilinearly sampled k-space data are interpolated onto a rectilinear grid. The use of k-space gridding allows the reconstruction of images in general non-uniform sampling schemes, and thus gives flexibility to designing various types of k-space trajectories. [0037] Conventional gridding algorithms have been proposed that are robust to noise and do not require a significant computational burden; however, a profile distortion of the reconstructed image often appears, unless the density compensation function (DCF) is sufficiently optimized. Other proposed algorithms with improved DCF's are often complicated, and it is still difficult to calculate the Optimal DCF' in general non-uniform sampling schemes.

[0038] Another type of conventional gridding algorithm, 'Block Uniform Resampling (BURS)', requires neither a pre- nor a post-compensation step, and the reconstructed image is usually of high quality. Although the originally proposed BURS algorithm is sensitive to noise, it has demonstrated that SVD regularization techniques can avoid amplification of these data imperfections. However, it is often difficult to determine the regularization parameters in advance, as the k-space data signal-to-noise ratio (SNR) is usually unknown before reconstruction. [0039] Another image reconstruction algorithm, known as 'Next Neighbor re- Gridding' (NNG), obviates the complicated procedures of the gridding algorithms described above. The NNG algorithm consists of the following four steps: 1) Density compensation of k-space data; 2) Distribution of the k-space data into a large rescaled matrix by a factor of s (= 2m), where m is a small positive integer. The location of each datum in the large rescaled matrix is determined by rounding off the original k-space coordinate in the target rectilinear grid after multiplying it by the same factor s. If more than one datum share the same matrix coordinate, the mean value is stored; 3) IFT of the large matrix; and 4) Extraction of the original-sized matrix at the center. In brief, each acquired k-space datum is simply shifted to the closest grid point of a finer rectilinear grid than the original grid, in order to directly perform IFT on a non-uniformly sampled k-space in the NNG algorithm. The errors caused by the data shifts are usually quite small in the reconstructed image if the scaling factor s is sufficiently large. Specifically, s = 4 or 8 is sufficient in practice. [0040] As is the case with the conventional gridding algorithm, the image quality of the NNG algorithm depends on the DCF used in step 1. In other words, non- negligible profile distortions of the reconstructed image are often observed if the DCF is not well optimized.

[0041] In one embodiment, the systems and methods described herein are directed at a new image reconstruction algorithm from non-rectilinearly sampled k-space data. The newly proposed algorithm is an extension of the NNG algorithm described above and will be referred to hereinafter as the 'Iterative Next-Neighbor re-Gridding (INNG) algorithm' as it includes an iterative approach. Although the algorithm requires a number of Fast Fourier Transforms (FFTs) of re-scaled matrices larger than the original-sized rectilinear grid matrix, no pre-calculated DCFs are required in the INNG algorithm, and the reconstructed image is of high quality. When the size of the rescaled matrices is significantly large, it is often impractical to perform FFTs on them. To overcome this, a 'Block INNG (BINNG) algorithm' has been developed. In the BINNG algorithm, k-space is partitioned into several blocks and the INNG algorithm is applied to each block. It will be shown that if data imperfections are non- negligible, e.g., low data SNR and/or a small scaling factor, the background noise level in the reconstructed image is increased as the iteration progresses in the INNG/BINNG algorithms. However, the rate of the increase is usually not significant unless the data imperfections are substantial. Hence, an adequate choice of stopping criteria can reconstruct a high-quality image given non-uniformly sampled k- space data. The INNG/BINNG algorithms are a simple new approach to accurate image reconstruction and an alternative to the previously-proposed optimized gridding algorithms that does not require DCFs or SVD regularization parameter adjustments.

Basic INNG Method

[0042] Referring now to FIG. 301 , the basic procedures of the INNG algorithm are presented as a flow chart. Suppose that the originally-designed rectilinear grid size is N x N. The initial image of the INNG algorithm can be obtained by steps 2) and 3) in the Next-Neighbor re-Gridding (NNG) algorithm described above. In other words, a 2D-IFT is performed on a large rescaled matrix after k-space data are distributed without density compensation. These steps are shown by the process (a) to (b) in FIG. 301. After the initial image is reconstructed, all the matrix elements except the central N x N region are replaced by zeros. In FIG. 301 , the process (b) to (c) represents this step. This procedure is equivalent to multiplication of the matrix (b) with a 2D-rect window function of amplitude 1 in the central N x N matrix and 0 elsewhere in the image. Therefore, if the matrix (c) is 2D-Fourier transformed, the obtained matrix (d) is the result of convolution of the matrix (a) with a 2D sine function (which is 2D-FT of the 2D-rect function used in the previous process). After the matrix (d) is obtained, at the matrix coordinates where the original data exist in the rescaled matrix (a), the data are replaced by the original data values, as shown in the process (d) to (e) in FIG. 301. Other matrix elements are left unchanged in this process. Then, 2D-IFT is performed on the matrix (e) leading to the updated reconstructed image (b). The procedures (b) -» (c) - (d) -» (e) → (b) (surrounded by dashed lines in FIG. 301) are repeated until the difference between the updated reconstructed image (b) and the image (b) at the previous iteration becomes sufficiently small.

[0043] In the basic INNG algorithm, it is assumed that the Nyquist criterion is satisfied for the entire k-space region which spans from -kmax to +km_ax along both k_x and k_y directions. In other words, at least one datum must exist in any s x s matrix region in the sNx sN rescaled matrix. In a practical implementation of the basic INNG algorithm, if there are non-sampled regions in the k-space, the corresponding regions in the rescaled matrix are set to zeros. For example, in spiral trajectories, k- space regions outside of the circle with a radius |km_ax| are usually not sampled. Correspondingly, the regions outside of the circle with a radius sN/2 are set to zeros in the sNx sN rescaled matrix, when the original data are inserted at each iteration. This procedure is also performed in the facilitated INNG algorithm and in the BINNG algorithm introduced in the following subsections.

[0044] The INNG algorithm described above can be classified as a well-known optimization method 'Projections Onto Convex Sets (POCS)'. In MRI, the POCS method has been used in half-Fourier reconstruction, motion correction and parallel imaging reconstruction. In the POCS method, each constraint can be formulated as a 'convex set', which is known in the art. In the INNG algorithm, two constraints are imposed on the data (or the image) at each iteration, that is, (i) the finite-support constraint and (ii) the data-consistency constraint. The constraints (i) and (ii) correspond to the process (b) to (c) and the process (d) to (e), respectively in FIG. 301. The constraints (i) and (ii) can be expressed as the following two convex sets Ωi and Ω₂₎ respectively:

Ω₁ = {I(x) \ I(x_out) = 0} , [301]

where l(x) is the image matrix of a large FOV (sN x sN) and x_out represents all the matrix elements except the central N x N matrix.

Ω₂ = {I(x) | I(x) = F{D(n)} , D(n_orig ) = D_orig (n_orig )} , [302]

where F is the Fourier Transform operator, D(n) is the Fourier data matrix (sN x sN) of l(x), n_orig represents all the elements in the larger scaled matrix where the original data exist, and D_orig are the original data values at these coordinates.

[0045] The constraint (i) is based on the signal sampling theory in which all the sampled signals must be expressed as the summation of rectilinearly located sine functions. If all the data values in the large rescaled matrix can be expressed as the summation of the 2D sine functions (each of which is the FT of the 2D-rect function with amplitudes 1 in the central N x N matrix and zero elsewhere), all the image matrix elements except the central N x N region must be zeros. The need for the constraint (ii) is to keep the original data values at the original data locations for each iteration.

[0046] Suppose that the operators which project an image matrix l(x) onto Ωi and Ω₂ are Pi and P₂, respectively. The image reconstructed using the INNG algorithm can then be expressed as:

ΛH.ι(*) = ^jyϊ #. (*)} ^■ [303]

where the subscript of l(x) denotes the iteration number.

[0047] It can be shown that P1 and P2 satisfy the following relations:

∑W_m (*)>- W-i W}]² ≤ ∑tf.(* .»ι(*)]^a [304]

where ∑ denotes summation of all the elements in the sN x sN image matrix. P

and P₂ are called non-expansive operators. The composite operator P₂Pι is also non-expansive, that is,

∑[ ₊₁( ) - ₊₂( )]²

-P2P i_mA^χ)}f ≤ ∑[ W - „,₊₁ω]² • [306]

X X X

Note that Eqs. [304,305,306] hold whether or not the data are ideal. The algorithms with non-expansive operators have certain convergence properties. If the data distributed in the larger rescaled matrix are ideal, then the above iterative algorithm has a unique convergence point. However, if the errors contained in the data are non-negligible, a unique convergence point may not exist. Since both Pi and P₂ are linear operators, P2P₁ is also a linear operator. Thus, the reconstructed image at the A77-th iteration can be expressed as the summation of the image values that originate from the ideal signal components, i.e. the signal components which satisfy the condition (i) I_ideal>m(χ) and the image values that originate from the residual imperfect signal components n_m (x) :

/„, (*) = I_idmhm (x) + n_m (x) (= P₂P_X {I_M^ (x)} + P₂P_λ {n_m__t (x)}) [307] As the iteration continues, the first term in Eq. [307] leads to the ideal reconstructed image. The second term in Eq. [307] usually manifests itself as background noise in the reconstructed image, and the noise level is increased with iterations as will be seen in the ensuing section. However, the increased rate of the noise level is reduced as the iteration progresses. This can be understood as Eq.[306] also holds for n_m(x) . Therefore, it is expected that if the data SNR is within a practical range, and the scaling factor s is sufficiently large, the magnitude of I _aUm(^χ) '^s sti" predominant over that of n_m(x) after a certain number of iterations.

[0048] Exemplary stopping criteria used with this algorithm are described below, although other stopping criteria can also be used. As Eq.[306] indicates, the quantity U _m(^χ) ~ I _m+\(^χ) monotonically decreases with iteration number m. In the present

embodiment, the sum of the squared difference [I_m(x) -I_m+l(x)]² is calculated within the central N x N image matrix instead of the entire sN x sN image matrix to facilitate the computation. Hence, the following quantity d is measured to determine where to stop the iteration:

NxN where ∑ denotes summation of the central N x N image matrix elements. The

iteration is stopped if d becomes lower than a predetermined value. The predetermined stopping criterion will be denoted by d_s where the scaling factor is s in the following sections.

[0049] The basic INNG algorithm described above requires a number of FFTs. Furthermore, the amount of computation for each FFT is usually demanding when a rescaled matrix is large. In the following section, a 'facilitated INNG' algorithm is described which reduces the computational load of the described INNG algorithm above.

Facilitated INNG Method

[0050] Referring now to FIG. 302, the facilitated INNG algorithm is shown in a flow chart. As seen from FIG. 302, the facilitated INNG algorithm modifies the basic INNG algorithm by employing consecutively increasing scaling factors. The first step of the facilitated INNG algorithm is the basic INNG algorithm with a small scaling factor, e.g., s = 2. The image reconstructed from the basic INNG algorithm with s = 2 ((a) in FIG. 302) is usually affected by noise because the errors caused by data shifts are significant. However, the image (a) is roughly close to the image reconstructed using the basic INNG algorithm with a larger scaling factor. Thus, the image (a) will be used for the basic INNG algorithm with s = 4 in the next step. A 4/V x 4/V zero matrix with the central N x N matrix replaced by the image (a) ((b) in FIG. 302) is a starting image matrix for the basic INNG algorithm with s = 4. In general, this basic INNG algorithm with s = 4 requires significantly fewer iterations than the same algorithm that starts with FIG. 301 (a) to satisfy the same stopping criterion. To further reduce the errors caused by data shifts, the basic INNG algorithm with s = 8 will be performed next, in a similar manner. In other words, the image (c) is transferred to the center of an 8/V x 8/V zero matrix to form the matrix (d) as shown in FIG. 302 and it will be used as a starting image for the basic INNG algorithm with s = 8. A scaling factor of 8 is usually sufficient in practice to reduce data shift errors.

[0051] In the facilitated INNG algorithm, an intermediate image reconstructed using one basic INNG algorithm is used as a starting image for the next basic INNG algorithm with a larger scaling factor. Although the final basic INNG algorithm must satisfy a rigorous stopping criterion, i.e., a small value of d in Eq.[308], in order to reconstruct a high-quality image, intermediate images do not have to satisfy a small d because they are merely 'estimate images' in the next basic INNG algorithm. Therefore, relaxed stopping criteria, i.e., relatively large d, can be used for all the basic INNG algorithms, with the exception of the last, in order to further improve the computational efficiency.

[0052] When the size of a rescaled matrix is significantly large, it is often impractical to perform an FFT on such a large matrix, and hence the INNG algorithms described above are difficult to implement. To address this problem, a Block INNG (BINNG) algorithm has been developed.

BINNG Method

[0053] FIG. 303 shows an exemplary partition scheme of the BINNG algorithm. In the BINNG algorithm, the acquired k-space region is partitioned into several blocks, and the basic or facilitated INNG algorithm is applied to each block. The sampled k- space is partitioned into, for example, 3 x 3 blocks. All blocks do not need to be exactly the same size. In FIG. 303, the acquired k-space region is denoted as a square with its side length 2 |k_maχ|. Suppose that the basic INNG algorithm is applied to the shadowed block at the upper left corner in FIG. 303. When the scaling factor is s, the k-space data that are within the square with bold lines (the side length |k_maχ|) are distributed to an sN/2 x sN/2 matrix. Zero data values are assumed for the non- sampled k-space region within the bold square.

[0054] The basic INNG algorithm is applied to the data within the bold square < region using, an sN/2 x sN/2 matrix as though the original target grid matrix size is N/2 x N/2. In other words, 2D-IFT is first performed on the sN/2 x sN/2 k-space data matrix (corresponding to (a)-»(b) in FIG. 301 ), while zeros are set outside of the central N/2 x N/2 region (corresponding to (b)→(c) in FIG. 301 ). A 2D-FT is subsequently performed (corresponding to (c)->(d) in FIG. 301), and the original k- space data within the bold square region are inserted into the updated sN/2 x sN/2 data matrix (corresponding to (d)-»(e) in FIG. 301). A 2D-IFT is then performed on the sN/2 x sN/2 data matrix (corresponding to (e)→(b) in FIG. 301 ). The above procedures are repeated until the difference between the updated matrix (b) and the matrix (b) at the previous iteration becomes sufficiently small. It is evident that an incomplete image appears in the central N/2 x N/2 region in the above iterations. However, both constraints (i) and (ii) of the INNG algorithm are effectively imposed on the sN/2 x sN/2 matrix at each iteration.

[0055] As is understood from the above procedures, the facilitated INNG algorithm can also be applied to the selected k-space data set by successively increasing the scaling factor. In this case, the extracted N/2 x N/2 matrix is transferred to the center of the next larger rescaled matrix of zeros after each basic INNG algorithm is performed.

[0056] The obtained sN/2 x sN/2 data matrix may contain non-negligible errors in the regions close to the edges as the k-space data are abruptly truncated when they are selected. Therefore, only the part of the matrix that corresponds to the originally determined block (the shadowed region in FIG. 303) may be kept from the obtained sN/2 x sN/2 data matrix.

[0057] After all the 3 x 3 blocks are processed in a similar manner, an sN x sN k- space data matrix can be formed. It is expected that this data matrix satisfies both conditions (i) and (ii) for the entire region. In order to reconstruct an image, as applying 2D-FFT to the sN x sN data matrix is computationally impractical (original assumption), a 2D-FFT is performed on the N x N data matrix obtained by s-fold decimation of the sN x sN data matrix.

[0058] In the above example of the BINNG algorithm, the sampled k-space region is partitioned into the exemplary 3 x 3 blocks, and the maximum size of the rescaled matrix is reduced to sN/2 x sN/2 from sN x sN required for the INNG algorithms. Other partition schemes and block sizes are also possible. For example, when the acquired k-space region is partitioned into 5 x 5 blocks, the maximum size of the rescaled matrix can be reduced to s/V/4 x sV/4.

[0059] The image reconstruction algorithms described above have been tested by reconstruction of a computer-simulated image and from physically acquired spiral data.

Computer simulation

[0060] Referring now to FIG. 304, a 128 x128 numerical phantom was constructed and ten interleaved spiral trajectories were designed. Each trajectory sampled 1765 points consisting of simulated data and shared the central point. Noise-corrupted data were also simulated by adding Gaussian white noise to the ideal data. The mean of the noise was 0, and the standard deviation (SD) of the noise was equal to 20% of the average magnitude of the original ideal data.

[0061] Five reconstruction algorithms were used: basic INNG; facilitated INNG; BINNG; conventional gridding; and BURS (Block Uniform Re-Sampling) for both ideal and noisy data. [0062] In the basic INNG algorithm, three different scaling factors were tested: s = 2, 4 and 8. For each s, 101 iterations were performed. In the facilitated INNG algorithm, as shown in FIG. 302, the scaling factor was increased as follows: s = 2 - 4 → 8. When s = 2 and 4, the predetermined stopping criteria (d₂, d_A) were set to (0.01 , 0.01 ). However, when s = 8, no predetermined stopping criterion was set, and the iterations were performed until the total number of iterations were 101. In the BINNG algorithm, the simulated k-space region was partitioned into 5 x 5 blocks, and the facilitated INNG algorithm was performed to process each block using sN/4 x sV/4 rescaled matrices. Two stopping criteria were used: (cfe, o , d_&) = (0.01 , 0.01 , 0.005) and (0.01 , 0.01 , 0.001 ).

[0063] In the conventional gridding algorithm, the data were convolved with a Kaiser-Bessel window (width 2.5). Two DCFs were used: the Area Density Function (ADF) and the Voronoi DCF which been shown to be well-optimized DCF. Gridding was performed onto a 2/V x 2/V matrix to avoid amplification of aliased side lobes, and the central N x N image matrix was cropped after post-compensation.

[0064] In the BURS algorithm, truncated SVD was used for SVD regularization for both the ideal and noisy data. As truncation schemes, Lp (the total number of singular values)] the largest singular values were used in our computation, where p was a constant ranging from 0 to 1 and L«J denoted the maximum integer less than •. Values of p ranging from 0.2 to 1.0 with a 0.1 interval were tested for both the ideal and noisy data (i.e. nine images were reconstructed for each data). The radii of the data window (3κ) and the rectilinear grid window (Δk) were set to 2.0 and 4.5 (in Cartesian step), respectively.

[0065] For all the reconstructed images from both ideal and noisy data, the root mean square (RMS) errors from the original numerical phantom were measured. For all the reconstructed images from the noisy data, the image SNR (Sj_m/N_im) was also measured where S_/m and Λ/_/m were the mean signal amplitudes of the 25 x 35 pixels in the upper right white rectangle and the SD of 15 x 10 pixels in the upper middle black rectangle, respectively. In the basic and facilitated INNG algorithms, the RMS error and the image SNR (for noisy data only) were measured for the image at each iteration. [0066] FIG. 305 shows the profiles at the 34th row and the 102nd column of the images (indicated by the dashed lines in FIG. 304) with ideal data reconstructed using the basic INNG algorithm with s = 8 after 30 iterations (a) and that after 101 iterations (b), the facilitated INNG algorithm after 30 total iterations (c), the BINNG algorithm with (d₂, d₄, d₈) = (0.01 , 0.01 , 0.001 ) (d), the conventional gridding algorithm with Voronoi DCF (e), and the BURS algorithm with p = 0.8 (f). As can be seen, FIGs. 305 (b)-(e) show little deviation from the original profiles. In FIG. 305 (a), the high spatial frequency components are not sufficiently reproduced. The BURS algorithm reconstructed the image with the least RMS error when p = 0.8 of all the nine tested values of p. Although no significant profile distortions were observed in FIG. 305 (f), the deviations from the original profiles are larger than those of FIGs. 305 (b)-(e). FIG. 306 shows the RMS error of the image with ideal data at each iteration in the basic INNG algorithm (the solid, dashed and dotted lines correspond to s = 2, 4 and 8, respectively) and facilitated INNG algorithm (the dash-dot line). The same types of lines are used in Fig. 307. Note that Fig. 306 is log-scaled. As Fig. 306 shows, in the basic INNG algorithm with s = 2, the RMS error decreases for the first few iterations and then increases although the rate of increase is decreased with iterations. Similarly, when s = 4, the RMS error is decreased to approximately 30 iterations and then remains almost unchanged even though iteration progresses. When s = 8, the RMS error is monotonically decreased at least until the 101st iteration. As expected, the larger s is, the smaller the RMS error becomes after a sufficient number of iterations. The RMS error of the facilitated INNG algorithm is decreased rapidly compared with that of the basic INNG algorithm when s = 8 and becomes almost constant after 30 iterations. [0067] The RMS error of each algorithm with ideal data is summarized in Table 301 where the number after the # sign represents the total number of iterations performed (also in Tables 302 and 303). For example, s = 2 (#7) denotes that seven iterations were performed when s = 2. The number of iterations indicated in the BINNG algorithm is the average number of iterations for all 25 blocks for the specified s (this is also the case in Tables 302 and 303). As seen in Table 301 , the INNG/BINNG algorithms involving s = 8 reconstruct the images with comparable or even smaller RMS errors than the conventional gridding algorithm with Voronoi DCF and the BURS algorithm. The RMS energy of the original numerical phantom was 48.581 x 10^'2.

[0068] FIG. 307 shows the measured image SNR at each iteration in the basic and facilitated INNG algorithms. The image SNR is transiently increased for initial parts of the iterations when s = 4 and 8. However, after certain points, the image SNR is monotonically decreased with iterations for all the basic INNG algorithms. Furthermore, the rate of decrease becomes more insignificant with iterations. In the facilitated INNG algorithm, the image SNR becomes almost constant after approximately 30 iterations. [0069] Table 302 summarizes the RMS error and the measured image SNR of each algorithm with noisy data. For noisy data, the BURS algorithm reconstructed the image with the least RMS error when p = 0.5 of all the nine tested values of p. In Table 302, the images reconstructed using the facilitated INNG and BINNG algorithms involving s = 8, show equivalent RMS errors and image SNR to that of the conventional gridding algorithm with Voronoi DCF.

[0070] The image reconstructed using the basic INNG algorithm with s = 8 after 30 iterations has the smallest RMS error and the largest image SNR in Table 302. However, as is evident from Fig. 307, this image is still in a transient stage of reconstruction in the basic INNG algorithm with s = 8. [0071] The INNG/BINNG algorithms described above lead to accurate image reconstruction when s is sufficiently large. The RMS errors of the facilitated INNG and BINNG algorithms involving s = 8 are comparable to those of conventional gridding with Voronoi DCF, which has been shown to reconstruct images with quite small errors. [0072] Although the BURS algorithm has been reported to be quite accurate, the deviation observed in FIG. 305 (f) is more pronounced than those seen in Figs. 305 (b)-(e). Since the BURS algorithm employs SVD matrix inversion, it is often unstable for over-sampled regions, even if the k-space data are ideal (This is also true for under-sampled regions. However, the sampling schemes employed in the context of the disclosed inventive methods satisfy the Nyquist criterion for the entire k-space.). Spiral trajectories usually sample the central k-space region more densely than the outer regions, which can result in relatively large gridding errors in the lower frequency regions. Although SVD regularization is effective for reducing these artifacts, it is difficult to determine appropriate regularization schemes before the algorithm is started. For this reason, the quality of the BURS algorithm is not as good as that of other comparable algorithms for the examples presented herein. [0073] Equation [7] suggests that the image reconstructed using the basic INNG algorithm can be regarded as the summation of the image originated from the ideal signal components l_m,ideai (x) and the image from the residual imperfections n_m(x). In other words, as iteration progresses, the RMS energy of [l_m,i_deai (x) - 1 (x)] is continuously decreased toward zero, and that of n_m(x) is increased. Thus, the initial descending part of each RMS error curve in Fig. 306 predominantly reflects the RMS energy of [l_m deai (x) - I (x)], and the RMS energy of n_m(x) becomes increasingly dominant in the RMS error with iterations. However, as mentioned earlier, the increase rate of the RMS energy of n_m(x) is reduced as the iteration progresses, since Eq.[306] holds for n_m(x) as well. This fact is also suggested from Fig. 306. As is noted in Fig. 307, the measured image SNR is almost unchanged for s = 2 and 4 after a certain number of iterations. hese results suggest that the image quality does not significantly change even if iterations unnecessarily continue from initial noise-corrupted data.

[0074] Table 302 shows that both RMS error and image SNR of the basic INNG algorithm with s = 4 are equivalent to those of the facilitated INNG/BINNG algorithms involving s = 8. A comparison between Figs. 305(a) and 305(b) suggests that a larger number of iterations are necessary for the higher frequency components to appear in the image for the basic INNG algorithm. This is also the case with noisy data. The image reconstructed using the basic INNG algorithm with s = 8 after 30 iterations has the least RMS error and the highest SNR in Table 302. However, fine structures may not be sufficiently reproduced. The expression of the fine structures in the image may require more iterations, although the image SNR becomes lower, as shown in Fig. 307. It should be noted that the RMS error is one of many metrics used to evaluate the image quality when the reference image is available. [0075] As can be seen in Fig. 306, the facilitated INNG algorithm substantially reduces the number of iterations from the basic INNG algorithm for the same target scaling factor s. In the facilitated INNG algorithm, intermediate images are used as starting images for the basic INNG algorithm with the next larger scaling factor, as shown in Fig. 302. As these images are only rough estimates for the next basic INNG algorithm, the number of iterations for the intermediate INNG algorithms could be reduced to further improve the computational efficiency. The employed stopping criteria (d₂, d ) = (0.01 , 0.01 ) are relatively relaxed. The numbers of iterations performed under this condition were 7 and 5 for s = 2 and 4, respectively, for the ideal data. However, Fig. 306 suggests that seven iterations for s = 2 may be unnecessary as the RMS error is minimized before 7th iteration when s = 2. These examples suggest that the stopping criteria used to obtain intermediate images should be as relaxed as possible to avoid unnecessary computations, although it is difficult to determine the optimal stopping criteria in advance.

[0076] The BINNG algorithm is useful when the target rescaled matrix is quite large and hence it is computationally impractical to perform a 2D-FFT. The BINNG algorithm was applied to the same simulated data used for the INNG algorithms in order to compare both algorithms. The image quality of the BINNG algorithm is comparable to that of the INNG algorithms, as seen in Fig. 305(d). The basic concept behind the BINNG algorithm is to recover the data that satisfy the constraints (i) and (ii) in each k-space block region. Therefore, the size of the blocks can be determined without any restriction. However, it would be desirable to have the central region of k-space covered by one complete block to avoid the data discontinuity between the blocks near the k-space center. In our examples, we set 5 x 5 blocks, and each block was processed using s/V/4 x s/V/4 matrices. As this example indicates, it is suggested that odd-by-odd blocks be used in the BINNG algorithm for the reason mentioned above. Physical scan

[0077] Images were reconstructed from physically acquired spiral data using the methods described above. In one exemplary experiment, an axial brain image was scanned from an asymptomatic volunteer. The acquisition was performed using a 1.5-Tesla Siemens Sonata scanner (Siemens Medical Solutions, Erlangen, Germany). All procedures were performed under an institutional review board- approved protocol for volunteer scanning. [0078] The followings acquisition parameters were used: 20 spiral interleaves; FOV 240 x 240mm; slice thickness 10mm; spiral readout time 16.0ms; and TE/TR=6.0/33.0ms. 1-4-6-4-1 binomial pulses were used for spatial-spectral excitation. The total flip angle for on-resonance spins was 16°. [0079] The images were reconstructed using basic INNG, facilitated INNG, BINNG, a conventional gridding algorithm with ADF and with Voronoi DCF, and BURS with p = 0.8. In the basic INNG algorithm, s = 2 and d₂ = 0.001. In the facilitated INNG algorithm, s = 2 → 4 and (d₂, d ) = (0.01 , 0.001 ). In the BINNG algorithm, the acquired k-space was partitioned into 5 x 5 blocks. Each block was processed by the facilitated INNG algorithm using s/V/4 x s/V/4 rescaled matrices with s = 2 → 4 -> 8 and (d₂, d_A, d₈) = (0.01 , 0.01 , 0.003).

[0080] For each reconstructed image, the image SNR (Si N_im) was measured where Sj_m and N_ιm were the mean signal amplitude of the 15 x 15 pixels in the white matter and the SD of 50 x 50 pixels in the background, respectively. Their locations are indicated by the squares in Fig. 308 (a).

[0081] FIG. 308 shows the in-vivo brain images reconstructed using the facilitated INNG algorithm (a), the BINNG algorithm (b), the conventional gridding algorithm with Voronoi DCF (c) and the BURS algorithm with p = 0.8 (d). There are visually no significant differences among these images. Table 303 summarizes the measured image SNR for each reconstructed image. There is no significant image SNR difference among the images shown in Fig. 308. The SNR of the basic INNG algorithm with s = 2 is equivalent to that of the conventional gridding with ADF, and the SNR of these two algorithms are approximately 60% of those of the other four algorithms shown in Table 303. [0082] The maximum scaling factors used to reconstruct the images of Figs. 308(a) and 308(b) were 4 and 8, respectively. However, there are no visually significant differences between these two images. This example suggests that scaling factor 4 would be sufficient for typical image contents as was also indicated by the results of the noisy simulated data shown in Table 302. [0083] The conventional gridding is robust to noise because convolution essentially averages the local k-space data. This is generally true independently of the type of DCF. However, as indicated in Table 303, the image SNR of the conventional gridding with ADF is lower than those of other algorithms whose images are shown in Fig. 308. This resulted from the non-negligible profile distortions of the image reconstructed using the conventional gridding with ADF. In other words, the signal amplitudes of the scanned object regions were overcompensated in the post- compensation step; hence the signal S_/m was measured low in this image relative to that in other images.

[0084] The INNG/BINNG method of the invention is quite simple and does not require complicated procedures to compute DCFs, while reconstructing images with small degrees of error. One primary drawback of the INNG/BINNG algorithms is that it is difficult to set appropriate stopping criteria for iterations. Exemplary stopping criteria employed herein are based on the metric d in Eq.[308]. Although the value of d is monotonically decreased for a specific s, it is quite difficult to determine in advance the stopping criteria for d. For example, when an image was reconstructed from the in-vivo data using the facilitated INNG algorithm, we set the stopping criteria (d , d ) = (0.01 , 0.001). These values were determined based on the simulation experiments. As FIG. 307 shows, the rate of image SNR change becomes smaller as iteration progresses. This fact suggests that stopping criteria can be alternatively devised by utilizing the rate of the image SNR change. However, it would be still difficult to set appropriate values to stop the iterations before the algorithms are started. In practice, the number of iterations to be performed may need to be optimized empirically.

[0085] In the conventional gridding and the BURS algorithms, interpolation coefficients can be calculated once the coordinates of the sampled points are given. In other words, the interpolation coefficients can be pre-calculated and stored and subsequently multiplied by k-space data after the data are acquired. This process facilitates the reconstruction speed after data acquisition. Generally, the procedures of the INNG/BINNG algorithms are performed after the k-space data are acquired, since the they employ iterations. [0086] According to another practice, partial Fourier reconstruction techniques can be employed to reduce scan time in spiral MR sampling schemes. This technique employs variable-density spiral (VDS) trajectories so that the Nyquist criterion is satisfied in the central region of the k-space, whereas the outer regions of k-space are undersampled. The projections onto convex sets (POCS) method is used in the reconstruction. The disclosed partial Fourier spiral reconstruction (PFSR) technique permits reduced scan time when compared with the conventional spiral imaging. [0087] Spiral imaging is a rapid MRI data acquisition technique that has gained considerable popularity, particularly in dynamic imaging, in the past few years. Although its scan time is usually a fraction of a second, artifacts due to sufficiently rapid motion are often observable. In the technique disclosed herein, a rectilinear partial Fourier reconstruction technique has been extended to image reconstruction from undersampled spiral k-space data sets, wherein the k-space data are incompletely sampled to further reduce the scan time of spiral imaging.

[0100] The PFSR technique applies the projection onto convex sets (POCS) method developed in rectilinear sampling schemes (3), to spiral sampling schemes. There are two constraints that are imposed on the data set at each iteration in the POCS method: (i) phase constraint in the image domain, (ii) data- consistency constraint in the k-space domain. To apply this method to spiral sampling schemes, the rescaling matrix reconstruction algorithm (the equivalent algorithm was proposed as the 'next neighbor re-gridding algorithm') has been modified. The first step of this algorithm is to create an estimated image phase map Φ_e from the low-resolution image reconstructed from the central k-space data. The next step is to perform iterative procedures to impose the two constraints on the acquired data set.

[0101] The proposed PFSR algorithm follows essentially the flow of the basic INNG algorithm described above with reference to FIG. 301, except that a phase constraint is imposed on the image (c) of FIG. 301. The original target grid is an N x N matrix. K-space data are distributed into a larger rescaled matrix by a factor of s (= 2m), where m is a small positive integer (FIG. 301(a)). The location of each datum in the large rescaled matrix is determined by multiplying the original k-space coordinate by s and then rounding the rescaled coordinate off to the nearest target rectilinear grid location. If more than one datum share the same matrix coordinate, the mean value is stored. An Inverse Fourier Transform (IFT) is performed on matrix (a), leading to image matrix (b). The intermediate reconstructed image appears in the central N x N matrix in (b). The phase constraint is imposed on the central N x N image region. That is,

I_new = * o iΦ_e) , [309]

where l₀i_d and l_new represent the image values at each pixel in the central N x N region of (b) before and after the phase constraint, respectively. The region outside of the central N x N matrix is set to zeros, resulting in (c). An FT is performed on (c), leading to (d) which is an estimate of the phase-constrained raw data. Then, the data-consistency constraint is imposed on this data matrix, i.e., the data where the original data exist are replaced by the original data values, as shown in (e). An IFT is performed on (e). The updated reconstructed image again appears in the central N x N matrix (b). The procedures (b) -> (c) - (d) - (e) - (b) (surrounded by dashed lines in FIG. 301) are repeated until the difference between the updated image and the image at the previous iteration becomes sufficiently small.

[0102] In-vivo spiral imaging experiments were performed in the manner described above to demonstrate the performance of the proposed algorithm. Head images were scanned from a healthy volunteer using a 10 interleaved VDS. Each spiral trajectory has 4.5 turns in the central portion (radius 30% of k_max: the Nyquist criterion is satisfied in this region) and 5.25 turns in the outer regions of k-space (the Nyquist criterion is violated). TE/TR=6.0/31.0ms, the slice thickness 10mm. For image reconstruction N=256 and s=4. The PFSR algorithm was performed with 50 iterations. A conventional matrix rescaling algorithm was also used as a standard reconstruction method for comparison. A comparative data set was also acquired with 20 interleaved constant density spiral trajectories (7.5 turns each). In the latter trajectory, the Nyquist criterion is satisfied for the entire k-space.

[0103] Fig. 309 shows the reconstructed images a, b: via reduced data sets ((a) via the matrix rescaling algorithm, (b) via the new PFSR algorithm, (c) via a full data set). As observed, (a) is affected by aliasing artifacts. The aliasing artifacts are reduced in (b) without loss of spatial resolution.

[0104] The reconstructed image quality in the conventional rectilinear partial Fourier reconstruction with POCS has been shown to depend on the estimated phase, which is also the case with the PFSR algorithm discussed above. The variable-density spiral can sample the central region of k-space with little additional acquisition time as compared with a constant-density spiral. In other words, the estimated phase map can be efficiently obtained with the use of a VDS in the PFSR technique. Constraint (ii) is difficult to apply when k-space data are sampled non-uniformly. However, the PFSR algorithm can overcome this difficulty since it uses large rescaled (i.e., rectilinear K-space) matrices. Both constraints (i) and (ii) can be readily imposed on the data set at each of the iterations depicted in FIG. 301 with the phase constraint imposed in (c). The PFSR technique permits image reconstruction with reduced artifacts from undersampled spiral data sets, thereby enabling further reductions in scan time in spiral imaging.

[0105] According to another practice, a Block Regional Off-Resonance

Correction (BRORC) can be employed as a fast and effective deblurring method for spiral imaging. Spiral acquisition techniques have advantages over other k- space trajectories because of their short scan time and insensitivity to flow artifacts, but suffer from blurring artifacts due to off-resonance effects. A frequency-segmented off-resonance correction (FSORC) method is commonly used to combat off-resonance effects and reconstruct a deblurred image. In this algorithm, several k-space data sets are first created by demodulating the original data by several different frequencies; separate images are reconstructed from each demodulated k-space data set via 2D inverse Fourier Transform (IFT). Deblurred image regions are selected from the reconstructed images under guidance of a frequency field map. The final reconstructed image with off- resonance correction is created by combining all deblurred regions selected from the appropriate demodulated image. The computational burden of FSORC is proportional to the number of demodulation frequencies used since the fast Fourier transform (FFT) is performed on each demodulated k-space data set.

Hence, FSORC is often computationally intensive, particularly when a wide range of off-resonance frequencies exists across a scanned object.

[0106] Other off-resonance correction algorithms with improved computational efficiency use, for example, a linear field map which, however, can only correct for linear components of off-resonance frequency variation. Therefore, residual frequency variations that deviate from the linear variation must be corrected with FSORC; hence several FFTs are usually required. In another conventional method called multi-frequency interpolation (MFI), images are reconstructed using a reduced number of demodulation frequencies. Images requiring other demodulation frequencies are estimated from the limited set of demodulated/reconstructed images via interpolation. In MFI, the interpolation coefficients need to be pre-calculated. The total number of demodulation frequencies used in MFI is typically one-fourth to one-third that of the conventional FSORC. Image domain deconvolution methods approximate the spiral time evolution function as a quadratic function with respect to a k-space radius. This enables correction via one-dimensional deconvolution (along the x and y directions) in the image domain since separable demodulation functions along the x and y directions can be formed. However, image quality degradations beyond those associated with FSORC may result when the difference between the actual spiral time evolution function and the approximated quadratic function cannot be ignored.

[0107] Accordingly, to improve image reconstruction, a novel fast off-resonance correction method, (a.k.a., 'Block regional off-resonance correction (BRORC)') is presented. In this method, off-resonance correction proceeds block-by-block through the reconstructed image, and FFTs are performed on matrices that are smaller than the full image matrix. Although the computational cost of BRORC relative to that of FSORC depends on the selection of the parameter values in these algorithms, the BRORC is usually computationally more efficient than FSORC. Furthermore, greater reduction of the computational costs can be expected in BRORC if only particular regions of the image need to be deblurred.

[0108] Referring now to FIG. 310, a block diagram of BRORC is shown having an original image matrix size ot Nx N (e.g., 256 x 256). The first step of the BRORC is to extract a small block region Mx M. For convenience, M is typically chosen to be a number expressed as a power of 2 (e.g., 16, 32), though this need not be the case. A 2D-FFT is performed on the Mx M image matrix. The obtained Mx M Fourier data is to be frequency demodulated. The demodulation function matrix for the Mx M data must also be M x M in size. This matrix can be obtained by /W-fold decimation of the original Nx N demodulation function matrix. Regions near the four corners of the Mx M demodulation function matrix should be handled carefully. Normally, after the acquired spiral k-space data are gridded onto an Λ/x /V grid, there are no data outside of the inscribed circle (radius N/2 in Cartesian step). These regions are usually set to zeros in the Nx N data matrix before frequency demodulation is performed. However, in the M x M Fourier data matrix, all the M M matrix elements usually have non-zero data values. If the corresponding M M demodulation frequency matrix has zero values in the regions near the four corners, artifacts originating from the inaccurately demodulated high spatial frequency components may appear after demodulation. Therefore, when the Mx M demodulation function matrix is formed, the regions outside the inscribed circle are filled with the maximum readout time values, thereby effectively performing /V/M-fold decimation without introducing such artifacts.

[0109] The demodulation frequency ('f indicated in Fig. 310) is determined from the central region of the Mx M sub-image matrix in the frequency field map. In practice, the mean off-resonance frequency of the central rM x rM pixels (0 < r ≤ 1 , and r is typically 0.5.) in the Mx M phase image matrix, is used as the demodulation frequency 'f. After frequency demodulation, the Mx M k-space data is subsequently 2D-inverse Fourier transformed. Since the outer regions of the obtained Mx M image matrix may exhibit artifacts, only the central rM x rM pixels of the Mx M deblurred image matrix are kept for the final reconstructed image. This procedure is repeated until the entire scanned object is deblurred. However, as is evident from the BRORC block diagram, it is also possible to only deblur particular regions of the image. This is not possible with the conventional FSORC

[0110] In-vivo spiral images were acquired to facilitate comparison of FSORC and BRORC. All acquisitions were performed using a 1.5-Tesla Siemens Sonata scanner (Siemens Medical Solutions, Erlangen, Germany). Axial brain images and cardiac images were acquired from asymptomatic volunteers using a quadrature head coil and four-element phased array surface coils, respectively. All procedures were performed under an institutional review board-approved protocol for volunteer scanning.

[0111] For the brain image acquisition, 20 spiral interleaves were used with a field of view (FOV) set to 240 x 240mm (resolution 0.94mm); slice thickness 10mm; spiral readout time 16.0ms; and TE/TR=6.0/1000.0ms.

[0112] ECG gating was used during cardiac image acquisition. There were 16 spiral interleaves; FOV 310 x 310mm (resolution 1.21mm); slice thickness 7mm; spiral readout time 14.0ms; and TE 6.0ms. [0113] For both brain and cardiac image acquisitions, T2 prep pulses were used prior to each acquisition and 1-4-6-4-1 binomial pulses were used for spatial- spectral excitation. The total flip angle for on-resonance spins was 64°. For each of the experiments, another image was acquired using the same sequence but with TE = 8.0ms. The frequency field map was obtained from the phase difference between these two reconstructed images.

[0114] Head images were reconstructed by gridding k-space data onto a 256 x 256 Cartesian grid using a modified Block Uniform Resampling (BURS) algorithm. Cardiac images were reconstructed via the matrix rescaling algorithm to facilitate the reconstruction from multiple coils' data (The NNG algorithm described above represents an equivalent algorithm). The reconstructed image matrix was 256 x 256 in size. The cardiac images were reconstructed via the sum-of-squares method from data acquired in each element of the phased-array torso/body coil. [0115] Both FSORC and BRORC were performed on each image data set for comparison. For FSORC, the total number of demodulation frequencies L is normally set to satisfy:

L > ^4Aω-^T , π [10]

where Δ<y_max is the absolute value of the maximum off-resonance frequency (in radians) and T is the spiral readout time. Eq.[310] was derived under the assumption that off-resonance frequency frequencies ranged from -Δω_max to

Δβ>_m3X . In practice, the range of off-resonance frequencies is often asymmetric with respect to the on-resonance frequency (i.e., 0 Hz off-resonance frequency). Thus, L was determined from the following modified equation:

L > ^{2 ' 2 /max ~ min ) ' r} = 4( _max -f_mia)T , [311] π

where _max and ^ represent the maximum and minimum off-resonance frequencies (in Hz) indicated in the frequency field map. L was set to the minimum integer that fulfilled Eq.[311].

[0116] Three combinations of (M, ή were tested for the BRORC algorithm in head imaging trials: (M, r) = (8, 0.25), (16, 0.5), and (32, 0.5). The BRORC algorithm was applied only to object regions in the head image after the background was properly removed by thresholding.

[0117] For the BRORC in cardiac imaging trials, (M, r) was set to (32, 0.5). In the cardiac image, only a 160 x 160 matrix centered on the heart was processed with BRORC.

[0118] A comparison of the computational costs of BRORC and FSORC shows that the total number of complex multiplications required for BRORC can be expressed as:

B = 2N² log₂ N+s - 2 -2M² log₂ M + s -M² , [312]

where s is the total number of rMx rM blocks that cover the scanned object regions. For example, if the entire Nx N image matrix is processed with BRORC,

N s = [313] rM

[0119] In Eq. [312], the first, second, and last terms represent the number of complex multiplications required for an Nx N 2D-FFT, those tor Mx M 2D-FFTs, and those necessary for frequency demodulations for M x M Fourier data, respectively.

[0120] On the other hand, the total number of complex multiplications necessary for FSORC can be expressed as:

S = L - 2N² log₂ N + L -N² = LN² (2log₂ N + 1) , [314]

where the first and second terms express the total number of complex multiplications required for an Nx N 2D-FFT and those for frequency demodulations tor Nx N k-space data.

[0121] Fig. 311 shows the axial brain images (a: the image before off-resonance correction; b: the image with FSORC; and c, d: the images with BRORC (c: (M, r) - (8, 0.25); and d: (M, r) = (32, 0.5)). The off-resonance frequencies ranged from -128.2 Hz to 46.7 Hz. Therefore, L was set to 12 in the FSORC according to Eq.[311]. Blurring artifacts can be observed for the anterior parts of the temporal lobes and at the contours of the brain stem in (a). These blurring artifacts are reduced in all images after off-resonance corrections (b)-(d). However, in (c), grid-like artifacts can be seen in some regions, including areas near the eyes (indicated by arrows) and parts of the brain. The quality of the image using BRORC with (M, r) = (16, 0.5) (this image is not shown) is almost equivalent to that of (b) although grid-like artifacts are still apparent in the eyes. The grid-like artifacts are not apparent in (d). Image (d) is visually identical to (b).

[0122] Fig. 312 presents the results of the off-resonance correction strategies on cardiac images (a: the image before off-resonance correction; b: the image with FSORC; and c: the image with BRORC (M, r) = (32, 0.5)). The off- resonance frequencies ranged from -139.1 Hz to 124.5 Hz. Hence, L was set to 15 in FSORC according to Eq.[311]. The regions indicated by arrows in (a), the definitions of a papillary muscle and an aorta, are improved in (b) and (c) when compared to the uncorrected image (a). There are no observable difference between (b) and (c).

[0123] Table 304 summarizes the total number of complex multiplications required for off-resonance correction in our experiments. Note that BRORC was applied only to the scanned object regions and not to the background in the brain images. The numbers in parentheses in the fourth column in Table 304 indicate the total number of complex multiplications if the entire 256 x 256 image matrix has to be processed with the same parameters (M, r) for each BRORC algorithm. Also note that in Table 304 the values for Figs. 312(b) and 312(c) are those required to process the data from a single coil. The variable s in Eq.[312] (the total number of rM x rM blocks processed using BRORC) is also shown in Table 304. The computational costs of BRORC with (M, r) = (16, 0.5) and BRORC with (A , r) = (32, 0.5) (FIG. 311 (d)) relative to that of FSORC with L = 12 (FIG. 311 (b)) are 21.2 % and 27.3 %, respectively. The computational cost of BRORC with (M, r) = (32, 0.5) (FIG. 312 (c)) relative to that of FSORC with L = 15 (FIG. 312 (b)) is 19.2 %.

[0124] The BRORC algorithm is usually computationally more efficient than FSORC even though the comparison depends on the parameter values in Eqs. [312,314]. For example, if r is small, a significant number of Mx M 2D-FFT's must be performed with the BRORC. Under these conditions, the BRORC may be computationally more intensive than FSORC. Also, if the range of off- resonance frequency across a scanned object is relatively narrow, i.e. the object is almost on-resonance, the total number of demodulation frequencies L would be small in FSORC. Under these conditions, the relative computational efficiency of BRORC to FSORC may be less than in the examples discussed above. However, the values shown in Table 304 represent the typical computational costs of BRORC for a 256 x 256 matrix image, and these values are independent of the range of the off-resonance frequency. In other words, when the regions of interest are approximately 40% of the entire FOV as seen in our images, BRORC always requires reduced computational demand than

FSORC when L is greater than 4 (because ^B^^M'^r'^s^ ⁼ (³²>⁰-5,121)) _ ₃__65x10 ^β ^a S(L = Ϊ)

/ [0125] One way to reduce the computational cost of BRORC is to decrease M, as suggested in Eq.[312]. However, as FIG. 311 (c) shows, when M is 8, grid- pattern artifacts can be observed although rM was set to the minimum possible value of 2. Therefore, a better value for M may be at least 16. Another way to improve the computational efficiency in BRORC for a given M is to increase ras it leads to the reduction of the total number of blocks to be processed (s in

Eq.[312]). However, r cannot be set too large as grid-like artifacts often remain, as explained below. In the examples described herein, rwas set to 0.5 in FIG. 311 (d) and FIG. 312 (c) when M = 32. Sufficient computational advantages were obtained in these BRORC algorithms compared to FSORC, as seen in Table 304.

[0126] When an Mx M image matrix is extracted to be Fourier transformed, the matrix is abruptly truncated. Thus, artifacts often appear in the outer regions in the Mx M image matrix after frequency demodulation. These artifacts are enhanced when the extracted matrix is located in regions with large off- resonance frequencies (in absolute values). To avoid these artifacts in the final

M — rM reconstructed image, the width of the outer regions to be discarded ( )

should be reviewed. For example, three pixels were discarded from the edges of the 8 x 8 matrices in FIG. 311 (c). However, grid-like artifacts appeared in the image. This example suggests that three pixels are insufficient to avoid apparent grid-like artifacts. The artifacts can be significantly reduced when four pixels are discarded from the edges of 16 x 16 matrices (This image is not shown). In FIG. 311 (d), eight pixels were discarded from the edges of 32 x 32 matrices and no evident grid-like artifacts can be observed. Since the width of the spiral blurring PSF is roughly given by (4*ADC time (in sec)*off-resonance frequency in Hz), the maximum blurring extent was estimated as 4 • 0.016 -128 = 8.2 pixels in the < depicted brain images. Although it slightly exceeds eight pixels, the residual artifacts are visually not apparent in Fig. 311 (d). Note that Fig. 312 (c) in which the same parameter was used as in Fig. 311 (d), also shows no grid-like artifacts. In general, the width of the discarded outer, region must be increased in areas of larger off-resonance effects to avoid grid-like artifacts.

[0127] In BRORC, a constant off-resonance frequency is assumed within each rMx rM block. Since a frequency field map is usually smoothly varying across the FOV, this assumption is, in general, valid. If abrupt frequency transitions occur in some regions, these regions may not be accurately demodulated. For this reason, rM may need to be set small to achieve accurate demodulation. Or, r may need to vary over the image based on the rate of change of the frequency field map. However, as seen in Fig. 311 (d) and Fig. 312 (c), it appears that rM = 16 is an acceptable size in these particular examples to achieve effective off- resonance correction.

[0128] Off-resonance effects are proportional to the product of the spiral readout time and the off-resonance frequency. Although the spiral readout time can be varied for each experiment and each clinical application, the spiral readout times used in the experiments described above were typical values for many spiral imaging experiments. The off-resonance frequency ranges are also comparable to those expected in clinical in-vivo imaging experiments. Therefore, the above discussions would be useful when choosing appropriate parameters in the

BRORC algorithm. In summary: for a 256 x 256 image, (i) To achieve effective off-resonance correction, rM= 16 is acceptable in practice even though a smaller rM is desirable, (ii) To reduce the computational costs of BRORC, r should not be

M — rM set too small, (iii) To avoid grid-like artifacts, > 4 is suggested. Based

Σ on the above results, (M, r) = (32, 0.5) is likely to be appropriate for BRORC for a 256 x 256 image in many clinical settings.

[0086] Accordingly, the BRORC algorithm is quite simple and it produces reconstructed image quality comparable to that using the FSORC. BRORC is typically computationally several times more efficient than the FSORC with no perceptual difference between the images. Moreover, BRORC can be applied to particular regions of interest to further reduce computational requirements.

[0087] While the invention has been disclosed in connection with the preferred embodiments shown and described in detail, various modifications and improvements thereon will become readily apparent to those skilled in the art.

[0088] What is claimed is:

Table 301. Algorithm comparison with ideal simulated data

a: The RMS energy of the original numerical phantom is about 48.581 x 10^"

Table 302. Algorithm comparison with noisy simulated data

a: The RMS energy of the original numerical phantom is 48.581 x 10^"

Table 303. Algorithm comparison with in-vivo data

Table 304. Comparison of computational costs of FSORC and BRORC

Claims

CLAIMS:

1. An MR coil for intravascular imaging, comprising:

a. a first coil wound in a first direction;

b. a second coil at a distance from, and coaxially arranged with respect to, the first coil, the second coil wound in a direction substantially opposite the first direction; and

c. at least one wire connecting the first coil and the second coil,

wherein the at least one wire extends substantially longitudinally, and is disposed substantially in the center, between the first and * the second coils.

2. The MR coil of claim 1 , wherein the distance is between about two times and about three times a winding diameter of the coil.

3. An MR coil for intravascular imaging, comprising:

a. a first coil wound in a first direction;

c. wires connecting the first coil and the second coil,

wherein the connecting wires extend mutually proximally, substantially collinearly, and essentially longitudinally along an imaginary line connecting a periphery of the first and second coils.

4. An MR coil for intravascular imaging, comprising: a. a first coil wound in a first direction; and

b. a second coil at a distance from, and coaxially arranged with respect to, the first coil, the second coil wound in a direction substantially opposite the first direction,

wherein each coil is connected to a distinct receive channel that is independently accessible during an MR imaging procedure.

5. The MR coil of claim 4, wherein the distance is between about two times and about three times a winding diameter of the coil.

6. The MR coil of claim 4, wherein each coil is individually frequency- tuned.

7. The MR coil of claim 4, wherein respective signals from the first and second coils are processed separately.

8. The MR coil of claim 7, wherein at least a portion of the processed signals are used for localizing the coil.

9. The MR coil of claim 7, wherein at least a portion of the processed signals are used for vascular imaging.

10. A wireless MRI tracking probe comprising a first material having an MR resonance frequency distinct from a resonance frequency of a second material adjacent to the first material.

11. The probe of claim 10, wherein the first material comprises a liquid contained in a lumen.

12. The probe of claim 11 , wherein the liquid comprises at least one of acetic acid and boric acid.

13. The probe of claim 11 , wherein the liquid comprises a contrast material.

14. The probe of claim 13, wherein the contrast material comprises a rare earth compound.

15. The probe of claim 11 , wherein the liquid is labeled with an isotope.

16. The probe of claim 10, wherein the first material comprises a biocompatible material.

17. The probe of claim 16, wherein the biocompatible material includes a material selected from the group consisting of: fluorinated compound, fluorocarbons, azo-compounds, amino acids, glucose, and a combination thereof.

18. The probe of claim 16, wherein biocompatible material is labeled with an isotope.

19. The probe of claim 10, wherein the second material includes at least one of water and fat.

20. The probe of claim 10, wherein the first material has a selectivity greater than about five times a selectivity of the second material.

21. The probe of claim 10, wherein the selectivity of the first material is greater than about ten times the selectivity of the second material.

22. The probe of claim 10, wherein the selectivity of the first material is greater than about twenty times the selectivity of the second material.

23. A catheter having a wireless MRI tracking probe, wherein the probe comprises a first material having an MR resonance frequency distinct from a resonance frequency of a second material adjacent to the first material.

24. The catheter of claim 23, wherein the first material comprises a liquid contained in a lumen.

25. The catheter of claim 24, wherein the liquid comprises at least one of acetic acid and boric acid.

26. The catheter of claim 24, wherein the liquid comprises a contrast material.

27. The catheter of claim 26, wherein the contrast material comprises a rare earth compound.

28. The catheter of claim 24, wherein the liquid is labeled with an isotope.

29. The catheter of claim 23, wherein the first material comprises a biocompatible material.

30. The probe of claim 29, wherein the biocompatible material includes a material selected from the group consisting of: fluorinated compound, fluorocarbons, azo-compounds, amino acids, glucose, and a combination thereof.

31. The probe of claim 29, wherein biocompatible material is labeled with an isotope.

32. A method of tracking a catheter having a wireless MRI tracking probe, wherein the probe comprises a first material having an MR resonance frequency distinct from a resonance frequency of a second material adjacent to the first material, the method comprising:

a. applying a conventional FLASH sequence to excite the MR frequency, to acquire an anatomic image of the second material, b. applying a CHESS pulse to excite the MR frequency of the first material, to acquire a catheter-selective image; and

c. alternating between the anatomic image and the catheter- selective image to track the catheter.

33. A method of reconstructing a magnetic resonance image from non- rectilinearly-sampled k-space data, comprising:

a. distributing the sampled k-space data on a rectilinear k-space grid;

b. inverse Fourier transforming the distributed data;

c. setting to zero a selected portion of the inverse-transformed data;

d. Fourier transforming the zeroed and remaining portions of the inverse-transformed data;

e. at grid points associated with the selected portion, replacing the Fourier-transformed data with the distributed k-space data at corresponding points of the rectilinear k-space grid, thereby producing a grid of updated data;

f. inverse Fourier transforming the updated data; and

g. until a difference between the inverse Fourier-transformed updated data and the inverse Fourier-transformed distributed data is sufficiently small, applying an iteration of steps b through f to the inverse Fourier-transformed updated data.

34. The method of claim 33, wherein the distributing includes scaling the grid by a scaling factor, thereby producing a scaled grid of distributed k-space data.

35. The method of claim 34, wherein the scaling factor is a positive integer power of 2.

36. The method of claim 34, including increasing the scaling factor after step g, prior to the iteration of steps b through f.

37. A method of reconstructing a magnetic resonance image from non- rectilinearly-sampled k-space data, comprising:

) a. distributing the sampled k-space data on a rectilinear k-space grid;

b. convolving the distributed data with a sine function;

c. at least partially based on a characteristic of the sine function, replacing a portion of the convolved data with a corresponding portion of the k-space data distributed on the rectilinear k-space grid, thereby producing a grid of updated data; and

d. until a difference between the updated data and the distributed data is sufficiently small, applying an iteration of steps b through c to the updated data.

38. The method of claim 37, wherein the distributing includes scaling the grid by a scaling factor, thereby producing a scaled grid of distributed k-space data.

39. The method of claim 38, wherein the scaling factor is a positive integer power of 2.

40. The method of claim 38, including increasing the scaling factor after step g, prior to the iteration of steps b through f.

41. A method of reconstructing a magnetic resonance image from non- rectilinearly-sampled k-space data, comprising:

a. distributing the sampled k-space data on a rectilinear k-space grid; b. partitioning the k-space grid into blocks;

c. inverse Fourier transforming distributed data of at least one of the blocks;

d. setting to zero a selected portion of the inverse-transformed data in the at least one of the blocks;

e. Fourier transforming the zeroed and remaining portions of the inverse Fourier-transformed data;

f. at grid points associated with the selected portion, replacing the Fourier-transformed data with the distributed k-space data at corresponding points of the rectilinear k-space grid, thereby producing a grid of updated block data;

g. inverse Fourier transforming the updated block data; and

h. until a difference between the inverse Fourier-transformed updated block data and corresponding inverse Fourier-transformed distributed block data is sufficiently small, applying an iteration of steps b through g to the inverse Fourier-transformed updated block data.

42. A method of reconstructing a magnetic resonance image from non- rectilinearly-sampled k-space data, comprising:

a. distributing the sampled k-space data on a rectilinear k-space grid;

b. partitioning the k-space grid into blocks;

c. convolving distributed data of at least one of the blocks with a sine function;

d. at least partially based on a characteristic of the sine function, replacing a portion of the convolved block data with a corresponding portion of the k-space data distributed on the rectilinear k-space grid, thereby producing a grid of updated block data; and

e. until a difference between the updated block data and corresponding distributed block data is sufficiently small, applying an iteration of steps b through d to the updated block data.

43. A method of adaptively adjusting at least one MR imaging parameter for an interventional procedure, comprising: a. adaptively tracking an MR-guided probe inserted into an object by i. locating the probe by acquiring first probe coordinates with respect to a reference coordinate system; and ii. calculating a velocity of the probe relative to the object by acquiring second probe coordinates with respect to the reference coordinate system; and b. based at least partially on the calculated velocity, adjusting a subset of the at least one MR imaging parameter to adaptively perform at least one of tracking of the probe and imaging of a target region of the object.

44. The method of claim 43, wherein acquiring the first probe coordinates includes acquiring a plurality of one-dimensional frequency- encoded projections to determine at least one of a three-dimensional position of the probe and an orientation of the probe.

45. The method of claim 43, wherein the at least one imaging parameter is selected from the group consisting of: field of view, image spatial resolution, image scan plane position, scan plane orientation, temporal resolution, bandwidth, slice thickness, imaging pulse sequence, image contrast, TR, TE, active receiver channels, k-space trajectory, excitation flip angle, MR scanner table position (e.g., to keep the probe proximal to an isocenter of the table), and a combination thereof.

46. The method of claim 43, wherein the adjusting is at least partially based on an auxiliary parameter.

47. The method of claim 46, wherein the parameter includes an element selected from the group consisting of: position of the probe relative to a target region of the object, position of the probe relative to the MR imager's receive coils, probe orientation, a physiological parameter associated with the object, and a combination thereof.

48. The method of claim 43, wherein adjusting the at least one parameter is at least partially based on a step function of the probe velocity.

49. The method of claim 48, wherein the step function includes a binary function of the probe velocity.

50. The method of claim 43, wherein adjusting the at least one parameter is at least partially based on a smoothly-varying function of the probe velocity.

51. The method of claim 50, wherein the at least one imaging parameter varies between an upper asymptotic value and a lower asymptotic value.

52. The method of claim 51 , wherein at least one of the upper asymptotic value and the lower asymptotic value is at least partially defined based on a hardware constraint of an MRI machine.

53. The method of claim 43, wherein the probe includes a stent.

54. The method of claim 43, wherein the probe includes a catheter.

55. The method of claim 43, wherein the probe includes a tuned resonant circuit having a resonance frequency substantially the same as a resonance frequency of tissue surrounding the probe.

56. The method of claim 43, wherein the object includes an anatomic tissue.

57. The method of claim 56, including using the probe to treat at least a portion of the anatomic tissue as part of the interventional procedure.

58. The method of claim 43, including providing a user interface to a user so the user can adjust at least one of: the at least one imaging parameter, the lower asymptotic value, the upper asymptotic value, the binary function, the smoothly-varying function, and a combination thereof.

59. The method of claim 58, wherein the user interface includes a graphical display for showing the user at least one of the probe and the object.

60. A method of adaptively adjusting at least one MR imaging parameter for an interventional procedure, comprising: a. locating an MR-guided probe inserted into an object by acquiring first probe coordinates with respect to a reference coordinate system; and . b. based at least partially on a parameter associated with the located probe, adjusting a subset of the at least one MR imaging parameter to adaptively perform at least one of tracking of the probe and imaging of a target region of the object.

61. The method of claim 60, wherein the adjusting is at least partially based on an auxiliary parameter.

62. The method of claim 61 , wherein the auxiliary parameter includes an element selected from the group consisting of: position of the probe relative to a target region of the object, position of the probe relative to the MR imager's receive coils, probe orientation, a physiological parameter associated with the object, and a combination thereof.

63. The method of claim 60, including acquiring second probe coordinates with respect to the reference coordinate system to determine a velocity of the probe.

64. The method of claim 62, wherein the adjusting is at least partially based on the determined velocity of the probe.