How to run a Fourier transform and PSD on a runway profile and read it at take-off speeds

A power spectral density (PSD) of a runway profile shows how much of the surface roughness sits at each wavelength, and a speed transform re-labels that wavelength axis as the frequency a wheel meets at each ground speed. This guide runs both steps on a synthetic profile, with every computed number reproduced by a script, and states what a spectrum cannot do: say where a bump is, or what peak load it causes (NLR, 2009) (SciPy, 2026).

Key facts

• A one-sided displacement PSD of a runway profile is reported in m³, which is m² per cycle/m, and its integral over a band is the variance in that band.
• The temporal frequency a wheel meets is f = v/λ, so a 30 m wavelength gives 0.67 Hz at 20 m/s and 2.50 Hz at 75 m/s in the synthetic example.
• The speed transform S_z(f) = G_d(n)/v assumes constant speed and a stationary profile, so an accelerating take-off roll is only approximated by snapshots at fixed speeds.
• NLR-TP-2009-190 notes that a PSD does not distinguish many small bumps from a few large ones at a given wavelength.
• In the synthetic case, two inserted bumps raised the fitted PSD level by about 30 % and moved single bins up or down (x0.95 and x0.82), but nothing in the spectrum says how many bumps there are or where; a rolling chord found both.

What are a Fourier transform and a PSD of a runway profile?

A runway profile is the surface elevation along a line, sampled at a fixed interval. A Fourier transform splits it into sine waves of different wavelengths; a PSD reports how much variance each carries. The spatial frequency n = 1/λ (cycles/m) is waves per metre, so the PSD G_d(n) is in m² per cycle/m, which is m³ (SciPy, 2026).

The variance in a band is the area under the PSD across it. The PSD used here is one-sided, with negative frequencies folded in, the reporting form of the road-profile standard ISO 8608 (Mechanical vibration, Road surface profiles, 2016). That standard's scope names roads, streets, highways and off-road terrain, not runways, and this guide prints none of its class boundaries (ISO, 2016). For the road index, see IRI, rut depth and texture depth and TMH 13.

How does wavelength become frequency as the aircraft accelerates?

FAA AC 150/5380-9 (Guidelines and Procedures for Measuring Airfield Pavement Roughness, 30 September 2009), paragraph 2.1, says the critical bump height depends on the relationship between the wavelength and the natural frequency of the aircraft (FAA, 2009). NASA Memo 2-21-59L (1959), PDF page 3, found gear loads greatest at impressed frequencies of 1.5 to 2 cycles per second. That is a 1959 finding, not a limit (NASA, 1959).

The formula, term by term: wavelength to frequency

Equation 1: f = v / λ = v · n

symbol	what it is in words	unit	where its value comes from	typical range
f	temporal frequency a wheel meets	Hz	computed	0.67 to 12.50 Hz in the table below
v	ground speed	m/s	roll speed at the chainage	20, 50, 75 m/s (about 39, 97, 146 kt), synthetic
λ	wavelength in the profile	m	the PSD axis, 1/n	30 m and 6 m, synthetic
n	spatial frequency	cycles/m	1/λ	0.033 and 0.167 cycles/m

Where it comes from: kinematics, not an empirical fit. A wheel moving at v covers one wave of length λ in a time λ/v. Frequency is the inverse of that period, so f = v/λ, and since n = 1/λ, f = v·n.

Worked numbers: λ = 30 m at v = 50 m/s gives f = 50 ÷ 30 = 1.67 Hz; λ = 6 m at 75 m/s gives 12.50 Hz. The author's calculation for the 1959 band: 2 Hz at 61.7 m/s (120 kt) gives λ = 61.7 ÷ 2 = 30.9 m.

What it means: the same wavelength excites a rising frequency as the aircraft accelerates. It says nothing about the aircraft's natural frequencies, which the aircraft data supplies.

The formula, term by term: where in the roll a wavelength meets a frequency

Equation 2: x = (f₀ λ)² / (2a)

symbol	what it is in words	unit	where its value comes from	typical range
f₀	frequency of interest	Hz	chosen by the analyst	1.5 Hz, the low end of the 1959 NASA band
λ	wavelength in the profile	m	the PSD axis	12, 30, 60 m, synthetic
a	constant acceleration from rest	m/s²	take-off roll model	2.0 m/s², synthetic
x	roll distance where λ meets f₀	m	computed	81, 506, 2025 m in the example

Where it comes from: derived from kinematics. Under constant acceleration from rest, v² = 2ax. Equation 1 gives the speed at which λ excites f₀: v = f₀λ. Substituting, x = (f₀λ)²/(2a).

Worked numbers: with f₀ = 1.5 Hz and a = 2.0 m/s², λ = 30 m needs v = 45.0 m/s, reached at x = (1.5 × 30)² ÷ (2 × 2.0) = 506 m; λ = 12 m gives 18.0 m/s at 81 m; λ = 60 m needs 90.0 m/s at 2025 m, beyond the synthetic 75 m/s top speed of the roll.

What it means and where it stops: it shows which wavelengths matter where in the roll. Real acceleration is not constant and the roll does not start from rest at chainage zero, so use the measured or planned speed against chainage.

Narrated film · 2:48

From runway profile to PSD at take-off speeds

Wavelength to frequency at 20, 50 and 75 m/s on a synthetic profile.

Read the film transcript

Does a runway's roughness spectrum tell you where the bumps are, and what they do to a wheel at take-off speed? A runway profile is the height of the surface along one line, sampled at a fixed step. This synthetic profile is three thousand metres long, sampled every quarter metre. A power spectral density, the PSD, says how much variance each wavelength carries. Equation three: the variance sigma squared in a band is the integral of G sub d of n, from n one to n two. G sub d is the one-sided displacement PSD, in metres cubed, which is metres squared per cycle per metre. The limits are the band you choose. The band zero point zero one to zero point five cycles per metre holds four point zero square millimetres, a root mean square of one point nine nine millimetres. Welch's method recovers the input slope: minus two point zero one, against minus two put in. Equation one turns wavelength into frequency: f equals v divided by lambda. Here f is the frequency in hertz, v is ground speed in metres per second, and lambda is wavelength in metres. A thirty metre wave at fifty metres per second gives one point six seven hertz. Equation two finds where in the roll a wavelength meets a frequency of interest, f zero. With constant acceleration a, x equals f zero times lambda, all squared, divided by two a. Taking one point five hertz, thirty metres and a synthetic two metres per second squared, x is five hundred and six metres. Equation four: the temporal PSD, S sub z of f, equals G sub d of n divided by v, with n equal to f over v. S sub z is in metres squared per hertz. At one point five hertz it is zero point four one, one point zero five and one point six four, times ten to the minus six, at twenty, fifty and seventy-five metres per second. The same wavelengths land at higher frequencies as speed rises: six metres is three point three three hertz at twenty metres per second. These are constant-speed snapshots, not a take-off spectrum. Two bumps were inserted: six metres by twenty millimetres, and thirty metres by twenty-five millimetres. They raised the fitted level by about thirty per cent, yet the spectrum cannot say how many bumps there are, or where. A rolling chord finds them: nineteen point nine four millimetres at four hundred and fifty metres, and twenty-four point zero eight at nine hundred. Load enters through a gear model: a sprung mass, a tyre spring, and a response H from profile to tyre force. In this synthetic two-mass model, the spectral and time-domain root mean square loads agree within about three per cent. That checks the arithmetic, not an aircraft, and a root mean square is not a peak load. A PSD also misses reach: highway laser profilometers catch wavelengths up to forty-five metres, an aircraft up to about one hundred and twenty. Keep the profile, the Welch settings, the speeds chosen, and the script that produced them.

How do you run the PSD, step by step?

1. Fix the profile and its sampling interval Δx. The synthetic example is 3,000 m at Δx = 0.25 m, so the shortest wavelength is 2Δx = 0.50 m (2 cycles/m).
2. Remove the mean and any grade, and state any high-pass cut-off; Welch's constant detrend removes only each segment's mean.
3. Choose the segment length. Here 2,048 samples is 512 m, so the bin spacing is Δn = 1/(512 m) = 0.00195 cycles/m and the longest wavelength is about 512 m. The lowest bins are attenuated by the window and by mean removal, so the level fit uses 0.03 to 0.5 cycles/m (SciPy, 2026).
4. Run Welch's method with a Hann window and 50 % overlap, one-sided, density scaling.
5. Check the level against a known input.
6. Choose the ground speeds, apply Equation 1 and Equation 4, and state the assumptions.

The formula, term by term: one-sided displacement PSD

Equation 3: σ²(n₁, n₂) = ∫ G_d(n) dn from n₁ to n₂

symbol	what it is in words	unit	where its value comes from	typical range
G_d(n)	one-sided displacement PSD	m³ (m² per cycle/m)	Welch estimate of the profile	4.0e-6 m³ at 0.1 cycles/m, synthetic
n₁, n₂	band limits	cycles/m	analyst's choice	0.01 to 0.5 cycles/m
σ²	variance of the profile in the band	m²	integral of G_d	4.0 mm² in the example

Where it comes from: derived. The squared Fourier magnitude, averaged over overlapping Hann-windowed segments, estimates the variance density; folding negative frequencies gives the one-sided form (SciPy, 2026).

Worked numbers: the synthetic input was G_d(n₀) = 4.0e-6 m³ at n₀ = 0.1 cycles/m with slope −2.0; the Welch fit on 0.03 to 0.5 cycles/m recovered slope −2.01 and G_d(n₀) = 4.07e-6 m³. The band 0.01 to 0.5 cycles/m holds 4.0 mm², an RMS of 1.99 mm, while the whole profile's RMS is 14.05 mm in this one random realisation (expected about 11 mm for this input), dominated by wavelengths up to the 3,000 m length that a stated high-pass cut-off removes.

What it means: the RMS depends on the band, so state it. Do not mix this form with a two-sided PSD, which is half the one-sided value, or with a PSD per rad/m, which is the per-cycle/m value divided by 2π.

How does the spectrum change with ground speed?

Variance does not change with speed; the frequency scale stretches.

The formula, term by term: the speed transform

Equation 4: S_z(f) = G_d(n) / v, with n = f / v

symbol	what it is in words	unit	where its value comes from	typical range
S_z(f)	temporal one-sided PSD of the wheel input	m²/Hz	G_d divided by v	0.41e-6 to 1.64e-6 m²/Hz, synthetic
G_d(n)	spatial PSD at n = f/v	m³	Welch bin nearest n	Equation 3 block
v	ground speed	m/s	as above	20, 50, 75 m/s

Where it comes from: derived from conserving variance. Because ∫ G_d dn = ∫ S_z df and df = v dn, S_z = G_d/v. This is standard random-vibration theory, and G_d(n) follows the ISO 8608 one-sided convention (ISO, 2016).

Worked numbers at f = 1.5 Hz, using the nearest Welch bin to n = f/v (exact f/v: 0.0750, 0.0300, 0.0200):

v (m/s)	nearest Welch bin n (cycles/m)	G_d(n) (m³)	S_z(f) (m²/Hz)
20	0.0742	8.183e-6	0.4091e-6
50	0.0293	52.384e-6	1.0477e-6
75	0.0195	122.970e-6	1.6396e-6

The same wavelengths land at higher frequencies as speed rises.

λ (m)	f at 20 m/s (Hz)	f at 50 m/s (Hz)	f at 75 m/s (Hz)
30	0.67	1.67	2.50
6	3.33	8.33	12.50

Where it stops: the transform assumes constant speed, and its use as a PSD assumes a stationary profile. An accelerating roll is non-stationary, so the 20, 50 and 75 m/s results are snapshots, not a take-off spectrum. DOT/FAA/TC-23/09 (January 2023) notes ProFAA's equations of motion did not allow variable speed (FAA, 2023).

What can a spectrum say about spikes and loads?

Very little about a single bump: a PSD cannot locate a single bump. NLR-TP-2009-190 (June 2009), Section 1, notes that a PSD does not distinguish many small bumps from a few large ones at a given wavelength (NLR, 2009).

Two bumps were inserted in the synthetic case: 6 m × 20 mm at 450 m and 30 m × 25 mm at 900 m. The bumps are half-sines, so a bump of length L is half of a wave of length 2L, and one pulse spreads its energy up to about 1.5/L cycles/m. Adding them raised the fitted level on 0.03 to 0.5 cycles/m from 4.07e-6 to 5.32e-6 m³ (about 30 %, slope −2.01 to −2.19) and the RMS from 14.05 to 14.16 mm. Single bins moved up or down through the cross-term with the random-phase profile (×0.95 and ×0.82 at two illustrative bins near 6 m and 30 m). Nothing in the spectrum says how many bumps there are or where. A rolling chord (a straightedge-like detector, not the Boeing Bump algorithm) found them: the 6 m chord read 19.94 mm at 450 m against a clean maximum of 4.03 mm, and the 30 m chord read 24.08 mm at 900 m against 8.96 mm. See finding the bumps that load the gear.

Load enters through a linear gear model, S_F(f) = |H(f)|² S_z(f), where H is the response from profile to tyre force. For one synthetic two-mass model the RMS from this integral (analytic input PSD) agreed with a time-domain run (one realisation) within about 3 % (3.0 %, −0.2 % and −1.0 % at 20, 50 and 75 m/s). The model has sprung mass 17,000 kg, unsprung mass 400 kg, tyre stiffness 4.0e6 N/m and a damped sprung mode of 1.29 Hz. This is internal consistency, not aircraft validation. An RMS is not a peak: a peak needs a Gaussian assumption that discrete bumps break.

What does a PSD leave out?

• Reach. NLR notes most laser profilometers developed for highway roughness catch wavelengths up to 45 m, whereas an aircraft responds to wavelengths up to about 120 m (NLR, 2009).
• Thresholds. No ICAO, FAA or EASA document sets a PSD limit for a runway profile; see the hub on runway surface thresholds.
• Space-domain index. See Boeing Bump Index, IRI and straightedge compared.
• Operations. This is engineering background, not operational guidance; decisions to operate, restrict or close belong to the aerodrome operator under its State's rules.

What should you keep after running a PSD?

Keep:

• the profile file, sampling interval, line position and survey date;
• the Welch settings, high-pass cut-off and any grade removal;
• the ground speeds chosen and the reason;
• the script and its output.

The procedure, step by step

1. Fix the profile and sampling interval

   State the line, the length and the sampling interval; the shortest resolvable wavelength is twice the interval.

2. Remove the mean and grade

   Detrend the profile and state any high-pass cut-off, because Welch's constant detrend removes only the segment mean.

3. Choose the segment length

   Pick a segment length whose bin spacing 1/(segment length) resolves the wavelengths of interest.

4. Run Welch's method

   Use a Hann window, 50 % overlap, one-sided density scaling and record every setting.

5. Check the level against a known input

   Run a synthetic profile with a known PSD and confirm the recovered level and slope.

6. Choose the ground speeds

   Pick the speeds to examine and say that each is a constant-speed snapshot.

7. Convert wavelength to frequency

   Apply f = v/λ to relabel the wavelength axis in hertz at each speed.

8. Convert the PSD to the temporal form

   Divide the spatial PSD by the speed to get S_z(f) in m²/Hz, keeping both PSDs one-sided.

9. Test for spikes in the space domain

   Locate discrete bumps with a chord or straightedge, because the spectrum averages position away.

Frequently asked

Why is a runway profile PSD in m³?

A one-sided displacement PSD reports variance per unit of spatial frequency. The variance is in m² and the spatial frequency is in cycles/m, so the density is in m² per cycle/m, which simplifies to m³. Its integral over a band of spatial frequencies gives the profile variance in that band, and the square root is the band RMS.

How do you convert runway wavelength to frequency at take-off speed?

Divide the ground speed by the wavelength: f = v/λ. In the synthetic example a 30 m wavelength gives 0.67 Hz at 20 m/s, 1.67 Hz at 50 m/s and 2.50 Hz at 75 m/s. The conversion assumes a constant speed, so during an accelerating roll each result is a snapshot, not a continuous spectrum.

Can a PSD find a bump on a runway?

No. A PSD averages over the whole profile, so it discards chainage, and NLR-TP-2009-190 notes it does not distinguish many small bumps from a few large ones at a given wavelength. In the synthetic case the bumps raised the fitted level by about 30 % and moved single bins up or down, yet the spectrum gave no count or position, while a rolling chord located both.

Where this connects

AC 150/5380-9 explained shows the Boeing Bump method in space, and surveying roughness on the centreline and wheel paths covers collecting the profile. For road profile practice, see from a profiler run to a register value.

References

1. FAA (2009). AC 150/5380-9, Guidelines and Procedures for Measuring Airfield Pavement Roughness. 30 Sep 2009, no revision letter, Active. Original source.Foreword; para 1.2b-g; para 2; 2.1; 2.2a-g; 2.3a-b; 2.4a-d; 2.5a-b; Fig 2-3; Fig 2-4; Ch 3 (3.1, 3.2)
2. NLR (2009). Evaluation Methods for longitudinal evenness of runway pavements: An overview of recent Dutch experiences, NLR-TP-2009-190 (P.A. van Gelder, M.J.A. Stet). June 2009. Original source.Sec 1 pp. 5-6; Sec 2.1 p. 7 and Table 1; 2.2; Elaborate Survey pp. 9-10; 3.1 p. 11; 3.3 p. 13; 3.5 p. 16; 4 p. 18
3. SciPy (2026). scipy.signal.welch. SciPy 1.17. Original source.Parameters: window, nperseg, noverlap, detrend, scaling
4. ISO (2016). ISO 8608:2016, Mechanical vibration, Road surface profiles, Reporting of measured data. 2016 (Ed. 2), confirmed 2021; systematic review may be in progress. Original source.Scope (roads, streets, highways, off-road; excludes rail; runways not mentioned)
5. NASA (1959). Study of Taxiing Problems Associated with Runway Roughness, NASA Memo 2-21-59L (B. Milwitzky). March 1959. Original source.PDF p. 3 (summary), p. 6
6. FAA (2023). Method to Identify and Quantify Multiple-Event Roughness, DOT/FAA/TC-23/09 (T. Gerardi, M. Gerardi, APR Consultants). Final Report, January 2023 (does not constitute FAA policy). Original source.Exec summary pp. ix-x; Sec 2 p. 2; Sec 2.1-2.3 pp. 2-3; Sec 3 p. 4; Sec 3.1 p. 9; Table 1 p. 7; App E pp. E-3/E-4, E-8/E-9

Use the manual edition and project requirements that apply to your work. Worked examples use synthetic values and are not site-specific design advice.