From 0404a1ad77f66c044d1cc7d6fbbb92b0dbecabf6 Mon Sep 17 00:00:00 2001
From: David Rotermund <54365609+davrot@users.noreply.github.com>
Date: Thu, 16 Nov 2023 17:25:38 +0100
Subject: [PATCH] Update README.md

Signed-off-by: David Rotermund <54365609+davrot@users.noreply.github.com>
---
 numpy_fft_1/README.md | 19 +++++++++++++++++++
 1 file changed, 19 insertions(+)

diff --git a/numpy_fft_1/README.md b/numpy_fft_1/README.md
index 1c7e5e5..4793318 100644
--- a/numpy_fft_1/README.md
+++ b/numpy_fft_1/README.md
@@ -54,3 +54,22 @@ Since we deal with a 1 dimensional time series
 y_fft: np.ndarray = np.fft.rfft(y)
 frequency_axis: np.ndarray = np.fft.rfftfreq(y.shape[0]) * sampling_frequency
 ```
+
+## Calculating a [normalized](https://de.mathworks.com/help/signal/ug/power-spectral-density-estimates-using-fft.html) power spectral density
+
+The goal is to produce a power spectral density that is compatible with the [Parseval's identity](https://en.wikipedia.org/wiki/Parseval%27s_identity). Or in other words: the sum over the power spectrum without the zero frequency has the same value as the variance of the time series.  
+
+```python
+y_power: np.ndarray = (1 / (sampling_frequency * y.shape[0])) * np.abs(y_fft) ** 2
+y_power[1:-2] *= 2
+
+if frequency_axis[-1] != (sampling_frequency / 2.0):
+    y_power[-1] *= 2
+```
+
+Check of the normalization:
+
+```python
+print(y_power[1:].sum())  # -> 0.5
+print(np.var(y))  # -> 0.5
+```