matlab程序代写 MRI 图像处理

In MRI, the spatial localisation of internal tissues are realised by gradient encoding, which innovation was awarded 2013 Nobel Prize. The theory behind it is as simple as Fourier transformation. The MR signal is received by RF coils and stored in k-space (frequency domain), which are later transformed into an image by applying inverse Fourier transform. Namely, the forward Fourier encoding produces k-space data, and the image reconstruction is realised by using inverse Fourier transform. However, in practice, acquiring full k-space data are time consuming; therefore, fast MR imaging are always preferable. One common method is to use phased-array coils acquiring partial k-space and then reconstruct the image with Sensitivity encoding (SENSE) method or GeneRalized Autocalibrating Partial Parallel Acquisition (GRAPPA).

Task1:

  1. (1)  Please code Fourier transform and inverse Fourier transform with MATLAB in two methods. You can use “for loops” as one of your method, and the other method of your own selection, such as “repmat” (exclude MATLAB build-in functions “fft”, “fft2”). Apply your Fourier transform on the given brain image to generate k-space data in the form of Ax = b. The brain image is 512×512, resize it to a 64×64 image to satisfy your computer memory when performing Fourier transform.
  2. (2)  Reconstruct the produced k-space data into the brain image with your coded inverse Fourier transform. Compare the computational time of your methods (two methods in Task 1-1) with the MATALB given function ‘ifft2’. Use the “fft2” generated k-space data as the reference and generate k-space error maps of your method to prove your method is correct. Is the Fourier transform matrix invertible? What is its inverse matrix?

Task2: From task2, you can use MATLAB functions fft2 and ifft2. In task 2, resize the image to 256×256.

(1). Generate and display individual coil images with provided sensitivity profiles (from 8- element phased-array coils) and then produce a combined image (256×256) with using the sum-

2

coili

(2). Generate an undersampled k-space with the provided undersampling pattern (reduction rate = 4, evenly undersampling), and then compare the reconstructed image with the original image. Why do we see aliasing artefacts when applying inverse Fourier transform? Produce an undersampling pattern that only sample central 64 lines. Compared the reconstructed image from this undersampling pattern with the image reconstructed from evenly undersampled k- space, why two images are differed so much even the same amount of data are sampled?

Task3: In GRAPPA, the central 32 k-space lines are normally fully sampled, which are referred as to auto-calibration (AC) lines. The reconstruction of missing sample in the i-th coil at an unacquired position r is simply given by :

x(r)(PR y)H g (1) i rr ri

Where gri are sets of reconstruction weights to be calculated. Rr represents the operation of choosing a block of k-space (from all the coils) out of the entire grid around the k-space positions indexed by r. The operators Pr represent local sampling patterns that choose only acquired samples from a block of k-space. Let y be a multi-coil k-space grid concatenated into a vector in which unacquired data are zero filled. So, the product PrRry is a vector containing only the acquired k-space neighborhood around the k-space position r.

Thegri canbeobtainedbysolvingtherelationinEq.(1)atdifferentpositionsink-space,where the xi are known. Typically, this is done in a fully acquired region in the centre of k-space, e.g., autocalibration (AC) region. To perform the calibration, it is useful to construct a so-called calibration matrix, denoted by A, from the AC portion of the acquired data. It is constructed by sliding a window throughout the AC data, taking each block (Rry)H inside the AC region to be a row in the matrix. The columns of A are shifted versions of the AC area, leading to a matrix structure known as Block-Hankel.

To obtain conditions for the weights gri , Eq. (1) is rewritten using the calibration matrix and applied to all locations inside the AC region. This yields a set of ideal conditions for the reconstruction weights:

yAC APHg (2) i rri

where y AC are data from the i-th coil inside the AC region. By construction, one of the columns i

of-square (SOS) method. Hint: element-wise production; Imsos 

 i1

(Im

)

Ae  yAC
of A is yAC . We can write this as i i , where, ei is a vector with “1” in the appropriate

i
position that chooses the i-th coil data, and “0” elsewhere.

Please try to:
(1). With the given data, generate a k-space undersampling pattern where the reduction rate is 4 but with central 32 lines fully sampled, then form the calibration matrix A.

8

(2). Explain the relationship between PH g e and calibration matrix A through SVD r ri i

analysis.