Add another goertzels example

README.md

@@ -96,6 +96,13 @@ A lot of times, you will want to use a sinewave, or need to do DFT for specific
 And sin/cos may not be fast.  Well, good news, goertzels is to the rescue.  You can
 generate a sin/cos without any sin or cos!
 
+Note that in the below sample, the power will be 0, because the sample is 1.0, so it
+is treated as a "DC offset" and thus gets nulled out.
+
+Also, goertzels does not need to compute an I and a Q to detect power.
+
+Note that you can store the "coeff" in a table for lookup, it is fixed per the frequency you are observing.
+
 ```c
 #include <stdio.h>
 #include <math.h>
@@ -114,17 +121,69 @@ int main()
 	float sprev2 = 0;
 	for( i = 0; i < numSamples; i++ )
 	{
-		float SAMPLE = 1.0; // If you wanted to do a DFT, set SAMPLE to your incoming sample.
-		float s = SAMPLE * coeff * sprev - sprev2; // Here is where the magic happens.
+		// If you wanted to do a DFT, set SAMPLE to your incoming sample.
+		float SAMPLE = 1.0;
+
+		// Here is where the magic happens.
+		float s = SAMPLE + coeff * sprev - sprev2;
 		sprev2 = sprev;
 		sprev = s;
 		printf( "%f\n", s ); 
 	}
+
+	// For DFT, your power will be:
+	float power = sprev*sprev + sprev2*sprev2 - (coeff * sprev * sprev2);
+	printf( "%f\n", power );
 }
 ```
 
-![goertzels](https://github.com/cnlohr/cnlohr_tricks/blob/master/media/goertzels0.png?raw=true)
+![goertzels0](https://github.com/cnlohr/cnlohr_tricks/blob/master/media/goertzels0.png?raw=true)
+
+Below is an example using Goertzel's Algorithm to extract the phase and magnitude of a sinewave of an incoming signal, against a specific target frequency. If you think this is magical, that makes two of us.
+
+```c
+#include <stdio.h>
+#include <math.h>
+
+int main()
+{
+	const float omegaPerSample = 0.015708; // pi / 200
+	const int numSamples = 400; // enough to go from 0 to 2pi
+
+
+	for( float phaseOffsetOfSignal = 0; phaseOffsetOfSignal < 3.1415926*2.0; phaseOffsetOfSignal += 0.01 )
+	{
+		float coeff = 2 * cos( omegaPerSample );
+		int i;
+
+		// TRICKY: When you want a sinewave, initialize with omegaPerSample.  This
+		// is crucial.  The initial state will have massive consequences.
+		float sprev = omegaPerSample;
+		float sprev2 = 0;
+		for( i = 0; i < numSamples; i++ )
+		{
+			// If you wanted to do a DFT, set SAMPLE to your incoming sample.
+			float SAMPLE = sin( phaseOffsetOfSignal + i * omegaPerSample );
+
+			// Here is where the magic happens.
+			float s = SAMPLE + coeff * sprev - sprev2;
+			sprev2 = sprev;
+			sprev = s;
+		}
+
+		// For DFT, your power will be:
+		float power = sprev*sprev + sprev2*sprev2 - (coeff * sprev * sprev2);
+		//printf( "Power: %f\n", power );
+
+		float coeff_s = 2 * sin( omegaPerSample );
+		double rR = 0.5 * coeff   * sprev - sprev2;
+		double rI = 0.5 * coeff_s * sprev;
+		printf( "%f,%f\n", rR, rI);
+	}
+}
+```
 
+![goertzels1](https://github.com/cnlohr/cnlohr_tricks/blob/master/media/goertzels1.png?raw=true)
 
 ## Linux