Algorithm Challenge

[ivconic]

This is only for phase, magnitude will need different method.
I used this for my lab instrument for comparing two signals with the same frequency and same amplitude. (audio application)
Zero Crossing method is the closest term that can descibe it.
Iterate through both signals and find the points where their values cross zero.
For each zero crossing in one signal, find the nearest zero crossing in the other signal and calculate the time difference between them.
Convert time difference to phase shift with:
phase_shift = 2 * pi * time_difference * frequency;

float zero_crossing_phase_shift(float *signal1, float *signal2, int n, float frequency)
{
  int i, j;
  float time_diff = 0.0f;
  int found_zero_crossing = 0;

  for (i = 0; i < n - 1; i++)
{
  if (signal1[i] * signal1[i + 1] < 0)
{
  found_zero_crossing = 1;
  for (j = i + 1; j < n - 1; j++)
{
  if (signal2[j] * signal2[j + 1] < 0)
{
  time_diff = (float)(j - i) / frequency;
  break;
}
}
  break;
}
}
  float phase_shift = 2.0f * M_PI * time_diff * frequency;
  return phase_shift;
}​

[Teleno]

This is only for phase, magnitude will need different method.
I used this for my lab instrument for comparing two signals with the same frequency and same amplitude. (audio application)
Zero Crossing method is the closest term that can descibe it.
Iterate through both signals and find the points where their values cross zero.
For each zero crossing in one signal, find the nearest zero crossing in the other signal and calculate the time difference between them.
Convert time difference to phase shift with:
phase_shift = 2 * pi * time_difference * frequency;

float zero_crossing_phase_shift(float *signal1, float *signal2, int n, float frequency)
{
int i, j;
float time_diff = 0.0f;
int found_zero_crossing = 0;

for (i = 0; i < n - 1; i++)
{
if (signal1[i] * signal1[i + 1] < 0)
{
found_zero_crossing = 1;
for (j = i + 1; j < n - 1; j++)
{
if (signal2[j] * signal2[j + 1] < 0)
{
time_diff = (float)(j - i) / frequency;
break;
}
}
break;
}
}
float phase_shift = 2.0f * M_PI * time_diff * frequency;
return phase_shift;
}​

I've asked ChatGPT and it returned quite a similar code, but using the peaks rather zeroes.

To calculate the phase shift between two sinusoidal signals sampled at interval T, you can use the following C code. This assumes that you have arrays (signal1 and signal2) representing the sampled values of the two signals.

#include <stdio.h>
#include <math.h>

#define PI 3.14159265358979323846

double calculatePhaseShift(double *signal1, double *signal2, int numSamples, double T) {
// Find the index corresponding to the peak of the first signal
int index1 = 0;
double max1 = signal1[0];
for (int i = 1; i < numSamples; i++) {
if (signal1[i] > max1) {
max1 = signal1[i];
index1 = i;
}
}

// Find the index corresponding to the peak of the second signal
int index2 = 0;
double max2 = signal2[0];
for (int i = 1; i < numSamples; i++) {
if (signal2[i] > max2) {
max2 = signal2[i];
index2 = i;
}
}

// Calculate the phase shift in radians
double phaseShift = 2 * PI * (index1 - index2) * T / numSamples;

return phaseShift;
}

int main() {
// Example usage
int numSamples = 100; // Replace with the actual number of samples
double T = 0.01; // Replace with the actual sampling interval

// Replace these arrays with your actual sampled signals
double signal1[numSamples];
double signal2[numSamples];

// Fill in your signals with actual data...

double phaseShift = calculatePhaseShift(signal1, signal2, numSamples, T);

printf("Phase shift between the two signals: %f radians\n", phaseShift);

return 0;
}
​

[ivconic]

Yes I am also using AI to make things simpler.
My original code is spaggetttiii and it will need me to spend hours to extract what I wanted to point here.
Anyway; that method is not what Carl asks here. Not good enough.
...
Previous one was just illustrative example.
Not perfect for direct implementation into md code because of the given conditions.
Conditions would be the same frequency (alright) and the same amplitude (will need additional hardware, bit tricky, not that trivial as it may seems).
I can't use Car's fancy symobols (Carl what you use to write math symbols into post?)
So it can be also expressed like this:
magnitude = sqrt(X^2 + R^2)
phase = arctan(R / X)
We don't have much of a choice here. 3 good methods to skin a cat only, though.
Direct Calculation: simple and efficient approach. (This is my favorite, because is "lightweight" for implementation and even can be tried with Atmega328P using IRQ)
FFT: requires more processing power and memory compared to direct calculation.
DSP: Several dedicated DSP libraries offer pre-built functions for calculating magnitude, phase, and other signal processing tasks.
These libraries can simplify implementation but may require additional hardware or software resources.
Or use dedicated DSP chip... ha! Opening another Pandora's box here!

Since I will rather try it on Atmega 328P and try to keep it as simple as possible; I vote for easiest one, the direct calculation.

// Function to calculate magnitude
float calculate_magnitude(float x, float r)
{
  return sqrt(pow(x, 2) + pow(r, 2));
}

// Function to calculate phase
float calculate_phase(float x, float r)
{
  return atan2(r, x);
}

int main()
{
  // Define X and R signal values
  float x_signal = 1.234;
  float r_signal = 5.678;

  // Calculate magnitude
  float magnitude = calculate_magnitude(x_signal, r_signal);
  printf("Magnitude: %.3f\n", magnitude);

  // Calculate phase in radians
  float phase_radians = calculate_phase(x_signal, r_signal);
  printf("Phase (radians): %.3f\n", phase_radians);

  // Convert phase to degrees (optional)
  float phase_degrees = phase_radians * 180.0f / M_PI;
  printf("Phase (degrees): %.3f\n", phase_degrees);

  return 0;
}​

[ivconic]

Why did I mention IRQ?
So if we do calculations all the time; even a much stronger processor than the Atmega328P will be "choked" and slowed down very quickly.
Maybe multi core will not.
We don't really need to do calculations all the time.
I suggest thinking of this in two steps. The first step is to simply recognize whether the target is ferrous or non-ferrous.
Only if the target is non-ferrous to "call" the IRQ and start the calculation.
All we have to do is come up with the simplest and easiest (for mcu) method to instantly recognize whether the target is ferrous or non-ferrous. (Zero Crossing method now makes sence, right?)
(I have something else in mind with this, who read my last post on another topic; has a chance to understand my intentions)

P.S.
And there! You've got me thinking about CPU usage again! We've slipped on the AMX road again!

[Carl-NC]

A zero-crossing algorithm only works for large-enough signals. Imagine the detector has a 100mV coil null signal and the target signal is 50mV, it won't work. That's why we demodulate and use high pass filters. The algorithm needs to work using the X & R signals.

[Carl-NC]

another solution ...
fast artan calculation using vector rotation.

So far I like this the best. It's called a CORDIC method and is basically a successive approximation algorithm. No multiplies or divides. It can be written for other common trig functions as well.

[ivconic]

A zero-crossing algorithm only works for large-enough signals. Imagine the detector has a 100mV coil null signal and the target signal is 50mV, it won't work. That's why we demodulate and use high pass filters. The algorithm needs to work using the X & R signals.

When zeroing the coils, we always tried to achieve a balance as close to zero as possible.
With hand held meters I have always found 7-14mV to be a well balanced coil.
Using the same measurement method on the factory coil, I would measure 0v, which corresponds to the logic that the DIY coil is always worse balanced.
Of course... how much to trust those measuring instruments.
But 100mV!
You probably mean the detectors from 50-60 years ago?
...
I was hoping for analog solutions... but somehow everything slipped back into the domain of programming.
Ok, that's what you called the topic.
Is there anyone interested in analog solutions?
If this all boils down to software solutions; it's not my cup of tea.
I am exclusively interested in an auto zero solution with analog circuits.
But to be advanced and efficient.
I won't bother you anymore, feel free to continue.

[moodz]

...
I was hoping for analog solutions... but somehow everything slipped back into the domain of programming.

Its been done by GARRETT .. patent US9989663B1 .. expires in 2036.

https://patents.google.com/patent/US9989663B1/en

[ivconic]

Its been done by GARRETT .. patent US9989663B1 .. expires in 2036.

https://patents.google.com/patent/US9989663B1/en

Oh sheeeee..ttt!
What a misunderstanding!!!
If you didn't mention that patent; I wouldn't have understood that all the time it was a misunderstanding.
I didn't mean auto balancing at coil!
How did this go that way??
All the time I think about the "auto-zero" circuit in non-motion detectors, whether it is IB or PI (probably a different morphology).
Ok, let's go slowly, so we have a non motion detector. It usually has a button next to which it often says "tune" or "retune" or "reset"...
When the detector starts to drift; you press that button and it calms down.
What is the name of the circuit that does that?
I'm interested in that circuit, to design it better, to be more perfect, to constantly keep the detector in "balance", whenever it drifts; for that circuit to return it.
So that we don't have to press that button often.
It would be cool to make a concept of such a circuit for both PI and IB detectors. Let's say Delta Pulse knows to drift and PS2 as well,
I remember that Cscope 1220B was drifting so I often pressed that button. And there are countless other examples.
Every time I say "auto-zero"; I mean that, to such a stage in the detector, ok.

[moodz]

So far I like this the best. It's called a CORDIC method and is basically a successive approximation algorithm. No multiplies or divides. It can be written for other common trig functions as well.

Its accurate to .06 degrees approx error .... ( calculation not VDI ) .... The Theta value array are the rotation steps direclty ( eg 4500 is 45.00 degrees, 2657 is 26.57 degrees, etc )

The X and Y values are 8 bit signed values.

[moodz]

When zeroing the coils, we always tried to achieve a balance as close to zero as possible.
With hand held meters I have always found 7-14mV to be a well balanced coil.
Using the same measurement method on the factory coil, I would measure 0v, which corresponds to the logic that the DIY coil is always worse balanced.
Of course... how much to trust those measuring instruments.
But 100mV!
You probably mean the detectors from 50-60 years ago?
...
I was hoping for analog solutions... but somehow everything slipped back into the domain of programming.
Ok, that's what you called the topic.
Is there anyone interested in analog solutions?
If this all boils down to software solutions; it's not my cup of tea.
I am exclusively interested in an auto zero solution with analog circuits.
But to be advanced and efficient.
I won't bother you anymore, feel free to continue.

So you want an auto zero function
....I was pretty sure you were talking about the coil balancing here ... sorry

[ivconic]

So you want an auto zero function
....I was pretty sure you were talking about the coil balancing here ... sorry

No, Carl intentionally started this fuzz!



This is an attempt, for full-scale range input, obviously it will be corrected and adapted to real x, y values, also depending on hardware and what levels coming from it, yet just an example.
I am tempted to see how Atmega328P will deal with this, will it choke or slow down?
Don't have prepared hardware though...

#include <Arduino.h>
#include <math.h>

// Sensor and ADC information
const float InputRange = 5.0f; // Full-scale range of the input
const int analogResolution = 10; // ADC resolution

// Scaling factor (optional)
const float scalingFactor = 10.0f;

// Define analog pins for x and y readings
const int xPin = A0;
const int yPin = A1;
// Define the K array for CORDIC algorithm
static const float K[8] = {
  0.6072529350088812f,
  0.4142135623730950f,
  0.2500000000000000f,
  0.1591549430918954f,
  0.1023381623379956f,
  0.0645161290312500f,
  0.0408407044921875f,
  0.0255017857524998f
};
// Fast arctan approximation using CORDIC algorithm
float fast_atan(float y, float x) {
  // Initialize variables
  float z = 1.0f;
  float angle = 0.0f;

  // Check for special cases
  if (y == 0.0f) {
    if (x > 0.0f) {
      return M_PI_2;
    } else if (x < 0.0f) {
      return -M_PI_2;
    } else {
      return 0.0f;
    }
  }

  // Perform CORDIC iterations
  for (int i = 0; i < 8; i++) {
    float d = x - K[i] * z;
    if (d < 0.0f) {
      angle += K[i];
      x = d;
      z += K[i] * y;
    } else {
      angle -= K[i];
      x = d;
      z -= K[i] * y;
    }
  }

  // Final adjustment and correction for quadrant
  if (y < 0.0f) {
    angle = -angle;
  }
  if (x < 0.0f) {
    angle += M_PI;
  }
  return angle;
}

void setup() {
  // Configure analog pins as inputs
  pinMode(xPin, INPUT);
  pinMode(yPin, INPUT);
}

void loop() {
  // Read analog values
  int x_raw = analogRead(xPin);
  int y_raw = analogRead(yPin);

  // Convert raw values to voltage
  float x = (float)x_raw * (InputRange / pow(2, analogResolution));
  float y = (float)y_raw * (InputRange / pow(2, analogResolution));

  // Apply scaling factor (optional)
  x *= scalingFactor;
  y *= scalingFactor;

  // Calculate arctan
  float angle = fast_atan(y, x);

  // ... further processing with the calculated angle ...
  // Example: print the angle to the serial monitor
  Serial.print("Angle: ");
  Serial.print(angle * 180.0f / M_PI);
  Serial.println(" degrees");

  // Delay to slow down the loop execution
  delay(100);
}
​

[moodz]

There are auto tune implementations ... heres one from Garrett again ( patent expired ).

[ivconic]

I know that, not good enough!
There is slightly better one in White's Eagle series, forgot the exact schmatic, but I am sure I have it... also not good enough!
Dang!
That's why I proposed that stage as separate project here.

[Carl-NC]

When zeroing the coils, we always tried to achieve a balance as close to zero as possible.
With hand held meters I have always found 7-14mV to be a well balanced coil.

OK, suppose you have 100mV of ground signal and a 50mV target signal. You can't extract the target phase with a zero-crossing algorithm.

I was hoping for analog solutions... but somehow everything slipped back into the domain of programming.
Ok, that's what you called the topic.

Yes, I specifically mentioned an algorithm for a micro. If you want to talk about analog solutions, then post a challenge/competition for that. Same for the autotune circuit. Or electronically-nulled coils.