Image Processing & Antialiasing

Image Processing & Antialiasing Part II (Aliasing, Sampling, Convolution, and Filtering)

Outline • Overview • Example Applications • Jaggies & Aliasing • Sampling & Duals • Convolution • Filtering • Scaling • Reconstruction • Scaling, continued • Implementation 9/19/2013

Jaggies & Aliasing • “Jaggies” an informal name for artifacts from poorly representing continuous geometry by a discrete 2D grid of pixels – a form of aliasing (later…) • Jaggies are a manifestation of sampling error and loss of information (aliasing of high frequency components by low frequency ones) • Effect of jaggies can be reduced by anti-aliasing, which smoothes out the pixels around the jaggies by averaging • Shades of gray instead of sharp black/white transitions • Diminishes HVS’ response to sharp transitions (Mach banding) 9/19/2013

Representing lines: Point sampling, single pixel • Midpoint algorithm: in each column, pick the pixel with the closest center to the line • A form of point sampling: sample the line at each of the integer X values • Pick a single pixel to represent the line’s intensity, full on or full off • Doubling resolution in x and y only lessens the problem, but costs 4 times memory, bandwidth, and scan conversion time! • Note: This works for -1 < slope < 1, use rows instead of columns for the other case or there will be gaps in the line Line approximation using point sampling Approximating same line at 2x the resolution 9/19/2013

Representing lines: Area sampling • Represent the line as a unit width rectangle, use multiple pixels overlapping the rectangle (for now we think of pixels as squares) • Instead of full on/off, calculate each pixel intensity proportional to the area covered by the unit rectangle • A form of unweighted area sampling – stay tuned: • Only pixels covered by primitive can contribute • Distance of pixel center to line doesn’t matter • Typically have more than one pixel per column so can go gradually from dark for pixels covered by the line to white background; the more area of overlap, the darker the pixel 9/19/2013

Filtering as form of averaging: filter is a weighting function • Apply a unweighted or weighted average to each sample and it’s neighbors • Example: filtering for student height or weight – give each sample equal weight or diminish influence of samples further away • Weight function gives weight for incremental area centered at • The volume of a filter function needs to be 1 to preserve overall brightness, so W =1 over its “support”, the area it covers, 0 elsewhere Pre-Filtered samples Weight function Post-filtered samples (decreases contrast between adjacent pixels): red = added value, light blue is subtracted area 9/19/2013

“Box Filter” Represents Unweighted Area Sampling • Box filter is constant over all area and is a single pixel wide here, but could vary in width • For each pixel intersecting the line, intensity contributed by each sub-area of intersection is • Then total intensity of the pixel (between 0 and 1) integrated over area of overlap is: • Integral is volume over area of overlap (in above figure, a rectangular wedge) 9/19/2013

“Cone Filter” for Weighted Area Sampling • Area sampling, but the overlap between filter and primitive is weighted so that sub-areas with closer to center of pixel count more • Cone has: • Linear falloff • Circular symmetry • Width of 2 (called support) • Intensity of pixel is the “subvolume” inside the cone over the line (see picture) 9/19/2013

2-unit circular support of filter Weighted Area Sampling (Continued) Pixel center (+) W Area of overlap between support and primitive Differential area dA1 Differential area dA2 Primitive 9/19/2013 W is the weight which is multiplied with at ; normalize Wto make volume under cone = 1

Another Look at Point Sampling – Even Box Filter is Better! • This simplistic scan conversion algorithm only asks if a mathematical point is inside the primitive or not • Bad for sub-pixel detail which is very common in high-quality rendering where there may be many more micro-polygons than pixels! A C D B Point-sampling problems. Samples are shown as black dots. Object A and C are sampled, but corresponding objects B and D are not. 9/19/2013

Another Look at Unweighted Area Sampling: Box filter • Box filter • Support: 1 pixel • Sets intensity proportional to area of overlap • Creates “winking” of adjacent pixels as a small triangle translates (b) Unweighted area sampling. (a) All sub-areas in the pixel are weighted equally. (b) Changes in computed intensities as an object moves between pixels. 9/19/2013

Another Look at Weighted Area Sampling: Pyramid filter • Pyramid filter • Support: 1 pixel • Approximates circular cone to emphasize area of overlap close to center of pixel (b) Weighted area sampling. (a) sub-areas in the pixel are weighted differently as a function of distance to the center of the pixel. (b) Changes in computed intensities as an object moves between pixels. 9/19/2013

Another Look at Weighted Area Sampling: Cone filter • Cone filter • Support: 2 pixels • Greater smoothness in the changes of intensity (b) Weighted area sampling with overlap. (a) Typical weighting function. (b) Changes in computed intensities as an object moves between pixels. 9/19/2013

Pseudocode and Results Related phenomenon: Moire patterns for each sample point p: //p need not be integer! place filter centered over p for each pixel q under filter: weight = filter value over q p.intensity += weight * q.intensity “Crawlie Demo” 9/19/2013 Aliased Anti-aliased

Anti-Aliasing Example Close-up of original, aliased render Antialiasing Techniques: Blur filter – weighted average of neighboring pixels Supersampling - sample multiple points within a given pixel and average the result Checkerboard with Supersampling Supersampling and Blurring 9/19/2013

Aliasing and Anti-Aliasing: Intro to Signal Processing • Easier to discuss aliasing and anti-aliasing of images and image scaling than of rendering anti-aliased geometric primitives - return at end of this unit • Need to do a brief intro to signal processing, since an image can be thought of as a sequence of scan lines, each displaying continuous intensity signals (photograph on fine-grain film) or discrete signals (pixels captured by a digital camera) • Signals are composed of combinations of waves (acoustical waves varying in time for sound, light waves varying intensity in 2D space for images) • Signals are complex waveforms but they can be decomposed into combinations of simpler waveforms (sines and cosines) of varying frequency, amplitude and phase: Fourier synthesis 9/19/2013

Sampling and Scaling of Images(images courtesy of George Wolburg, Columbia University, in 2nd Edition Computer Graphics: P&P) • Scan converting an image is digitizing (sampling) a series of continuous intensity functions, one per scan line. • We will use single scan lines for simplicity, but everything still applies to full images Scan line from synthetic scene Scan line from natural scene 9/19/2013

The Sampling/Reconstruction/Display Pipeline – overview Original continuous signal: Sampled signal (lose information): Reconstructed signal (approximation to original): (many reconstruction methods) One reconstruction method (linear interpolation): 9/19/2013

Fourier Waveform Synthesis – Approximate a Continuous Signal Approximation of scan line signal from image improves with addition of more sine waves • Any signal can be approximated by summing sine (and cosine) waves of different frequencies, phases, and amplitudes • A signal has 2 representations. We’re familiar with the spatial domain, but every signal also can be represented in the frequencydomain 9/19/2013

Digression: Music • A vibrating string creates a pressure wave with multiple frequency components • Intensity/volume of sound is related to amplitudes of components of pressure wave • Lowest frequency is “fundamental frequency” Labels show 1/f -- the wavelength • Integer multiples of the fundamental frequency are called “harmonics” or “overtones” (non-integral multiples also exist) • Power of two multiples of the fundamental frequency are called “octaves” For more, applets to interact with these concepts: http://homepages.gac.edu/~huber/fourier/index.html http://www.phy.ntnu.edu.tw/ntnujava/index.php?topic=17 http://www.chem.uoa.gr/Applets/AppletFourier/Appl_Fourier2.html 9/19/2013

Analogous Operations, Duals • Easier to visualize signals and filtering/averaging in the frequency domain • Many problems are easier solved by transforming to another problem (an untransforming!) • Characterize a signal by its frequency spectrum showing amplitudes of component waveforms • Signal and frequency domains are duals and represent identical information. • Easier to filter in frequency domain than in signal domain, as we shall see. • Take a familiar problem: multiplication of numbers • Can take logarithm of number, perform operations on log(number), then take antilogarithm) • Dual of multiplication is addition in logarithmic “space” If ab=cthen log(a) + log(b) = log(c) • This “invertible transformation” makes slide-rules such effective tools for multiplication: manipulating sliders corresponds to manipulating numbers via their logs 9/19/2013

Frequency Spectrum of a Signal Signal Domain Frequency Domain • Sine wave is characterized by amplitude and frequency • Frequency of a sine wave is number of cycles per second for audio, or number of cycles per unit length (e.g., inter-pixel distance) for image scan lines • Can characterize any waveform by enumerating amplitude and frequency of all its component sine waves (Fourier transform – see chapter 18 in book) • This can be plotted as a “frequency spectrum”, a.k.a. power spectrum, (we ignore negative frequencies, but they are needed for mathematical correctness) To see spatial and frequency domains of simple signals: http://www.cs.brown.edu/exploratories/freeSoftware/repository/edu/brown/cs/exploratories/applets/fft1DApp/1d_fast_fourier_transform_guide.html 9/19/2013

Sine waves in 2D • Let white be max positive amplitude, black max negative amplitude Higher Frequency (Harmonic – integer multiple) Low Frequency http://cns-alumni.bu.edu/~slehar/fourier/fourier.html 9/19/2013 Images: Hays

Simple Fourier Synthesis example (1/2) • example : g(t) = sin(2πf t) + (1/3)sin(2π(3f) t) + = 9/19/2013 Slides: Hays

Simple Fourier Synthesis example (2/2) • Approximate signal using sines and cosines Spatial domain Original Function Partial approximation with 2 cosines Partial approximation with 3 cosines Partial approximation with 4 cosines Approximation with 5 cosines f(x) We are just summing waves to get bolded blue function to approximate the red signal x 9/19/2013

Sampling: The Nyquist Limit • In order to capture all frequencies in a signal, we must sample at a rate that is higher than 2 times the highest frequency (with a non-zero value from the Fourier transformation) in the signal (the Nyquist limit). • Here is an approximate analog sine wave: • The sine wave sampled at an acceptable rate (4 times the highest frequency): • Reconstructed wave based on these samples: 9/19/2013

Aliasing: Know Thine Enemy • Aliasing occurs when we sample a signal at less than twice maximum frequency • Here is our analog sine wave again: • Here is the sine wave sampled at too low a rate: • Here is the reconstructed wave based on these samples: • The reconstruction isn’t even close! 9/19/2013

Temporal Aliasing: Another Sampling Error • Ever seen tires spin in a movie? Have you ever noticed that sometimes, they seem to be spinning backwards? • It’s because the video frame-rate is lower than twice the frequency at which the wheels spin. This is temporal aliasing! • You see this a lot in movies because the effect is so striking.It’s known as the stage-coacheffect. • Porsche Dyno Test 9/19/2013

Sampling at the Nyquist Limit • Sampling right at the Nyquist limit can also be problematic • Here is our perfect analog sine wave again: • Here is the sine wave sampled at the Nyquist limit. This time it works fine: • Here is the same sine wave sampled at the Nyquist limit, with the sample points shifted. Now we get no signal: • For an applet dedicated to illustrating the Nyquist Limit: • http://www.cs.brown.edu/exploratories/freeSoftware/repository/edu/brown/cs/exploratories/applets/nyquist/nyquist_limit_guide.html 9/19/2013

The Enemy is Recognized • Aliasing is shown in bottom diagrams on previous slides – Signals that are sampled at too low a rate can reconstruct high frequencies as low frequencies • These low frequencies are “aliases” of high frequencies • The low sampling rate data could not adequately represent the high frequency components, so it represented them incorrectly, as low frequencies • So, we just sample above the Nyquist limit, right? • Regrettably, we can’t always do that • What about this square wave signal? just a test pattern of alternating black/white bars • Let’s try using Fourier Synthesis Number of sine waves: 9/19/2013 125 25 5

Infinite Frequencies • Square waves have infinite frequencies at the jumps, so how can we sample correctly at those edges? • We can’t. Pure and simple. , and we can’t sample at an infinite rate. • Unfortunately, infinite frequencies are the rule in synthetic computer graphics - discrete transitions between adjacent pixels • So, do inverse operation. Instead of increasing our sampling frequency to meet that of the signal, we: • Pre-filter out high frequencies we can’t show • The signal is now guaranteed to consist only of frequencies we can represent and reconstruct with reasonable accuracy • This isn’t the same signal that came in, but it’s better than an aliased version • Reconstructing the pre-filtered approximate signal beats reconstructing sampled original signal, with its corrupting aliases for the highest frequencies 9/19/2013

Infinite Frequencies • The more high frequencies we pre-filter out, the lower the sampling frequency needed but the less the filtered signal resembles original. • Note: pre-filtering is often just abbreviated as filtering, but the prefix “pre” helps remind us that post-filtering (i.e., another stage of filtering after image computation or transformation) is also practiced. If post-filtering is done on reconstruction of samples of original signal, it will blur in the aliases present in the corrupted reconstruction! 9/19/2013

Scale Aliasing, or “Why do we have to pre-filter?” This doesn’t look right at all. There are no stripes and the image now has a blacker average intensity Image with sample points marked Image scaled using point samples Original Image 9/19/2013

Scale Aliasing II, or “Close, but no cigar?” Better, but not perfect… Prefiltered image with samples marked Prefiltered image scaled Original Image 9/19/2013

Scale Aliasing III, or “Why is it still wrong?” • Good: The pre-filtered scaled image now has same overall brightness • Bad: still no stripes • Filter removed “high frequencies” from image • discontinuities that were stripes • Given number of points to represent image once scaled, not enough points to represent high frequencies; even higher sampling rate doesn’t help… • Nyquist limit says that we can’t represent frequencies higher than ½ our sampling rate, thus can’t do better than this blurred approximation 9/19/2013

Filtering as an Application of the Convolution Integral • Remember filtering is just (weighted) averaging of samples in a neighborhood • How is it done in practice? Just a nested for loop through all the sample points • either sample at integer pixel locations or in between (e.g., for image up-scaling) • based on the theory of signal convolution – the convolution integral • we’ll use a finite approximation to that integral • First we’ll low-pass pre-filter to remove (as much as possible) un-representable frequencies to minimize aliases (doable in most non-pathological cases for images, not for rendering, where we’ll use “supersampling” and post-filtering – lecture IP-IV) • Second, there will be another form of filtering when the samples are “reconstructed” as a continuous signal approximating the original signal • in the case of raster graphics, reconstruction combines • the h/w generation of pixels from pixel values as an image that feeds • the perceptual apparatus of the HVS (Human Visual System) which typically smoothes out discontinuities (spatially and temporally), a kind of integration 9/19/2013

Low-Pass Filteringto Eliminate HighFrequencies (shown for one scan line in Spatial Domain) 9/19/2013

Ideal Low-Pass Filtering (Frequency Domain) • Multiplying by the box function in the frequency domain • Frequencies under the box are kept, and the high frequencies are cut off • Corresponds to convolution with the sinc function in the spatial domain … … 9/19/2013

Discrete convolution – Function as an Array • This is the type of convolution we will be using once we begin coding • If we evenly sample a function, get an array of point samples, e.g., (R,G, B) triples as color intensity values • Some mathematicians think of an array of numbers as a function. Think this way for the next few slides! • The larger the array (the more samples/unit interval) the more accurate the representation of the function and the better our ability to anti-alias and improve the reconstruction [ 1.2, 0, -1.2, 0, 1.2, 0, -1.2 ] 9/19/2013

Discrete Convolution for Filtering • Have two arrays (functions), f and g, and we assign g as our “filter”. For the sake of convenience, force the filter to have an odd number of points so that there is a center • Once we are done with our convolution, we will produce a new array (function) • You’ll do this for blur and edge detection, as well as for image scaling • Take the filter, and line up the center point of the filter with each element, i, of f. Multiply each pair of elements and sum the products. Assign this value to element i of the output array • When we are scaling, this procedure will change slightly because where we center the filter and where our value goes will be different! • One hiccup: The output array will be slightly larger than our initial array because we start our multiplication/sum process at the edges, where it will generate a non-zero result • It is easier to see! 9/19/2013

Discrete Convolution – Visually (1/2) [ ⅓, ⅓, ⅓] [ 1.2, 0, -1.2, 0, 1.2, 0, -1.2 ] Function, f Function, g, a box filter 9/19/2013

Output Discrete Convolution – Visually (2/2) [ 1.2, 0, -1.2, 0, 1.2, 0, -1.2 ] .4 [ .4 _ _ _ _ _ _ _ _ ] [ ⅓, ⅓, ⅓ ] [ 1.2, 0, -1.2, 0, 1.2, 0, -1.2 ] [ .4 .4 _ _ _ _ _ _ _ ] [ ⅓, ⅓, ⅓] • [ .3, , .3] [ 1.2, 0, -1.2, 0, 1.2, 0, -1.2 ] [ .4 .4 0 _ _ _ _ _ _ ] [ ⅓, ⅓, ⅓ ] 0.0 [ 1.2, 0, -1.2, 0, 1.2, 0, -1.2 ] [ .4 .4 0 -.4 _ _ _ _ _ ] [ ⅓, ⅓, ⅓ ] -.4 9/19/2013

Discrete Convolution -- Review • Think of an array as a function • We take two arrays and generate a third • We “slide” the filter along the other array and at each element, calculate a value by multiplying the pairs and summing the products to do the (weighted) averaging • Our output array will be slightly larger than our input array if we start too far over the edges. We can (and will) ignore this for image processing and just take “inner” part of the array • Another hiccup to be aware of: Our process has an additional step if our filter is not symmetric. We will mention this more in the coming slides. Don’t worry too much about it though, in this course most filters will be symmetric 9/19/2013

What is (Continuous) Convolution? • Convolution is a process that takes two function, f and g, and produces a new function h. • When you create a function you need to define it's value for each input, i.e. h(x) = ? • For convolution, h(x) equals the integral of the product ofsignal f and filter g, centered at x. • This should sound very similar to our discrete notion of convolution (integrals are like sums) • You can think of convolution as sliding g over f to produce h because as we shift the center of g, we evaluate h at a new point, x • in theory, an infinite number of x’s; in practice for discrete convolution, a finite number of x’s. We’ll have to figure out where those should be in the original image, and that will differ for up-scaling and down-scaling, as a function of scale factor 9/19/2013

Convolution Mathematically (1/2) • Convolving signal with filter function • Changing variable of integration to be 𝜏, ranging from 0 to ±∞starting at x Center filter g at x Multiply the pairs Exactly what we were doing, just in a formula! Integrate (sum) the products 9/19/2013

Convolution Mathematically (2/2) • Convolving signal f(x) with filter function g(x): • At each point x, h(x) is integral of product f of and g, except g is flipped and translated so its origin is at x • In practice, f or g is often an even function (symmetric about the y-axis) and we don’t need to flip • If g(x) has finite support, it does a weighted average centered at x • f(x) (blue signal) convolved by g(x) (red filter) to get result h(x) (black signal) • f(x) and g(x) are box functions • Each point on black signal is result of an integral, represented by the area under the product of f(x) and g(x) (yellow area) • Note: in filtered signal h(x), discontinuities of f(x)are smoothed out and the base is broadened http://www.cs.brown.edu/exploratories/freeSoftware/repository/edu/brown/cs/exploratories/applets/convolution/convolution_guide.html 9/19/2013

Simple Convolution Example Convolution is a lot like multiplication Try the applet: http://www.cs.brown.edu/exploratories/freeSoftware/repository/edu/brown/cs/exploratories/applets/twoBoxConvolution/two_box_convolution_guide.html 1111 ×1111 1111 1111 1111 +1111 1234321 * 9/19/2013

Properties of Convolution • Commutative • Associative • Distributive • Identity • is the Dirac delta– a mathematical artifice • everywhere except at • Area under is • Can be thought of as the limit as a Gaussian narrows 9/19/2013

Image Processing & Antialiasing