Skip to content

This issue was moved to a discussion.

You can continue the conversation there. Go to discussion →

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

is it time to improve Circling Wind calculation? #1388

Closed
bomilkar opened this issue Mar 26, 2024 · 50 comments
Closed

is it time to improve Circling Wind calculation? #1388

bomilkar opened this issue Mar 26, 2024 · 50 comments

Comments

@bomilkar
Copy link
Contributor

XCSoar version

7.42

What should XCSoar do differently, what functionality should be added?

The current implementation of the circling wind calculation has some room for improvement. The comments in src/Computer/Wind/CirclingWind.cpp suggest this already. For instance, it assumes “the pilot flies in perfect circles” and “with constant airspeed”. Also, the calculation of the wind bearing can be improved.

What are you trying to do, what is the use case for the suggestion?

I have a prototype running which implements the following:

  1. it uses changes in turn rate to validate the assumption of perfect circles.

  2. it uses true airspeed to compensate for changes.

  3. it uses an iterative fitting to determine the bearing of the wind.

  4. to evaluate the level of perfection of the circle, it uses the track information from the GPS. This works reasonably well as long as the wind speed is small enough compared to TAS (maybe 1/3). GPS receivers are black boxes wrt their ability to cope with changes in track and ground speed: how quickly will the values in the GPRMC sentence follow the real movements of the ship? However, future implementation may use turn rate information from gyros of simple/inexpensive IMUs. That change will make sense as soon as IMUs have a better support in XCSoar and are more widely available.

  5. my assumption is that airspeed information is available in many installations. So it uses ground_speed – airspeed instead of just ground_speed from the GPS reciever. However, the prototype works even if TAS is not available (or derived). In that case it lacks that compensation.

  6. in an ideal case, the speed oscillates in a sine curve. So fitting the actual speed values of a given circle to a sine curve usually provides better results for the wind direction.

I’m running the prototype in Replay mode using files from the “NMEA logger” which I usually capture during my flights. So far results look promising.

How shall I proceed? Shall I post a PR? I would like to get your feedback, even though the implementation is far from final.

@lordfolken
Copy link
Contributor

In general this sounds interesting. Open a dratft PR.
Please follow xcsoars git rules. Small building block changes that build with well crafted commit messages.
https://xcsoar.readthedocs.io/en/latest/policy.html#code-style

@kobedegeest
Copy link
Contributor

Just my thoughts

on point one, to get change in turn rate from position is like third order derivative, taking derivatives increases noise level.

On point 2, I think the reason the current circling wind implementation does not use airspeed is because you don't need circling wind if you have airspeed.

Point 5, if you assume everyone has airspeed you might as well assume everyone has wind data already as well. As most varios who measure airspeed also calculate wind. But I think it is irrelevant for circling wind as it is a methode to get an estimate of the wind when you don't have this data.

Point 6 is completly correct i think.

avg(v_x) = w_x and avg(v_y) = w_y in general i think whether you fly circle or ellipse or any other closed shape. Only error there is pilot choosing to move the circle. And i don't think there really is a way of knowing whether pilot moves the circle or the wind moves the circle (with just gps data).

@bomilkar
Copy link
Contributor Author

Thanks for the feedback, @kobedegeest

I also have doubts about what a GPS receiver can reliably deliver. As you said, it is a third order derivative and the quality of the results can vary greatly depending on operating conditions, brand and type. For this reason, I have started to analyze the replay results of different flights. Although it looks promising, there are still cases where I get unexpected results from CalcWind.

I don't buy your assumption "avg(v_x) = w_x and avg(v_y) = w_y in general i think whether you fly circle or ellipse or any other closed shape". The effect is that a perfect circle in the air is an ellipse on the ground. The shape of the ellipse can be used to determine the wind. That's the underlying principle of circling wind, I guess.

I remain convinced that you can distinguish between a (near) perfect circle and an arbitrary pattern (to move the center point while circling) based on the ship's rate of turn. Certain GPS receivers may not provide sufficient track data to derive a useful rate of turn. But even simple IMUs with their gyros should suffice. Fortunately, my algorithm doesn't care where the turn rate comes from.

@bomilkar
Copy link
Contributor Author

One more point, which might also affect the ZigZag algorithm: ground speed and TAS might have different delays in the signal processing chain. In my case it seems to be around 2 sec, which I noticed looking at the InfoBoxes when I replayed a takeoff.
TAS comes 2 sec after VGND. While circling that translates to ~30° difference in heading.

@kobedegeest
Copy link
Contributor

Maybe i should clarify that when i said closed shape i mean that in the reference frame of the air of course. In which case your average ground speed vector will be zero if there is no wind and will be equal to wind speed if there is wind since you keep stationary w.r.t. the air but air moves w.r.t. ground. So my reasoning is that you try to keep stationary as you try to keep flying around center of thermal and the thermal drifts so you drift with it and in that case avg(v_x) = w_x and avg(v_y) = w_y is correct. But if the pilot decides to move its center by 1m west every turn there is no way to know this is the pilots chose and not the drift of the wind. So best you can do is throw out circles where big chance in center is found and hope the rest averages out.

Current implementation also weights first turns less i assume because at that point your still finding the center.

I remain convinced that you can distinguish between a (near) perfect circle and an arbitrary pattern (to move the center point while circling) based on the ship's rate of turn.

I am quit interested in knowing how you think this can be done or knowing how you do it.

Gyros are useful cus they provide info on change in inclination of the plane from which you can estimate current inclination from which you can get predicated turn rate and radius and from the error you can get the wind. But once you have gyros and airspeed you get algorithms with extended kalman filters used in like hawk, larus and anemoi. And once you have gyros i totaly agree you can differentiated between pilot and wind, because you can measure change in inclination and pitch not just because you get turn rate.

Fortunately, my algorithm doesn't care where the turn rate comes from.

It should cus if you have gyros you can do so much better than if you don't have them.

I don't want to critic to much or discourage the effort as i also find the circling wind to be not too good in flight. But the answer should not be use more data (which often is not present). And also like you mention there might be delay between different data especially if you collect data from external devices and then do the calculation in XCSoar.

@bomilkar
Copy link
Contributor Author

I'm trying to get away with GPS (and optional airspeed) only. But I may not succeed.
A good percentage of Android phones have built in IMUs. But I have no clue how to get (raw) data from it. It should be possible, though,
I don't want to stick a ton of data into a kalman filter (with 10 - 15 dimensions). I'd rather use my physical understanding of the problem instead of brute force compute power.
For me it's not just the first few circles where I try to find the best climb rate. I adjust quite often. And I found the change of turn rate (even if it's 3rd order derivative from GPS) to be quite good to identify "perfect" circles.
I'll post a draft PR in the coming days. I just to run a few more flights with it. Unfortunately speeding up the replay doesn't work.
BTW: I really appreciate the discussion! Are you the author of the current implementation of CirclingWind ?

@bomilkar
Copy link
Contributor Author

I'm submitting the code changes.
First draft PR #1390 is to compare the results of different wind algorithms. It is in a CSV format to load for processing after the replay of the nmea file.
Columns are:

  1. time in seconds since 00:00 UTC
  2. flight mode: 1 == circling, 0 == forward
  3. wind speed in user units
  4. wind.bearing in degrees
  5. baro_altitude in user units

This PR is for analyzing results, not intended to include in production.

@bomilkar
Copy link
Contributor Author

Got stuck with github (again).
I created a new branch bomilkar:new_circling_wind and added 4 commits:
ddcf81e
37044ee
cb56627
1fb6a91

The PRs build on top of each other.
After the first commit I created PR #1390
Now that I'm at this point I can't create any more PRs (for any of the remaining commits).
And I can't figure out what's wrong! Any idea?

@bomilkar
Copy link
Contributor Author

bomilkar commented Apr 1, 2024

Just complete the 2nd attempt to submit the stuff.
#1391
To be safe, I included the 2 commits in one PR. (Had trouble with that on the first attempt.)
Any feedback is highly appreciated!

@bomilkar
Copy link
Contributor Author

bomilkar commented Apr 1, 2024

Use environment variable to select between the new and the old algorithm:
export CIRCLING_WIND_EXPERIMENTAL=1
or equivalent on other platforms.

@lordfolken
Copy link
Contributor

You create a PullRequest to a specific branch of yours.
If you add commits to your branch, the PR will be updated automatically.
You cannot create mulitple PRs to the same branch.

If you need to fix stuff in the current commit:
git add
git commit --amend (will merge into the previous commit)

if you need edit something multiple git commit back:
git rebase -i HEAD~5 (where 5 is the ammount of previous commits)

You can then write edit or pick or merge infront of the commit, and it will dump you to that point in the history.
you can make changes and then commit --amend to that commit.
to continue you can do git rebase --continue

All of this is reflected in
git status

Important: If you edit your git history you might have to force push to the remote (replacing your current commits with the new ids)
so git push --force

@bomilkar
Copy link
Contributor Author

bomilkar commented Apr 2, 2024

Thank you very much, @lordfolken !
Is there a good (and concise) description on how to use git and github together?

@kobedegeest
Copy link
Contributor

I believe there is a mistake in your reasoning. You look at change in direction of track to determine quality of circle. You assume that this should be a constant. But track is Heading plus Wind. And indeed change in heading when flying nice circles is constant but that does not imply that change in track is constant just like your ground speed changes during a circle due to wind also the change in angle track changes. (Unless you do account for that somewhere and i missed it)

@bomilkar
Copy link
Contributor Author

bomilkar commented Apr 2, 2024

I believe there is a mistake in your reasoning. You look at change in direction of track to determine quality of circle. You assume that this should be a constant. But track is Heading plus Wind. And indeed change in heading when flying nice circles is constant but that does not imply that change in track is constant just like your ground speed changes during a circle due to wind also the change in angle track changes. (Unless you do account for that somewhere and i missed it)

I'm aware of this. The changes should be small, not 0. And it works as long as wind_speed is sufficiently low. That's why I wrote "wind_speed < TAS/3" . That's what I confirm in my replays.
One way to account for this is to use turn_rate from a gyro. But gyros aren't widely available.
One other way would be to look at the max difference between ground_speed and TAS. That gives an indication if method is still valid. (Not yet implemented.)
The nice thing about the new algorithm: it need just GPS. Airspeed sensor is good to have. And just hoping the "pilot flies perfect circles" and "keeps constant speed" is addressed and no longer just an invalid assumption.

@kobedegeest
Copy link
Contributor

So you have

int quality = 5; // top quality, perfect circle if (quality_metric > 0.2) quality = 4; if (quality_metric > 0.3) quality = 3; if (quality_metric > 0.4) quality = 2; if (quality_metric > 0.5) quality = 1; if (quality_metric > 0.7) return Result(0);

with quality_metric equal to (max_d_alpha - avg_d_alpha)/avg_d_alpha
that is your check on turn rate to see if you have good circle right?
you mention this wind_speed < TAS/3 how did you come up with that number? Did you check what the "quality metric" is at that condition? If you didn't I can tell you it's 0.5.

The flaw with this metric is that regardless of how good you fly this metric becomes worse with increasing wind speeds. Whereas the need for good calculation (i would say) increases with stronger winds. At wind_speed = TAS/2 this "quality_metric" becomes 1. Your 0.7 value is reached around 40%. So someone thermaling at 80km/h cant measure winds larger than 32km/h.

So I think everyone can agree that this quality metric has not much to do with how good your circle is.

Again feel free to correct me if i miss understand the code, i must admit i don't fully know where or how the quality is used.

@bomilkar
Copy link
Contributor Author

bomilkar commented Apr 2, 2024

You can calculate the quality metric on paper using geometry and trigonometry, But that doesn't take into account what happens inside the GPS receiver. We discussed that earlier.
I have no flight logs for wind speeds higher than 35 km/h. That's how I determine "wind_speed < TAS/3". However, it doesn't mean it won't work at higher wind speeds.
Best is to replay various flight logs and review the results. I encourage you to do that if you have nmea files with air speed.

The quality number which is fed into Store has 6 distinct levels: 0 - 5. I translate the "roundness" of the circle (circle_quality_metric) and the fit to the cosine curve (min_fit_metric) into these 6 level. I have to admit, there is no elaborate theory behind that. And there should be more experiments to validate or improve that. Food for thought!

@kobedegeest
Copy link
Contributor

What i currently did is just straigth up calculate/ simulate what the formula would give (input heading and airspeed but not use airspeed as data). Way more time efficient than replaying flight logs, and easier than working out the equations delta_angle_track is not that nice of a formula :D . Plus can actualy check if the resulting wind is equal to the input wind. Plus you can input whatever error you want to check on the data noise on airspeed noise on gps position etc and evaluate exaclty how the calculation changes based on that.

And currently i must say the formula you used for the magnitude of wind speed seems to be quite good actualy. And i think much better than the current methode (have not tested current methode). Btw did you look into just fitting everthing to the cosine now you just get angle from the cosine fit could also fit magnitude.

But having a qualty metric that becomes worse just because wind gets stronger (or flight speed goes down) is not a good design. Espacily since the reverse is true for example if you have random variation on air speed the mistake on calculated wind goes down as wind increases.

I understand how you calculate the 6 distinct levels. The part that i do not know is what you use the level for. Like if it just for some weighted averaging of multiple circles it might not matter as every turn will have the same bad quality if wind is high. But it also looks like you just dont do the calculation if it is in the bad qualty.

And note the values given before are for perfect circle with no noise on the angle. Having some small variation like 10% noise on d heading is like adding 0.1 to that number.

But yeah purely the calculation for the magnitude and direction do seem to be quite good.

@bomilkar
Copy link
Contributor Author

bomilkar commented Apr 3, 2024

Do you have a testbench ready to use? That was one of the next things on my to-do-list and would be happy to use yours if you want to share.
The testbench should instantiate the DUT, pretty much the same way as in Computer.cpp together with Store. And a testvector generator, which should be able to generate deterministic and random pattern, mostly to see how it behaves in corner cases. It should also be able to generate pattern from an NMEA file. (But much faster than replay!) However, deterministic pattern will most likely not reflect what happens inside the GPS receiver. I guess it derives track and ground speed from the fixes (3rd order derivative...). But I don't think the manufacturers will share a spec. Maybe high precision GNSS receivers use other methods (like doppler) which might deliver much better track and speed results. There is no choice than to use NMEA files to simulate realistic scenarios.

If nothing remains from my experiments, I think at least the method to calculate the magnitude of wind speed should go into CirclingWind.cpp . It is indeed much better. In fact, it is good enough so that an additional fitting loop will probably not improve the result (beyond the noise level). BTW: look at CirclingWind2.cpp line 223. There is the TO DO note on that subject.

@bomilkar
Copy link
Contributor Author

bomilkar commented Apr 3, 2024

Some thoughts on realistic scenarios:
Flights where (wind_speed/2 >= circling_TAS) are very rare. (I'm not talking about waves in the mountains.) With a Ventus you will not circle below 100 km/h IAS if the wind speed exceeds 50 km/h. It'll be so bumpy, that you better have good speed margin not to stall or fly in the bad section of the polar. And I don't think very many pilots will fly a Ka8 at 50 km/h wind. And if they do, their circling percentage will have to be VERY small otherwise they fly backwards.
Can you use your testbench to figure out how the algorithm performs at (wind_speed/3 < circling_TAS < wind_speed/2)? That'll be interesting.
Another point to bear in mind: The horizontal wind speed in the thermal is not the same as between thermals. This primarily has to do with momentum conservation. However, the wind speed between thermals is important for cross country performance.

@kobedegeest
Copy link
Contributor

kobedegeest commented Apr 3, 2024

I just have some matlab script actualy :D just cuz i was to lazy to work out the equations past d_angle_track which already had way to many sines and cosines in it. And initially that was what i wanted to check the effect on d_angle_track and whether or not this is constant. And like i said its not.

So i don't have a test bench setup like you discribe. Just looking at the interesting mathematical problem, my c++ and github skils are quit limited actually.

After i got d_angle_track i looked at the formula you used to calculate the wind speed. And indeed fitting a cosine to get the amplitude does not realy improve upon it. Where i think it might out perform is when the sampling rate is not constant and not uniform over the circle lets say in one half of the circle the reception is worse and you miss 50% of the gps fixes.

Can you use your testbench to figure out how the algorithm performs at (wind_speed/3 < circling_TAS < wind_speed/2)?

Yes this i can (but i assume you meant to swap TAS and wind speed?) do from what i have tested up to now i would say if that it is quit good if windspeed is larger than TAS/10 . Below that random noise on the speed or gps fix or turn rate gives really large error on wind speed. above that +-5% error on any or all parameters keeps the error within +-10% and will average out over multiple circles. But below measured wind seems to always be larger than the actual wind. This is also where fitting a cosine out performs the formula.

So in using the formula their seems to be an over estimation of windspeed at low windspeeds with noisy data. When fitting cosine this is not the case and the error is more centered around the correct value so will average to correct wind speed. Don't know how hard it is to find best fit in c++ have to fit 3 parameters in the cosinefunction.

So for now i mainly looked at just adding random noise to the data not really at conscious decision of pilot to move the center. That is something i still want to do, that is also the part that should be attempted to filter out. I hope this explanation already gave some extra insight.

@bomilkar
Copy link
Contributor Author

bomilkar commented Apr 3, 2024

"you miss 50% of the gps fixes" that's very unlikely in any scenario. So the sampling rate is constant at 1 per sec.
I also noticed the results getting worse at lower wind speed. Strange. I need to figure out what happens. It's clear the wind bearing is meaningless at low speed (below 1 m/s). But above that it should work. Maybe a mistake in the code?

@kobedegeest
Copy link
Contributor

Maybe a mistake in the code?

No i just looked at the math behind it. So i did not actualy check your implemantation. So if we both notice that i am pretty sure it is in the formula used, but also don't know yet why. But the error is quite "constant" don't know how to say. If i sweep wind speed from 0/TAS to 0.5/TAS the calculated wind at low speeds converges to a value scaling with noise on the data +-5% noise wind/TAS converges to 0.05, +-10% noise it converges to 0.1. And it looks to come from noise on the speed so either pilot not flying cst speed or noise due to gps data. Noise on the angle of heading does not influence wind at low speed so pilot not having cost turn rate.

So the sampling rate is constant at 1 per sec.

Yeah i noticed in your code that you don't seem to look at time between "samples" as you call it and thus assume it to be the same between all of them. But is it though? In real cases? I am looking maybe more at the worst case scenarios but if you have perfect 1sec sampling rate and fly nice looking circles any algorithm will give good results.

"you miss 50% of the gps fixes" that's very unlikely in any scenario.

What i meant was 50% in half of the circle so 25% in total but even if it is less than that. Jus having a non constant sampling rate over the circle will give error in the averaging if you have your cosine and you have more data points on the positive part of the cosine than on the negative side the average will be off. In my experience its not unlikely that gps reception is not equal good across the full 360° turn depending on your gps reciever position. Just pointing out potential pitfalls of the method.

@bomilkar
Copy link
Contributor Author

bomilkar commented Apr 4, 2024

There is/was a mistake in the code. And I fixed it. I'll validate a bit more and then push, maybe later today.

Line 153 checks the avg time between the samples. But you are right, I need to check if that changes too much between samples. So far I came to the conclusion the algorithm works for any step width as long as it's constant over the circle (and not to far apart). Agree?

@bomilkar
Copy link
Contributor Author

bomilkar commented Apr 4, 2024

A related question: what if there are 2 or more GPS receivers feeding into XCSoar? Any wind algorithm will most likely produces bogus numbers!
I don't understand the priority mechanism behind the device list enough to be sure. But I guess if it sees GPS fixes from one device, it'll ignore the fixes from devices lower on the list. @lordfolken Right?

@kobedegeest
Copy link
Contributor

I do know whats going on at low wind speeds, you take the average of absolute value of noise and than multiply by 1.5708. So with zero or low wind any variation in speed will contribute a positive effect to the wind.

My looking at a fitted cosine would give negative wind speed half of the time (which i now realise is not physicly meaningfull)
I assume for low wind speed you would find random directions. How is the wind averaged over multiple circles sum of the vectors? Than the problem at low speeds would be averaged out i guess. Also your cosine fit would give it a bad qualty metric.

So far I came to the conclusion the algorithm works for any step width as long as it's constant over the circle (and not to far apart). Agree?

Yes 10 samples over the circle seems to be sufficient. But i only checked the formula for speed for now.

A related question: what if there are 2 or more GPS receivers feeding into XCSoar?

First device in the list that provides info will be used so if device B and C both give gps only gps of B will be used.

@bomilkar
Copy link
Contributor Author

bomilkar commented Apr 4, 2024

Here is the essence of the fix:

double
CirclingWind2::FitCosine(const size_t n_samples, const double amplitude, const double offset, const Angle phase)
// fit the measured speed data to an inverse cosine
// uses the "sum of squares" algorithm
{
  double sum_of_squares_of_deltas = 0;

  for (size_t i = 0; i < n_samples; i++) {
    const double a = (samples[i].track - phase).AsBearing().Radians();
    double net_speed_diff = samples[i].ground_speed - samples[i].tas - offset;
    // to make the fit metric (somewhat) independent of the amplitude
    // wind speeds less than 1 m/s are not interesting anyways
    if (amplitude > 1.0) net_speed_diff /= amplitude;
    sum_of_squares_of_deltas += Square(-fastcosine(a)-net_speed_diff);
  }
  return sum_of_squares_of_deltas;
}

The push will follow later.

@kobedegeest
Copy link
Contributor

What's your opinion on w_x = mean(v_x), w_y = mean(v_y). I know you said before that you don't believe that but it is in fact true.
The error at low wind speeds is less as there is no absolute value taking before averaging and it is computationally easier

You now take two averages one of speed to get offset and then one of abs(speed - offset) and that just gives the strength.
By taking average of v_x and v_y you take two averages and find both strength and angle. And then a sum of squares of the expected cosine (no need to fit just check as we know everything) can give a quality metric. (v_x and v_y are ground speed)

@bomilkar
Copy link
Contributor Author

bomilkar commented Apr 4, 2024

Consider this: average of angle == 180° if it's a full circle. No matter what.
You can compute w_x = mean(v_x), w_y = mean(v_y) but there is no information on the associated angle. Therefor w_x / w_y has nothing to do with a direction of the wind.

@kobedegeest
Copy link
Contributor

average of angle does not matter. The reason w_x = mean(v_x), w_y = mean(v_y) works is because overtime you drift in the direction of the wind. w_y/w_x is tangens of the direction of the wind. I can write down a proof as well if you like? To show this to be true for full circle.

Based on your answer i assume that once i convince you of the math you'd be down for this approche?

@bomilkar
Copy link
Contributor Author

bomilkar commented Apr 4, 2024

What exactly do you call "mean()" ?? Average, or Median, or what?

"w_y/w_x is tangens" is true, but only if w_y and w_x are orthogonal. Hence w_y and w_x must be vectors, not scalars.

Average is a scalar, not a vector.

Good luck if you want to deliver a proof. I'll be curious to see that.

@kobedegeest
Copy link
Contributor

kobedegeest commented Apr 4, 2024

ok here it is

for easy of notation let:
w be magnitude of wind
a magnitude of airspeed (unknown to gps or xcsoar though)
g magnitude of ground speed (provided from gps data)
alpha is heading (direction the plane points in also unknown to gps or xcsoar)
beta ground track angle (also provide from gps)

index x and y refer to the x component of the vectors x and y axes are global and fixed (for example say x is east and y is north)

We can say pilot trys to fly constant airspeed and heading to change at fixed rate we thus have 
For ease of proof lets define angles to be measuresd with respect to x-axes counter clockwise positive as is usual convention in math
ax = a*cos(alpha)
ay = a*sin(alpha)

vec(g) being the vector ground speed is the sum of vec(w) and vec(a) to sum vectors we can sum their components thus have 
gx = ax + wx = a*cos(alpha) + wx
gy = ay + wy = a*sin(alpha) + wy

we also get g and beta from gps and have that 
gx = g*cos(beta)
gy = g*sin(beta)

over full circle alpha goes from 0 to 360 thus since delta_alpha is cst average over time (or data points) is same as average over alpha (here avg is average over time or data points).
avg(gx) = avg(a*cos(alpha) + wx) = a*avg(cos(alpha)) + wx = wx
avg(gyx) = avg(a*sin(alpha) + wy) = a*avg(sin(alpha)) + wy = wy

a, wx and wy are constant and can thus be brought out of the avg().

and we have that vec(w) = [wx, wy] thus the wind vector is uniquely defined.

QED

note1: unlike alpha, beta does not have a constant change we can thus not say avg(gx) = g as average over time is not same as average over beta, and g is also not constant. So don't be confused by that.

note2: choice of angles to be wrt x-axes will not change result of proof just the expressions of the ax, ay, gx and gy etc

note3: assumption of cst airspeed and turn rate are the basis of any circling wind algorithm including the current implementation and your proposal, error on these assumptions is what we have been discussing.

I hope this clears up any confusion in notation as well. And that every step is clear. But happy to clarify further if needed.

@bomilkar
Copy link
Contributor Author

bomilkar commented Apr 4, 2024

a = g + w
thus g = a - w

You can multiply vec * scalar.
But you can't add vec + scalar
As stated before avg() is a scalar.

avg(acos(alpha) + wx) = aavg(cos(alpha)) + wx
this is wrong.
Same for y.

Your calculation of avg(gx) and avg(gyx) are wrong. It should be a * x.
And what is gxy ??

Sorry, but your proof doesn't hold.

@kobedegeest
Copy link
Contributor

a = g + w
thus g = a - w

this comes down to how one defines w i guess has no impact on proof just sub w for -w

But you can't add vec + scalar

add which point do i do this? cus i don't

And what is gxy ??

clearly a typo lets not get confused by that shall we.

It should be a * x.

What the hack is x now all of a sudden? You claim a typo invalidates my proof but i should know what a undefined variable x is?

Sorry, but your proof doesn't hold.

at what point please enlighten me cuz all you say is it is wrong without proving it is wrong and saying I add scalers to vectors?
whenever i used a vector a stated vec(..) and never did i write vec(..) + something.

Please don't dismiss an idea cuz your to stubborn to admit a mistake.

average of a sum is sum of the averages and average of a cst times something is that cst times the average of that something if you are gone claim that's wrong please proof so i would like to see that one.

@kobedegeest
Copy link
Contributor

And unlike you i have already tryed it and it does give the wind speed so maybe give it a try before holding on to your believes

I proved you wrong before about change in angle track being constant. So sorry if I have more faith in my math skills than yours.

@kobedegeest
Copy link
Contributor

and its plus for all i care if my airspeed is 100km/h and their is a wind of 20 km/h in the direction i am flying my ground speed will be 120km/h and not 80

a = g + w
thus g = a - w

and who is violating vector operations here? not me thats for sure cuz i defined a, g, and w as bein a scalar of magnitude of the vector. So what you say not even correct anyways.

@kobedegeest
Copy link
Contributor

And you say avg() is a scalar, if i want to take the average of a vector i take the average of a vector which will still be a vector. But besides the point of the proof as all variables used in the equations are scalars.

@bomilkar
Copy link
Contributor Author

bomilkar commented Apr 4, 2024

@kobedegeest You're calling me stubborn, etc. I don't take that!
I stop here.

@kobedegeest
Copy link
Contributor

etc? that's realy the only thing i called you. And in my defense you asked for a proof and dismissed it without disproving it. And don't stop a good idea just cuz one person challenged the idea.

Anyway here are some examples compairing current implemantation your proposal and what i just asked you if you had or would consider. Just to give an idea noise is +-5%.

image
image

and if you don't want any feedback fine but don't stop with something you want to do because some random guy. Was just trying to help.

@lordfolken
Copy link
Contributor

Why don't you guys plan a session together via some telefconferencing tool and start charting here on how the current or a new algorithm could be implemented.
I think the asynchronous format here, only leads to conflicts, while creating an idea together would work better.

@bomilkar
Copy link
Contributor Author

bomilkar commented Apr 4, 2024

Thanks @lordfolken for chiming in!
I agree that this platform is only partially suitable for facilitating a discussion like this. But I'm not convinced that a conference call would be better. I'm open to trying it, but ask for a certain amount of respect.
However, we should use github where it is good: when collaborating on actual code: If someone has an idea, they are welcome to implement it and put it up for discussion.

@bomilkar
Copy link
Contributor Author

@lordfolken : how can I proceed?

After extensive tests I’m now happy with the new CirclingWind algorithm. I have tested:

  1. the new algorithm on a testbench in various scenarios
  2. replay in XCSoar of 5 flights with a totals of ~25 hours flight time, in various scenarios comparing results vs the old CirclingWind algorithm and vs ZigZag, with and w/o TAS input

The bottom line:
The new CirclingWind algorithm behaves more robustly than the old one, which manifests itself in 2 symptoms: 1. there are fewer "erratic spikes" and 2. the result is less susceptible to the courses, regardless of whether it is a headwind or a tailwind. (2nd was one my main motivators to rewrite the algorithm. With the old wind calculation I noticed quite often that it underestimated the wind on a tailwind course and overestimated the wind on a headwind course. It didn’t always happen but it was very annoying.)
Although I haven’t seen any cases where it produced worse results, but of course this new algorithm needs a lot more testing. BTW: I’ve created my own build and will fly with it.
I have not replayed any flights with wind speed over 40 km/h. And I have only tested it with one particular GPS receiver (Power FLARM). These are probably the most important tests still missing.

If you are interested I would like to share and discuss the results of the replays in a teleconference. Just let me know.

On a related topic: running both Circling and ZigZag wind calculations simultaneously doesn’t seem to be a good idea. There is one gliding-average unit which is fed by either Circling or ZigZag depending on the flight phase. But these 2 calculations’ results are different since they measure different winds: one measures the speed at which the thermal moves the other measures the speed at which the air in between thermals is moving. The first is usually slower compared to the second. Hence, the gliding-average unit’s input switches back and forth between “different” winds. The result is a mix of “apples and oranges”. We should rethink this. For me I draw the conclusion that I will not switch on both Circling and ZigZag at the same time.

@bomilkar
Copy link
Contributor Author

Here is some insight related to the quality metrics discussion under #1391 . I copied some of the xcsoar.log files to a DropBox:
xcsoar.log files
In there you find LogString() outputs from src/Computer/ConditionMonitor/ConditionMonitorWind.cpp (after the commit bc40b2f ) which reflects what the pilot sees in the Wind InfoBox in 10 sec steps. (in addition to time and altitude)
Then there are the "CalcWind2 circling:" entries which come from the LogString() calls in src/Computer/Wind/CirclingWind2.cpp These give you some insight into the quality metrics and input into the gliding average unit.
The dates in the filenames reflect the date of the flight as posted on WeGlide.
The part "_noAS" stands for no TAS. I simply removed the lines containing airspeed data from the nmea file.
Hope that helps.
Let me know if you like to see anything else.
If you are interested in how the quality metrics work at higher wind speed, then look at 2023-07-17.

@kobedegeest
Copy link
Contributor

How come QM can be negative? So QM can be -1 and still give Q=5.
QM is the qualti metric that you check is below 0.6 i believe and than bin in Q: 0-5 but their it is assumed that QM is larger than 0
From reading the code i also believed this to be a positive number. left turn vs right turn maybe?

Also i would be happy to provide you with additional flight logs if that might help testing (no airspeed though).

@bomilkar
Copy link
Contributor Author

Good catch!
I misunderstood AbsoluteDegrees().
QM should never be negative.

@lordfolken
Copy link
Contributor

@lordfolken : how can I proceed?

First bring the PR up to shape. No codacy or other issues.

@bomilkar
Copy link
Contributor Author

bomilkar commented Apr 13, 2024

First bring the PR up to shape. No codacy or other issues.

Will do.
Can I run Codacy alone? I don't want to waste CPU-cycles on builds while I'm iterating on Codacy issues.

Can you give me a hint how to get rid of issues like "struct/class member 'xxx' is never used" in hpp files?
One example:
it complains about "class member 'CirclingWind2::active' is never used" in src/Computer/Wind/CirclingWind2.hpp
but it doesn't complain on CirclingWind::active
and I don't see the difference.

@bomilkar
Copy link
Contributor Author

No codacy or other issues.

There are 6 Codacy issues which I cannot address:

  1. „Clones added“. This one is expected. The intention of this draft PR is that eventually after some more validation the previous CirclingWind class be replaced. That'll resolve the clones issue.

  2. there are 5 issues of the type “struct member X is never used” or “class member X is never used”. I don’t know how to address this. There are some articles in the web which classify this as a “false positive”. But no hint how to address this. See also https://forum.xcsoar.org/viewtopic.php?t=4338

@lordfolken what next?

@GuillermoHazebrouck
Copy link

I am also interested in this topic, but for a different project. It is true that you don't need this algorithm when there is barometric altitude, heading and airspeed available. If you have all that, you can make the corrections for the TAS and just solve a triangle. If you lack one of those, then you still might want to fall back on circling however, and even try to enrich the algorithm with whatever data that is still available.
I have actually implemented an algorithm based on the rate of turn and tested it in a simulation a few times a few years a go, and it seemed to give fairly good results, even when the circles were not perfect. What I do is scanning the history to find a sample where the rate of turn remains within certain bounds. The rate of turn I use is computed by taking the average course variation. The course is the one given by the GNSS, falling back to an own computation if the GNSS is not providing it. Then, I use homologous points in the trajectory to evaluate the drift and take the mean value over the sample. The algorithm is quite simple and robust, but it has obviously some pitfalls.
Another solution I was considering is to collapse the cycloid to a circle using a collapsing function and minimum squares. This would in theory help to better understand the quality of the data and to dismiss trajectories that are not well suited.
Based on my general experience I can tell that increasing the complexity of a program will not necessarily yield better results. If you are dealing with noisy data and you don't understand the nature of the noise or the behavior of the system, then you will not be able to improve the predictions. In this case, the best is to evaluate the input and only to produce an output when the data looks clean. You could evaluate this statistically for example. An algorithm that throws a value without evaluating the quality of the data is an "arrogant" algorithm, and you should try to avoid that. It is better to return "I don't know" than "trust me, this is the answer".

@bomilkar
Copy link
Contributor Author

bomilkar commented May 6, 2024

@GuillermoHazebrouck thanks for the feedback!
Several of the implementation details you mention sound familiar and are in my current version of the circling algorithm.
I've been testing it in simulations (a lot) and in flight (not as much as I would like to): the results are good and address issues which I had with the old algorithm.
Taking the rate of turn from the GNSS is a bit obscure, for reasons that have been described above. Therefore I'm currently testing the algorithm with the rate of turn from a gyroscope. Surprisingly enough it makes nearly no difference for the testcases I have tried so far. These testcases have wind speeds up to 22 km/h and TAS around 100 km/h.
You mention a different project. Which one is this?
Are you interested in looking deeper into this project here?

@GuillermoHazebrouck
Copy link

GuillermoHazebrouck commented May 6, 2024

Take a look at my repositories for the other project... I won't mention it here or Max will censure my comments accusing me of unethical practices (yes, even when being also open source).
Indeed, everything coming from the GNSS can be a bit obscure, but on the other side, these things come usually with good filters. Unless you have the specifications of it, you will never know what they are doing. UBlox has very well documented spec sheets for example.
I don't know of the GNSS providing rate of turn directly however (unless it is using multiple sensors?). I only know it provides a calculated course. In my case I calculate the rate of turn manually from that. But the problem I see with doing this is that it will shift from the actual value the higher the wind, because there is a point in the cycloid where the ground curse tends to a discontinuity and the numerical algorithm will start having troubles much before than that... So I think using gyroscope data when available still makes sense! Maybe you can try to see what happens when the wind speed goes above 50% of the airspeed... (although that's not representative of an usual flight).
My tests were also at about the same speeds as you, but I was not using sensor data, just simulations from flight gear.
For me it is still important to be adaptive with the different sources of data, to be able to cope with different scenarios. You cannot trust sensors 100%.

@XCSoar XCSoar locked and limited conversation to collaborators May 21, 2024
@lordfolken lordfolken converted this issue into discussion #1449 May 21, 2024

This issue was moved to a discussion.

You can continue the conversation there. Go to discussion →

Labels
None yet
Projects
None yet
Development

No branches or pull requests

4 participants