Fouriertransformaties: van warmte tot beeldbewerking

De analyse van golffuncties en golfvergelijkingen kreeg historisch een concrete vorm in de 18^e eeuw. Grote namen uit die tijd, waaronder d’Alembert, Bernoulli, Euler, en Lagrange, bestudeerden het probleem van een vibrerende snaar (zoals een gitaarsnaar). De ééndimensionale golfvergelijking (voor de wiskundeliefhebbers: zie formule 1) werd uiteindelijk in 1746 door d’Alembert ontdekt, en de driedimensionale tegenhanger een decennium later door Euler.

\(\frac{\partial^2 u}{\partial t^2} = c^2 \frac{\partial^2 u}{\partial x^2}\)

Formule 1. De ééndimensionale golfvergelijking. De één dimensionale golfvergelijking zegt dat de ‘kromming’ van een snaar in de loop van de tijd (links van het =-teken) proportioneel is aan de kromming van die snaar in de ruimte (rechts van het =-teken).

De oplossingen voor de golfvergelijkingen – de beschrijving van hóé een snaar kan trillen, dus – waren echter lang niet eenvoudig te vinden. Enkele specifieke functies die aan de vergelijking voldeden waren wel gevonden, maar een algemene analyse van deze vergelijking was nog niet beschikbaar. Dit brengt ons naar de uitvinding die een Franse revolutionair zou vereeuwigen: fourieranalyse. Een wiskundige techniek, geïntroduceerd door Jean-Baptiste Joseph Fourier.

Fourier

Joseph Fourier werd geboren in het Franse Auxerre, in 1768. Hij was van gewone komaf, maar was vanaf zijn negende levensjaar een weeskind. Desondanks deed hij het goed op school, en specifiek in wiskunde blonk hij uit. Op zijn 17^e meldde hij zich aan bij het Ministerie van Oorlog, om bij de genie te mogen. Aangezien Fourier niet van adel was, kon hij echter niet toetreden tot deze elite. In plaats daarvan trad hij toe tot de kerk in Auxerre waar hij les ging geven.

De uitbraak van de Franse revolutie die niet veel later zou volgen, veranderde Fouriers leven dramatisch. Omdat hij een aanhanger was van de revolutie, werd hij enkele malen gevangen gezet. Hij overleefde zijn gevangenschappen en verhuisde naar Parijs, waar hij aan de École Normale Supérieure les zou gaan geven. Niet veel later, in 1795, wisselde hij naar de École Polytechnique waar hij de beroemde wiskundige Joseph-Louis Lagrange zou opvolgen. Interessant genoeg ging Fourier vlak na zijn aanstelling bij de École Polytechnique mee met een expeditie van niemand minder dan Napoleon Bonaparte, als wetenschappelijk adviseur. Bij deze expeditie werd onder andere door een Franse militair de bekende Rosetta-steen gevonden.

Na zijn terugkeer in Frankrijk werd Fourier door Napoleon gevraagd voor de positie van gouverneur van het department Isère in Grenoble. Hij wordt beschreven als een zeer competente ambtenaar, en overzag onder andere de constructie van wegen tussen Frankrijk en Italië. In diezelfde periode leverde hij een aanzienlijke bijdrage aan het epische werk “Description de l’Égypte” (De beschrijving van Egypte), aan de hand van zijn ervaringen tijdens Napoleons expeditie. Dit werk was toentertijd het grootste verslag over de Egyptische cultuur vanaf de prehistorie tot de 19^e eeuw.

Het was tijdens zijn aanstelling in Isère, rond zijn 35^e levensjaar, dat Fourier zijn werk aan zijn magnum opus begon: “Théorie analytique de la chaleur” (Analytische theorie van warmte).

De theorie van warmte

In zijn boek beschreef Joseph Fourier enkele belangrijke bevindingen. Zo vond hij de vergelijkingen voor de voortplanting van hitte, introduceerde hij nieuwe wiskundige technieken om deze vergelijkingen op te lossen en heeft hij die oplossingen kunnen toepassen om uiteindelijk zijn resultaten te vergelijken met metingen. Desondanks werd Fouriers werk aanvankelijk niet goed ontvangen. De beoordelaars van zijn manuscript waren niet tevreden over zijn methodes en twijfelden aan de universaliteit en rigorositeit van zijn voorgestelde theorieën. In 1822 werd het werk uiteindelijk toch gepubliceerd, en werd de wereld bekend met fourieranalyses.

Laten we beginnen met de warmtevergelijking zoals die werd geïntroduceerd door Fourier. Voor het gemak beschouwen we de ééndimensionale versie van de warmtevergelijking, maar de beschrijving kan vrij eenvoudig veralgemeniseerd worden naar meerdere dimensies. De warmtevergelijking beschrijft hoe de temperatuur in een bepaald medium, bijvoorbeeld een gas of een vloeistof, verandert als functie van de tijd en plaats

Laten we voor de duidelijkheid de volgende notatie introduceren: u(t,x). De functie u(t,x) beschrijft de temperatuur van het medium op een bepaald tijdstip en op een bepaalde positie. De snelheid waarmee deze temperatuur verandert is gegeven door de tijdsafgeleide van deze functie: \(\frac{\partial u(t,x)}{\partial t}\). Vervolgens maakte Fourier de volgende belangrijke observatie: de snelheid waarmee de temperatuur verandert op een positie (x,t), moet proportioneel zijn aan de verandering van de hoeveelheid warmte die toe- of wegstroomt van de positie (x,t). De verandering van de warmtestroom in de x-richting (de enige richting in ons ééndimensionale systeem) is de verandering van de zogeheten temperatuurgradiënt. Dit klinkt allemaal erg complex, maar je zou het kunnen vergelijken met water in een leiding: de snelheid waarmee water door een pijp stroomt wordt veroorzaakt door de verandering in de waterdruk op elke positie, en moet daardoor daarmee evenredig zijn. Hoe meer de druk verandert op een positie x, hoe sneller het water gaat stromen (de verandering in tijd). Voor warmte werkt dat alles precies hetzelfde. Als formule ziet het bovenstaande er als volgt uit:

\(\frac{\partial u(x,t)}{\partial t} = K \frac{\partial^2 u(x,t)}{\partial x^2}\)

Formule 3. De warmtevergelijking.Het symbool u geeft de temperatuur aan, die afhangt van zowel de plaats als de tijd. De formule zegt hoe de temperatuur in de loop van de tijd verandert (de linkerkant) als de afhankelijkheid van de plaats op een gegeven moment bekend is (de rechterkant).

Hier staat precies wat we hierboven al in woorden zeiden. De linkerkant van de vergelijking in formule 3 geeft de verandering van de temperatuur in de loop van de tijd weer. De rechterkant heeft een constante K, die afhangt van het materiaal waar warmte doorheen gaat, en verder een dubbele afgeleide van de temperatuur u naar de positie. Dit komt overeen met de toestroom van warmte op dat punt. (De eerste afgeleide is de temperatuurgradiënt.) De tweede afgeleide naar de plaats geeft dus de verandering van de temperatuurstroom weer.

Het is deze vergelijking die opgelost dient te worden als we willen begrijpen hoe warmte ‘stroomt’. Om dit te doen maken we gebruik van een eigenschap van de warmtevergelijking die ‘lineariteit’ heet: als we meerdere oplossingen kennen, dan is de som van deze oplossingen ook een oplossing. Zonder verder in te gaan op details: we kunnen daarmee uiteindelijk de formule vinden die ons vertelt hoe bijvoorbeeld een staaf van een meter lengte afkoelt als functie van tijd. Om dat te doen hoeven we slechts twee ‘basisoplossingen’ te vinden; lineaire combinaties (optellen en aftrekken) van deze oplossingen geven ons dan elke andere mogelijke oplossing.

De twee basisoplossingen van de warmtevergelijking worden gegeven door twee sinusvormige functies. Met behulp daarvan kunnen we dus een formule vinden waarmee bijvoorbeeld een lange, dunne staaf (of een ander ééndimensionaal object) afkoelt of opwarmt.

In alle wiskundige glorie ziet die oplossing eruit als hieronder, formule 4. We nemen daarbij voor het gemak aan dat we beginnen met een uniforme temperatuurverdeling over de hele staaf – de staaf is aan het begin bijvoorbeeld overal precies u₀=100 graden warm.

\(u(x,t) = \frac{4 u_0}{\pi} \sum_{n=1}^\infty \frac{\sin((2n-1)\pi x)}{2n-1} e^{-K(2n-1)^2 \pi^2 t}\)

Formule 4. De oplossing van de warmtevergelijking. Exacte wiskundige oplossing voor het afkoelproces van een staaf van een meter lengte (x varieert dus van 0 tot 1). Een belangrijk detail is dat de oplossing gegeven wordt als een gewogen som van sinusoïden. De gewichten worden gegeven door de exponenten helemaal aan het eind van elke term in de oneindige som. De constante K vertelt ons hoe snel de termen afnemen in waarde als functie van de tijd. Een grotere waarde van K komt overeen met een snel afkoelende staaf.

Omdat deze formule er wellicht wat cryptisch uit ziet, visualiseren we hetzelfde afkoelproces met behulp van de grafiek in afbeelding 5. Hierin is de hoogte (ook weergegeven door de kleur) overeenkomstig met de temperatuur. De breedte van de grafiek komt overeen met de positie x, en de diepte geeft de tijd t weer.

Het nieuwe idee dat Fourier introduceerde in zijn werk, en bij het vinden van deze oplossing, is het volgende: elke functie op een eindig interval, dus hier voor waarden van x tussen 0 en 1, blijkt te schrijven te zijn als een som van trigonometrische functies (sinusoïden). In formuleform luidt deze bewering als volgt: elke f(t) kun je met behulp van bepaalde constantes a_n en b_n schrijven als de volgende fourierreeks:

\(f(t) = \sum_{n=-\infty}^\infty a_n \sin \left(2 \pi \frac{n}{T} t \right) + b_n \cos \left(2 \pi \frac{n}{T} t \right)\)

Formule 6. Fourier’s theorema. Elke functie functie f(t) op een eindig interval is te schrijven als een som van sinussen en cosinussen, elk met een bepaalde amplitude (de a_n en de b_n), en een fase (de \(2\pi \frac{n}{T}\)).

Fourier’s theorema is ook terug te zien in formule 4. De oplossing is ook daar simpelweg een som van sinussen, met een amplitude a_n gegeven door de exponenten. In dit voorbeeld komen geen cosinussen voor, dus elke b_n is gelijk aan 0.

Hoe vind je een fourierreeks?

Het doen van fourieranalyses (het vinden van de coefficienten a_n en b_n) blijkt geweldig nuttig te zijn, zelfs voor onze alledaagse wereld. Voor we overgaan naar enkele real-life-toepassingen van Fouriers werk, geven we eerst twee voorbeelden van fourierreeksen – reeksen zoals in formule 6 – en beschrijven we hoe die reeksen bepaalde functies kunnen benaderen. Zoals gezegd: indien een functie aan de juiste voorwaarden voldoet (bijvoorbeeld dat die op een eindig interval is gedefinieerd, maar men kan hetzelfde doen voor periodieke functies), dan kan die geschreven worden als een oneindige som van sinusoïden. Wat we moeten zien te achterhalen zijn de coëfficiënten a_n en b_n die staan voor elke (co)sinus zoals in formule 6.

Laten we een concreet voorbeeld bekijken om te laten zien hoe dit proces in zijn werk gaat. Stel dat we een functie hebben, en die kunnen schrijven als de som van twee sinussen met een periode  van twee- en éénmaal π respectievelijk:

\(f(t) = a_1 \sin(t) + a_2 \sin (2t)\)

Wat we dus moeten zien te achterhalen zijn de waardes van de coëfficiënten a₁ en a₂. Nu is er een mooi trucje om deze coëfficiënten te achterhalen. Dit kan namelijk gedaan worden door gebruik te maken van het feit dat de lange tijdsgemiddeldes van producten van twee sinusoïden met verschillende frequenties altijd nul zijn. Alleen bij twee gelijke frequenties is het gemiddelde van zo’n product niet gelijk aan nul maar aan 2 pi. Dat klinkt ingewikkeld, maar is eenvoudig te visualiseren aan de hand van ons voorbeeld. We vermenigvuldigen beide kanten van de formule hierboven met sin(t). Vervolgens integreren we het resultaat over meerdere periodes (bijvoorbeeld over een tijdsinterval 4π) en delen door het aantal periodes (middelen!) In formulevorm ziet dat er zo uit:

\(\sin(t) f(t) = a_1 \sin^2(t) + a_2 \sin(t) \sin(2t)\)

\(\frac{1}{4\pi} \int_0^{4\pi} \sin(t) f(t) = \frac{a_1}{4\pi} \int_0^{4\pi} \sin^2(t) + \frac{a_2}{4\pi} \int_0^{4\pi} \sin(t) \sin(2t)\)

Wat we hierboven zeiden, is dat we nu verwachten dat de tweede term nul zal zijn. Wat in dat geval overblijft aan de rechterkant is de term met a₁, en die constante kunnen we vervolgens berekenen met de integraal aan de linkerkant. De reden waarom de tweede term met de verschillende frequenties nul oplevert kun je zien in afbeelding 6. Zoals te zien in de afbeelding heffen de oppervlaktes onder en boven de x-as elkaar precies op, wat betekent dat het gemiddelde gelijk is aan nul.

Met behulp van het proces dat hierboven beschreven is, is het relatief ‘eenvoudig’ om de coëfficiënten a_n en b_n te bepalen voor elke sinusoïde in de fourierreeks. Door deze coëfficiënten vervolgens in te vullen krijgen we het gewenste resultaat. Ter illustratie is dit proces grafisch weergeven in afbeelding 7 voor twee typische periodieke functies. Het aantal termen wat we meenemen in de reeks moet feitelijk oneindig groot zijn om de functie f(t) perfect te benaderen, maar dit is voor een computer natuurlijk niet praktisch. Daarom gebruiken we simpelweg een groot maar eindig aantal termen (sinusoïden) uit de fourierreeks. Het aantal termen is in de figuren aangegeven boven de afbeelding. Zoals te zien wordt de benadering steeds beter, naarmate we meer termen toevoegen aan de reeks.

Het nut van fourieranalyses: geluidmodulatie

Laten we nu naar enkele praktijkvoorbeelden kijken waar fourieranalyses (het beschrijven in termen van sinusoïden) een bijzonder nuttige rol spelen. Het allereerste voorbeeld is direct een van de bekendste toepassingen: geluid. Specifieker, het bewerken van geluid. Stel dat je bijvoorbeeld een geluid hebt opgenomen, maar om de een of andere reden zit er in de opname een hoge toon die het geluid vervelend maakt. De microfoon kan natuurlijk alleen de totale drukgolf van het geluid meten, niet de individuele componenten uit de bijbehorende fourierreeks. Op de een of andere manier moeten we de hoge toon toch zien weg te filteren uit het totale signaal, om het geluid ‘prettig’ te maken. Fourieranalyse geeft hierbij de uitkomst.

Wat we moeten doen is uiteindelijk vrij eenvoudig: we identificeren alle frequenties die bijdragen aan het geluid (schrijven de functie als een fourierreeks!). Vervolgens kunnen we kunstmatig de hoge ‘pieptoon’ wegfilteren, simpelweg door de bijdrage van de hoge frequentie te verwijderen – de bijbehorende coëfficiënt a_n zetten we dus in de computer op 0. Vervolgens doen we een ‘inverse’ fouriertransformatie waarbij we alle individuele componenten – behalve natuurlijk de component die we net weggefilterd hebben – weer bij elkaar optellen. De uitkomst is een ‘prettige toon’ zonder de hoge piep.

Het nut van fourieranalyses: Lord Kelvins getijdevoorspeller

Een ander illustratief voorbeeld van de kracht van fourieranalyses is te vinden in de techniek voor het voorspellen van de getijden. Zoals waarschijnlijk bekend worden eb en vloed ‘aangestuurd’ door gravitationele effecten van de zon en de maan. Deze effecten zijn elke dag net iets anders. Deze variaties komen voor het grootste deel door maar liefst vijf verschillende effecten: de rotatie van de aarde om zijn as en om de zon, de rotatie van de maan, en de precessie van de maan, zowel in zijn rotatievlak als de precessie van dat rotatievlak zelf. De precieze details zijn op zichzelf al bijzonder interessant, maar niet voor dit artikel van belang. Wat wel van belang is, is het feit dat we, aan de hand van deze astronomische effecten, het getijde (symbolisch weer weergegeven met u(t)) kunnen beschrijven als een som van ongeveer 120 sinusoïden. De frequenties van al deze sinusoïden zijn bekend, en worden bepaald door de duur van de hiervoor genoemde astronomische rotaties. In formulevorm ziet dat er weer als volgt uit:

\(u(t) = a_1 \sin(\nu_1 t) + b_1 \cos(\mu_1 t) + a_2 \sin(\nu_2 t) + \ldots\)

met bekende waarden voor \(\nu_n\) en \(\mu_n\). Het probleem is dus weer om de amplitudes a_n en b_n te vinden. Deze amplitudes zijn, in tegenstelling tot de periodes, echter niet eenvoudig uit eigenschappen van de natuur te achterhalen. Ze hangen namelijk af van de diepte van de zeebodem en andere lokale parameters die sterk afhankelijk zijn van de plek op aarde waar het getij voorspeld dient te worden.

Nu blijkt het dat het getij zich vrijwel exact elke 19 jaar herhaalt. Dat is namelijk de periode waarin alle astronomische fenomenen die het getij beïnvloeden zich samen gaan herhalen. Stel dus dat we gedurende 19 jaar metingen hebben gedaan op een bepaalde locatie, bijvoorbeeld in een haven. Dan kunnen we de coëfficiënten bepalen in de fourierreeks met behulp van de technieken die eerder in dit stuk zijn uitgelegd (formule 6). Vanaf dit moment is het voldoende om simpelweg de formule in te vullen om het getij te bepalen voor elk tijdstip voor de komende duizenden jaren! Dat is een stuk makkelijker en accurater dan duiken in oude archieven en oude metingen opzoeken, en maakt het bovendien mogelijk te interpoleren tussen gedane meetpunten.

De berekeningen die nodig zijn of de fouriercoëfficiënten te bepalen zijn redelijk complex, en deden mensen uiteraard niet met de hand. Dat is veel te veel werk om te doen voor elke dag voor elke haven. Om het proces te vereenvoudigen is in de 19^e eeuw een machine ontwikkeld die de fourieranalyse voor je doet. De uitvinder van deze machine was niemand minder dan Lord Kelvin, naar wie de gelijknamige temperatuurschaal is genoemd. In afbeelding 9 is een illustratie van zijn machine te zien.

De grootte van de verschillende tandwielen komt overeen met de frequenties van het getijde – die zijn gelijk voor elke plek op aarde. De lengtes van de touwen komen overeen met de coëfficiënten a_n en b_n, en die kunnen ingesteld worden voor elke plek op aarde. Als dit gedaan is, kan men simpelweg aan de hendel draaien, en zal een potlood aan het einde van de machine een grafiek tekenen die het getij voorspelt voor de plek waarop de machine is ingesteld!

De machines die Kelvin ontwierp werden gebruikt tot ver in de 20^e eeuw. Pas toen digitale computers krachtig (en goedkoop) genoeg werden om dit proces in te programmeren met software, werden deze machines uitgefaseerd.

Uiteraard is het niet zo dat er voor elke plek op aarde 19 jaar aan meetdata beschikbaar zijn. Een kortere meetperiode kan in principe ook gebruikt worden. De grootste effecten van het getij worden immers bepaald door de fenomenen met de hogere frequenties (de rotatie van de aarde om zijn as, elke 24 uur, en van de maan om de aarde, elke maand). Een maand aan metingen kan dus al gebruikt worden om heel behoorlijke voorspellingen te doen aan de hand van de resulterende coëfficiënten!

Het nut van fourieranalyses: afbeeldingsbewerking (1)

Als afsluiter van deze tekst heb ik mezelf een beetje laten gaan, en het Quantum Universe-logo onder handen genomen. Het is namelijk niet alleen mogelijk om een ‘lijn’ (dus: een functie van één variabele zoals de tijd) te fouriertransformeren. We kunnen ook een functie op een 2-dimensionaal oppervlak fouriertransformeren – een afbeelding, bijvoorbeeld! We maken hiervoor gebruik van een 2-dimensionale variant op de fouriertransformatie om de coëfficiënten te vinden van de 2-dimensionale golven die samen de afbeelding vormen. De figuren in afbeelding 10 geven een weergave van het resultaat. De linker twee figuren geven de informatie weer die komt uit de fouriertransformatie (een coefficient, en een frequentie). De rechter figuur is het resultaat van een inverse fouriertransformatie, waarbij alle informatie uit de eerste twee figuren gebruikt wordt om de tweedimensionale golven op te bouwen. Wanneer het resultaat van de ‘ruis’ in de linker twee figuren volgens Fouriers recept bij elkaar wordt opgeteld, resulteert dat dus in de oorspronkelijke afbeelding, het logo.

Het nut van fourieranalyses: afbeeldingsbewerking (2 – de leukste!)

Het vorige voorbeeld is nog vrij droog en misschien weinig illustratief. Daarom proberen we iets anders! Als we het logo van de website (afbeelding 10) namelijk goed bestuderen, zien we dat het bestaat uit twee ingekleurde losse randen. Deze randen zijn één dimensionale paden, die we, omdat ze periodiek (gesloten) zijn, volgens Fourier zouden kunnen beschrijven met sinussen en cosinussen. We willen dus de vorm van de binnenste rand en de buitenste rand schrijven als een som van sinusoïden. Hiervoor moeten we een paar stappen doorlopen. Eerst moeten we de beweging van de randen ontbinden in een beweging in de x-richting en y-richting. Het resultaat hiervan voor de binnenrand is weergegeven in afbeelding 11. Merk op dat de begin- en eindpunten van beide lijnen in deze grafiek op dezelfde hoogte liggen. Dit is uiteraard noodzakelijk voor een periodieke baan: aan het einde ben je weer terug waar je begon.

We kunnen de twee lijnen uit de rechter grafiek in afbeelding 11 nu fouriertransformeren. De fouriergetransformeerde van de blauwe lijn vertelt ons hoe we de x-coördinaat kunnen beschrijven met sinusoïden. Die van de oranje lijn vertelt ons hetzelfde over de y-coördinaat. We kunnen al deze informatie combineren om een ‘functie’ te maken, die ons vertelt hoe de lijn zal lopen in het vlak als een periodiek pad. Hoe meer termen uit de fourierreeks we meenemen in ons uiteindelijke functievoorschrift voor de paden, hoe accurater het resultaat zal zijn. In afbeelding 12 is het resultaat te zien voor verschillende aantallen termen in de fourierreeks. De lijn die we zien bewegen is het resultaat van het fouriertransformeren van het functievoorschrift voor de x– en y– coördinaten uit afbeelding 11. De gekleurde achtergrond is alleen voor de show gekozen, om de afbeelding er nóg mooier uit te laten zien.

Daarmee hebben we de nodige voorbeelden van het gebruik van fouriertransformaties gezien. Er valt echter nog veel meer over dit onderwerp te vertellen – bijvoorbeeld over de rol die fouriertransformaties spelen in de quantummechanica. Daarover veel meer in deel 2 van dit artikel!