Homework+3+(not+due)+2D+FFT

=This homework will use Fourier analysis to characterize multiscale structure.=

Let's use 2D scalar data (Black and white images, in other words) and do 2D Fourier transforms.
We can take any image (here .jpg format but others are easy) and play with its spectrum, to see how these spectral slope considerations correspond to the degree of detail/ smoothness in a visual space. [|spectrum_adjust_jpeg.pro] =spectrum_adjust_jpeg.mov = Here are the original are photos and spectra that the above came from (and code) [|spectralslope_from_jpg.pro] an IDL program with one input: a filename for a .jpg file [|spectralslopes_of_art.mov] shows a bunch of examples (full images at bottom of page)

Background on spectral slope in mesoscale/ multiscale meteorology:
1. Theory: If energy is put into a fluid at one scale (wavelength), and passed CONSERVATIVELY to NEIGHBORING (local) scales of motion, in an equilibrium situation in which dissipation at some other scale equals the energy input, then the energy spectrum will have a slope of -5/3 with scale. This -5/3 slope follows from purely dimensional arguments, recognized by Kolmogorov in the 1940s. (Vallis 2008 textbook has good discussion of this). It doesn't depend on the mechanism of scale interaction, although the classic paradigm is vortices advecting vorticity (rotational turbulence). But is it applicable to convection, whose net condensation pumps energy into large scales directly (through the positive definite heating)? Is this log-log slope estimation exercise just too floppy to be very discriminating? Let's look at some art as a baseline, where scales can be whatever somebody drew.



Observations: The Nastrom and Gage (1983) aircraft wind spectra (in the upper troposphere) are the starting point for many discussions. Here for example is a 2008 paper which reproduces the famous figure, shown below. []

Background: about 2D Fourier transforms:
2D fft examples: gathered as .mov files (frames ever 3 seconds, or you can just manually step among them): [|Filter_demo_AndyAnka.mov] shows filtering for one image. [|Filter_demo_jpgs.mov] shows spectra of a variety of art images (see below .jpg links for color originals).

__**THIS HOMEWORK**__: I want you to make the "total wavenumber" spectrum (Power Spectral Density) of one of these (or your own) 2D arrays (jpg images), and try to estimate the spectral slope using a log-log plot.

1. Make an ordinary spectrum plot: PSD = d(variance)/dm plotted vs. m, OR  d(variance)/d(log m) = m d(variance)/dm plotted vs. log(m). In both of these forms of PSD plot displays, area under the curve is proportional to variance.

2. The main point: estimate a power (brightness variance, in these images) spectral slope with scale. Make a log-log spectrum plot: log-power vs. log-wavenumber (or log-wavelength).

The challenge here is to average all the power in each m (total wavenumber) bin. m = sqrt(k^2 + l^2), so {k=3, l=6} and {k=6, l=3} have the same total wavenumber (they represent the same scale). So you need to make m bins and sum up all the power in each bin. This can be done with a call to histogram, but it is a power-weighted histogram: the simple histogram simply counts how many points in the 2D spectrum fall within a given m bin. Histogram with reverse_indices in IDL is what you need. Matlab must have something similar. Examples exist in http://mpo581-hw2.wikispaces.com/. If somebody can find the exact right lines of code, why not paste them right here. No need to repeat efforts. I will eventually if nobody else does.

__Matlab weighted histogram recipe:__ Sil has a solution using find to fill the histogram bins, one at a time. []

__IDL weighted histogram recipe:__ Just substitute normal histogram for hist_2d in this example code. The reverse_indices syntax is identical. http://mpo581-hw2.wikispaces.com/IDL+reverse_indices+trick+for+2D+histograms [|spectralslopes_of_art.mov] [|spectralslope_from_jpg.pro] [|GATE_GigaLES_cloudpic_spectrum.png]

Its Fourier transform is fhat(k,l)
x wavenumber is called k y wavenumber is called l (careful it looks like a 1) total wavenumber is m = sqrt(k^2 + l^2)

Matlab:
fft, ifft, fftshift is the 1D Fourer transform. fft2, ifft2, fftshift2 is for 2D transforms. Here is a code that reads the data, does the fft, applies filters (masking out regions of the spectrum), and transforms back to show the filtered image. Also here is the associated dataset (a photograph, bitmap format, .bmp). You can read any image format just as easily. [|fft_2D_playground.m] [|AndyAnka.BW.bmp]

IDL:
fft by default does both dimensions. For 1D transforms, you give it a dimension keyword. For an inverse transform, you give it the /inverse keyword. Shifting the spectrum so zero frequency is in the middle for display must be done by hand (or with the /shift keyword, available in IDL8 only). [|fft_2D_playground.pro]

Related homework including examples from MPO 581 homework 4:[| MPO581-HW4-2011]. Changheng Chen was lead homework "TA".

Data:
Let's work with images of a few hundred pixels (no sense having it finer that you can even display on a ~1000 pixel computer screen) //This permits us to see 2-3 decades of scale.//

Files:
[|GATE_GIGALES_180x180km2.jpg]
 * ATMOSPHERIC DATA**: Usually science works with wind speed arrays, so that the data value squared (variance or power, the quantity displayed in a Fourier spectrum) corresponds to turbulence-related -3 and -5/3 spectral theories of Kinetic Energy (wind variance). But time is short. Let's just work with a brightness (cloud) image from a cloud model. (Marat Khairoutdinov's "Giga-LES" simulation: 100m grid spacing, 180 x 180 km size). See [|http://rossby.msrc.sunysb.edu/~marat/SAM.html] for more details, animations, etc.), Let's just broadly characterize its spectrum of scales.


 * IMAGES:** These will be our "null" cases where we learn the properties of the Fourier transform. Art - how does the power spectrum relate to its visual impact? These are color, make B&W using the formula in the Matlab code above. Or one could analyze it channel by channel (color images are NX x NY x 3 arrays).

[|IMG_2024_abstract_fewlargescales.jpg] [|IMG_2051_blocks.jpg] [|IMG_2057_circles_fewsizes.jpg] [|IMG_2059_skateboards.jpg] [|IMG_2022_butterflymandala.jpg] [|IMG_2033_oilpainting_multiscale.jpg] [|IMG_2047_comic_withfinelines.jpg] [|IMG_2043_comic_withbigsquare.jpg] [|IMG_2064_rainsplashes_lightgradient.jpg] [|IMG_2063_fewlargescales_landscape.jpg] [|IMG_2061_fractal_leafveins.jpg]

GATE_GigaLES_cloudpic_spectrum.png: