Transcription
Hello everybody, and welcome back. This is part two in a three-part series looking at image reconstruction in CT. In the first part, we focused on simple back projection. Now, we're going to move on to a separate entity known as filtered back projection. And you'll notice there's a lot of overlap in the methodology here, but we're going to make some subtle changes that are going to account for some of the shortcomings that were noticed in simple back projection.
So, let's revisit the image that we created in simple back projection. Remember, we've got projection data at multiple different angles around the patient, and we're smearing that data back, that attenuation data back, across these different angles. And ultimately, that data is going to accumulate and form this image here. Now, you can see the major issue with this image is the blurring here. We've got poor edge definition, and we've got a lot of blurring, radial blurring, occurring from all the attenuating structures here. That's the main downfall of simple back projection. We've got good contrast. What we've got here is what's known as a low-pass filter. We're going to cover that later on. We want to now get rid of this blurring.
Remember, when we're creating this image in simple back projection, this is a mathematical calculation. We're overlapping data in a specific array, and then we're allocating each pixel within this image a specific value, an attenuation value, based on all these different projections here. This is a calculation, it's not an actual image. And the process of simple back projection is introducing blur into the image in a predictable and calculable way. What we want to do then is account for that blur with some mathematical function that's going to allow us to extract the true attenuation data that's underlying here. It's this process here of applying some mathematical function to this simple back projection data that is known as filtered back projection, and that's what we're going to focus on today.
Now, I think it's useful to revisit the formation of that simple back projection image and see how we're introducing blur into the image. And we're going to look at what a low-pass filter is. Notice that as we add each projection here, we're smearing attenuation, higher attenuation values, beyond the subject that's attenuating our X-rays. There's no attenuation occurring in these regions here, but the data suggests that there is, and that's because of that smearing of the data. As we add each layer, you can see how we lose edge definition and add blur to the edges here. The edges aren't well-defined. What also happens as we add each layer of projection data, the attenuating structures are getting brighter and brighter. We're adding more and more values to these attenuating structures in the matrix that we're forming here.
Now, importantly, when we're creating an image like this, we're again creating an array of numbers, an array of attenuation values that represent our calculated attenuation at each pixel in the image. If we were to look at areas in the image where there's very little change over distance to signal intensity, that's what's known as low-frequency parts of the image. You can see if we were crossing this part of the image here, there's very little change in true attenuation data. As we reach the edge here, where we get a large change in signal intensity, that's what's known as our high-frequency areas in our image. We're going to look at later, we're going to Fourier transform spatial images into frequency domain images, and we'll look closer at low and high-frequency data.
What a low-pass filter essentially is, it's allowing low-frequency data to pass through the filter whilst attenuating high-frequency data. We are reinforcing areas in the image that are low-frequency, where signal change doesn't occur very much, and we are suppressing the high-frequency data areas, causing blur at the edges of our image. Now, obviously, we're creating a very simple image here. You can imagine if this was a clinical image, how much more complex these structures would be, and the edges in that image are going to be very poorly defined. We're going to create an image that actually isn't usable to us clinically.
Now, we can think of this image that we've created as containing two separate data sets. It's got our true attenuation data, and it's got some blur function that we've introduced by the process of simple back projection. Now, if we were to perform some mathematical operation between these two data sets, we would create a simple back projection image, that blurred image that we saw. This mathematical operation is what's known as convolution. It's going to be outside the scope of this talk. You think about multiplying two data sets, or adding two data sets, or subtracting two data sets. Another way that we can combine two data sets is by this process of convolution. You might see this sign here, you might also see this. This also represents convolution. I'm going to leave some resources below this video if you're unsure as to what convolution is. Convolution is often used in image processing, and you'll see it when you are trying to blur images, or when you're trying to sharpen images. We use convolution a lot. I'd highly encourage you to go check out those resources if you want to dive more into this. The explanation that we're using today is more of a visual explanation, not a mathematical explanation.
But if we convolute this data, you can see the blur function here. If we look at it as a 3D function, how much higher the magnitude is at the center of this blur function, and then we get radial expansion of the blur function out here. This is showing us a symmetrical nature of blur in simple back projection. It's also what's known as a point spread function. This blur function is going to be applied to every pixel within this data set here. It's not only applied to the middle and spread out. So, we combine these two data sets, and we create our simple back projection image. This process is just defining how simple back projection goes about creating an image.
I want to show you, if we were to have an image with just two points on here, notice how the blur is radial from each pixel within this image here. We're not getting blur from the center out, where we get more blur on the outside of the circle. This is what's known as a point spread function. Every point in our calculated image has blur emanating out from each one of these points. That's a function of this blur function here. So, we take our simple back projection image. Now, we need a way of extracting that blur function out and calculating the true attenuation data. That's filtered back projection. Let's look at that a little bit more closely.
Now, in order to remove this blur function using filtered back projection, we first need to understand how blur is represented in an image. And up until now, we've been looking at images in the 2D spatial domain. We've got X and Y coordinates that show us what the attenuation value is at specific locations in this 2D image. And we said that blur is introduced into this 2D image by convoluting the blur function here. Now, we've also seen projection data that's one-dimensional data. We're looking at the attenuation data through one angle, through our patient. So, we can get 2D spatial domain, and then we can represent data in the one-dimensional spatial domain. This is for a specific angle through the patient, and we represent this X-axis now as R. It's known as the radial distance away from the center of the detector. Again, this is a one-dimensional spatial domain. We say that the attenuation data at a specific location is represented by this graph here.
Now, if we were to use this one-dimensional spatial data to back project and create a simple back projection image, the blur that would be introduced would be represented by this function here, a one-dimensional spatial domain blur function. Now, the formulas here aren't important. What's important to note is in the spatial domain, we convolute the blur function into the spatial domain that describes how blur is introduced into the image. What's also important to note is notice how this blur function is discrete data. When we measure data, we take an analog signal and we convert it into a digital signal. It's that digital signal that we're going to use to smear over in simple back projection. Therefore, our blur function corresponds to that discrete data.
Now, why am I going in a lot of depth here in the spatial domain? It's because we can also represent images in a separate domain known as the frequency domain. Now, when we're using the spatial domain, importantly, the blur function is convoluted. Now, why am I going over this? Well, convolution is computationally very demanding. It's not multiplying or adding two data sets. It's a very computationally demanding process because we have to take the blur function, rotate it a full 180 degrees, and then apply it to each and every location on the spatial domain image. It's a very demanding process.
Now, the frequency domain is different. When we introduce blur in the frequency domain, we introduce it through a multiplication, not a convolution. When we look at our image here, we say there are high-frequency and low-frequency regions within the spatial domain. Areas where signal doesn't change very much over a certain distance are what's known as low-frequency regions in the image, and areas where there's an abrupt change in attenuation values, like here on the edges of our image, is what's known as high-frequency regions in our image. So, high-frequency data corresponds to the edges, low-frequency data corresponds to regions where signal doesn't change very much over a given distance.
Now, whenever we have a graph like this one here, that's going over a set distance with different magnitudes, we can perform a function on it known as a Fourier transform. That's going to show us what frequency of sine waves need to be combined together to create this specific function here. And if we were to perform a one-dimensional Fourier transform on this projection here, we'd get something that looks like this. The center of this graph represents the low-frequency data. As we head out to the peripheries, we get the sine waves or the frequencies that represent the edges, the high-frequency data points in our image. This one-dimensional Fourier transform, we use our axis known as K here, which represents the different spatial frequencies that occur within this original one-dimensional spatial domain. You can see here that the lower frequencies predominate, and we get very little higher frequencies out at the periphery. We've got very few edges within this image here as well. The center of this image is what's known as the DC component, that represents the average intensity for our entire projection here.
Now, the same can be done in the two-dimensional spatial domain. We can perform what's known as a two-dimensional Fourier transform, and we create an image that looks like this. This again is a data set, and it's essentially these one-dimensional Fourier transforms placed at different angles on a specific matrix. The center of this image is what's known as the low-frequency region of the image. As we head out to the peripheries, we're talking about the high-frequency data from this two-dimensional image. The specifics here aren't hugely important, just knowing where the low-frequency and where the high-frequency data is. The intensity, the brightness of the signal here, represents the magnitude. We've got a lot of low-frequency data, very little high-frequency data. It turns out that for different projection angles, we can create different one-dimensional frequency domains for each projection angle, and essentially those projection angles correspond to the different angles running through this two-dimensional Fourier transform. They're not, we can't place them directly on here, but they do correspond. That's what's known as the central slice theorem, which again is outside the scope of this talk.
Now, just for completeness' sake, let's label the axes here. We've got U and V. V represents the spatial frequency changes in the vertical portion of our original image, and U represents the spatial frequency changes along the horizontal range of our image here. We can see we get combinations of U and V as we change that angle through our two-dimensional image. Now, why go through the bother of creating these frequency domains? Surely we're adding an extra step into the process. Well, the key here comes into figuring out how is blur represented in the frequency domain. And I said to you that in the spatial domain, we represent blur in the image by convoluting that data. In the frequency domain, we multiply blur into the image. This represents the low-pass filter that I was talking about. Notice how the frequencies in the center of this, the low frequencies, are amplified by our blur function, and the high frequencies are suppressed by this blur function. We can multiply this function directly with this one-dimensional frequency domain data set. The same occurs in the two-dimensional frequency domain. We can see here's the blur function. The center of the image, all of these values get really amplified. We get a low-pass filter. Low-frequency data is allowed to pass through the filter, high-frequency data is suppressed, and this is the crux of filtered back projection.
What we want to do is reverse this blur function. So, we've created different images in the spatial domain and the frequency domain. They all represent the data that used to create a simple back projection image. Now, let's look at how we go about removing that blur. And I'm going to start with the frequency domain because this is what's most commonly done in practice. We're removing blur in the frequency domain here. Here's our simple back projection image, and here's our two-dimensional and one-dimensional Fourier transform of these images. Now, we said that blur is a multiplication of the blur function and this one-D Fourier transform, as well as this blur function with a 2D Fourier transform. These are representing the same thing, but we're dealing with a two-dimensional frequency domain and a one-dimensional frequency domain. If we were to perform an inverse Fourier transform, we would create a projection that looks like this. We've seen this before. Again, this is not going to work. If we then back project at multiple different angles, we're going to create the blur because of this blur function. The same in a two-dimensional Fourier transform.
Now, this is where we go about removing blur from our image. We apply what's known as a ramp filter. This is the filter in filtered back projection. A ramp filter suppresses this low-frequency data that's been increased by our blur function, and we are increasing the high-frequency data in the blur function. We're counteracting this blur function. We want to bring about more high-frequency data and suppress some of that low-frequency data, and this ramp filter is what allows us to do that. We can see in the two-dimensional domain, the same thing happens. A ramp filter amplifies the edges of our blur function here, allows us to get more high-frequency data. Now, this blur function is innate within our one-dimensional Fourier transform, so we need to apply this ramp filter to the one-dimensional Fourier transform, and that's going to negate that blur function.
Now, if we were to perform an inverse Fourier transform, take frequency data and convert it into spatial data, we would get a data set that looks like this. Much, much different. Notice the change here. We've got these negative lobes at the edges of our image here, and we've got these positive lobes on the other side of each edge. The same mathematical function can be performed in the two-dimensional Fourier transform. We apply that ramp filter, we only multiply it. Very easy for a computer to multiply two arrays together instead of having to convolute data. Again, we amplify the peripheries here. The brighter signal here is amplifying the higher frequency and suppressing the lower frequency data. We take this data that has been filtered using a ramp filter, and we can back project our one-dimensional inverse Fourier transform, or we can perform a two-dimensional inverse Fourier transform here, and essentially what's happened is we've removed the blur function and created a filtered back projection image. It's a lot to go through, but you can see where the filter comes from. This ramp filter is what's counteracting the blur function that has been introduced through the process of spreading data at multiple different angles through an image.
Now, we can also remove blur in the spatial domain. It's only ever really done in the one-dimensional spatial domain. If it's ever done, it's so much easier to remove blur in the frequency domain because of that multiplication with the ramp filter, not needing to deconvolve. So, let's get our projection through here. This is our 1D projection. We know that blur is introduced through the convolution of this function and our projection data. In order to counteract this blur that we're going to introduce through the process of simple back projection, we use what's known as a deconvolution kernel. You may have heard of kernels before, you may have heard of deconvolution. Within filtered back projection, when we're talking in the frequency domain, we're adding a ramp filter. We're allowing the higher frequencies to pass through that ramp filter, and we're suppressing those lower frequencies. When we're talking in the spatial domain, we no longer really talk about a ramp filter, we talk about deconvolution. This is just convoluting this kernel, this data set, with our projection here and creating the same projection, filtered projection, that we created in our one-dimensional filtered back projection.
Now, the reason that this filtered projection is the same in the spatial domain as the one we created in the one-dimensional frequency domain is because that deconvolution kernel is essentially an inverse Fourier transform of the ramp filter that we applied in the frequency domain. Again, we can use this, back project all of these at different angles, these filtered projections, and create a filtered back projection image. Notice the better quality here.
Now, I know these processes that we've gone through are a little bit abstract. We haven't actually defined any mathematical functions, we've just been giving broad brushstrokes. What I want to do is give you a visual representation of how this is actually occurring. And we're going to take a very simple image here, just a circle. This is the only attenuating structure within this image. Now, if we were to pass X-rays at this angle through this image, we're going to create this projection data here. And I want to put that projection data to the side.
Now, let's go about the steps of filtered back projection. The first thing we need to do is perform a Fourier transform of this data. A low-frequency data at the middle, higher-frequency data is out to the periphery. We know that to remove blur from this frequency data, we need to apply a ramp function. Now, you'll notice I've been curving this ramp function out like this. I'm going to show you why that's important later on. Once we've applied a ramp function to this frequency data, we can then perform an inverse Fourier transform, and our initial projection data has now been converted to something that looks like this. Again, those negative lobes at the edges of our image and positive lobes on the other side of each edge.
Now, this data is a little bit complicated. What I want to do is draw over it with a solid orange line. Again, I want to give you a visual representation here. This data, I'm going to put back onto our image here. You notice that on the peripheries, the filter has given zero values. Then we get negative lobes, we've got negative values here, and positive values where the attenuation is going to occur. Again, negative and zero values. We now just apply the process of simple back projection. We smear that data over the angle that the X-rays came from. Positive values, and these are going to be represented by a light color here. Negative values are going to be represented by a yellow color on our image. Notice how as we project that data over, now we've created negative lobes on either side of this projection that we've averaged out here. The zeroed part of the image remains zero, remains black. There's zero attenuation occurring here.
Now, let's change the angle again. Our projection data is going to be the same, no matter why we go around this image, because we've created such a simple circle in the middle here. When negative lobes cross over previously positive values, those values are going to cancel each other out, and that's the beauty of this function here. It's very specific. The values that we have smeared over this region here are going to be canceled out by these negative values here. Watch how that happens here as I simple back project that function over a different angle. We get negative and positive lobes canceling each other out. This is the crux of filtered back projection. That's why we need this filtered back projected function here. If we change the angle again, notice how projecting that over, we again get cancellation. These negative lobes cancel out these positive lobes here, and you can imagine these being placed much tighter together. Remember, we're doing a thousand, 2,000 projections around a patient for a single slice. Notice now I'm not going to run it over each time, I'm just going to place these now at the different angles. Notice, let me go back. Notice how the black is increasing here, we're canceling out with these negative lobes, and we're reinforcing the center. The center is getting brighter and brighter. How edges are also getting more and more smooth.
If we were to continue this process for the 1,000 projections, you would see how that would ultimately end in an image that looks like this. It's this process of applying the ramp filter, either in the one-dimensional or the two-dimensional frequency domain, that's going to give us these crisp edges here. We're essentially amplifying the higher-frequency data and decreasing the lower-frequency data. We've got an image now with great contrast and great edge definition. Now, remember, a patient doesn't look like this, but the process stays exactly the same. We're going to get better edge definition with filtered back projection over simple back projection.
Now, this is all well and good, but why then don't we often use filtered back projection anymore? Why are we using the next talk's topic, iterative reconstruction? Well, the main issue comes when we're looking at noise. We know inherently there's going to be noise in our image, and later on, when we look at image quality, I'm going to show you how noise is introduced in the image. But say we were to create a sinogram here, but we were to add noise into the sinogram, either through the detectors or through patient factors here. We can see that there's background noise. How then does filtered back projection deal with this background noise? Well, it doesn't deal with it very well.
Let's take one line from our sinogram and plot the predicted projection that we would have through this patient. This is what we would predict. But when there's noise in the image, I'm going to show you what noise looks like with this green data overlying our predicted data. This is what our detectors actually register in the computer in the CT scanner. If we were to perform the Fourier transform, we would see our Fourier transform, our low-frequency and high-frequency data. But what's really important is to zoom in on this point. See the green line here? That's the noise contribution to the frequency data. If we zoom in on the peripheries, look how the noise is actually higher than some of that high-frequency signal that was creating our edge definition.
If we apply the ramp filter, we will increase these higher frequencies. We're going to increase noise, and we're going to increase the true higher frequencies from our image. This is the reason that the ramp filter curves up like this in an attempt to suppress some of that noise that's being added into the image. However, even with the curving of this ramp filter, not just perfectly straight lines like you might see in some textbooks, we're still going to amplify noise here in the one-dimensional frequency domain and create an image that looks like this. We've got some contrast, we can see where the attenuating and non-attenuating structures are, and we've got edge definition. We have amplified the higher-frequency data points here, so we've got good edge definition, but we've created a terrible image with lots of noise. It's going to inhibit our ability to create any sort of meaningful image.
Now, to get around this, what we can do is apply a low-pass filter, or get rid of the higher-frequency data here. I'm just going to draw over with an orange line here. We're getting rid of this, and then we can perform an inverse Fourier transform on this data set and back project that image, and we end up creating an image that looks like this. We're back to square one. We've got an image that's very blurry, got good contrast, but we've lost our edge definition. We've had to get rid of those higher frequencies. And this is the major issue with filtered back projection: is that if there's any noise within the image, we're going to amplify that noise with our ramp filters, and that's going to create an image that's incredibly noisy. If we get rid of that noise, it comes at the expense again of edge definition.
Now, why go through that whole process only to end up with essentially the same type of image that we started with? Well, in the next talk, in part three of this image reconstruction, we're going to look at a process known as iterative reconstruction, and it often starts with this as our initial image, a filtered back projection image, and it uses this data that we created here and iterates on it to go about removing noise, accountable mechanical factors, and ultimately create a much better image than we've created here. So, let's move on to iterative reconstruction.
I know there's a lot to cover in this talk, but I'd encourage you to just know the basic steps: how we go about creating frequency data, apply a ramp filter, an inverse Fourier transform of that frequency data that's had the ramp filter applied, and then create our filtered back projection image. So, I'll see you all in the next talk. Until then, goodbye everybody.