Analiza charakterystyki danych akcelerometru i zaprojektowanie filtra

Mam dane z przyspieszeniomierza o wartości około 32 sekund w podstawowym scenariuszu jazdy 25MPH normalnych dróg wraz z uderzeniem około 7 dziur i nierównego odcinka drogi. Akcelerometr jest zamontowany na desce rozdzielczej mojego samochodu za pomocą dwustronnej taśmy.

Problem: Mam wszystkie dane, które są zaszumione z akcelerometru, i muszę w prosty sposób wykryć, że wystąpiło zdarzenie w wyboju. Poniżej znajduje się kilka wykresów danych w dziedzinie czasu i FFT. Akcelerometr mierzy w GForce

Zasadniczo chcę, aby moje arduino wiedziało, że dziura pojawiła się z dość dużą dokładnością i nie używa matematyki i technik na poziomie magisterskim.

Akcelerometr próbkowany przy częstotliwości 100 Hz ma prosty filtr o niskiej częstotliwości 50 Hz RC NA OSI Z

Here is the CSV data for the 32 seconds of accelerometer readings TIME, GFORCE format:

http://hamiltoncomputer.us/50HZLPFDATA.CSV

AKTUALIZACJA: Jest to pełna przepustowość RAW akcelerometru próbkowanego 1000 Hz z najwyższą częstotliwością próbkowania, jaką mogłem uzyskać na Arduino. Bezpośrednie pobieranie pliku CSV: około 112 sekund danych

http://hamiltoncomputer.us/RAWUNFILTEREDFULLBANDWIDTH500HZ.csv

Czarny ślad jest niefiltrowanym RAW Dane akcelerometru: Niebieski ślad jest filtrowany przez filtr pasmowy na podstawie ekstremalnych częstotliwości występujących w FFT, Dominate 2HZ i 12HZ.

Wydarzenie wyboju wygląda tak w domenie czasu:

nie jestem pewien, jaki składnik 10 do 15 HZ znajduje się w FFT, czy jest to rzeczywisty otwór, czy jest to podskakiwanie kół do drogi, czy też jest to częstotliwość rezonansowa samochodu?

FFT:

wygląda na to, że są to rzeczywiste zdarzenia wybojowe, tutaj jest HPF @ 13HZ Dominujące cechy dziur wydają się ulepszone

Chcę być w stanie wykrywać i liczyć dziury w czasie rzeczywistym

Wydaje się to sprzeczne z intuicją, zawieszenie powinno poruszać się o wiele wolniej niż 10 do 13 Hz, co spowodowałoby chorobę lokomocyjną

AKTUALIZACJA:

Zgodnie z sugestiami AngryEE wykorzystałem pełną szerokość pasma akcelerometru 1000 Hz i maksymalną częstotliwość próbkowania, jaką mogłem uzyskać na arduino.

FFT:

oto przykładowy fragment danych o zdarzeniu w wyboju oraz niektóre nierówności i hałas drogowy:

Dodano obwód detektora obwiedni diody, wyjście wygląda tak samo ... Akcelerometr zawsze ma moc wyjściową od 0 do 3,3 V, a nie ujemną ...

AKTUALIZACJA:

Z wielu testów drogowych nigdy nie przekroczyłem wartości 1.6G do 45 MPH w moim samochodzie na osi Z, użyłem rand () do wygenerowania pseudolosowego przyspieszenia Gforce.

Mój pomysł polega na tym, że jeśli mogę spojrzeć na okna danych z 1 do 3 sekund, mogę obliczyć przesunięcie osi Z, ale martwiłem się o dryft akcelerometru i błędy w integracji. Nie muszę tu być nawet 90% dokładny,> 70% byłoby miło, ale jeśli patrzę na przesunięcie od jednej do trzech sekund na raz, czy byłoby to możliwe w czasie rzeczywistym? W ten sposób widzę, czy przemieszczenie jest większe niż jak 1 cal, 2 cale, 5 cali. Im większe przesunięcie, tym bardziej nierówny był guz lub dziura:

Czy możesz sprawdzić, czy robię to dobrze, po prostu skonfigurowałem na pulpicie, używając rand () do generowania losowego przyspieszenia od -1,6 do 1,6 G, rejestrując 3 sekundy danych przy symulowanej częstotliwości próbkowania 50 Hz

Jeśli chcesz uruchomić * nix, używam Sleep () z Windows.h, aby opóźnić 20 ms, częstotliwość próbkowania 50 Hz

Chciałem tylko sprawdzić, czy kod wygląda dobrze dla ciebie, nie zrobiłem jeszcze bufora, jestem trochę zdezorientowany, jak go zaimplementować: skomentowany kod pochodzi z klasy, nad którą pracuję , ale jeszcze tego nie rozumiem w 100%. Okrągły bufor pozwoliłby na ciągłe przesuwanie okien danych, prawda?

#include <cstdlib>
#include <iostream>
#include <iomanip>
#include <ctime> // USED BY RAND
#include <windows.h> // Used for delay


using namespace std;



#define SAMPLE_RATE   0.020 // Sample rate in Milliseconds
#define GRAVITYFT_SEC 32 // Gravity velocity 32 feet/sec
#define INCH_FOOT     12 // 12 inches in foot, from velocity to inch displacement calculation










int main(int argc, char *argv[])
{
    srand((unsigned)time(0)); // SEED RAND() for simulation of Geforce Readings

    // SIMULATING ACCELERATION READINGS INTO A CIRCULAR BUFFER

   // circular_buffer Acceleration; // Create a new Circular buffer for Acceleration

   // cb_init(&Acceleration, 150, 4); // Sampling @ 50HZ, 3 seconds of data = 150, size is float data of 4 bytes

    //Simulate a sample run of Acceleration data using Rand()

    // WE WILL BE SIMULATING "RANDOM" GEFORCE RATINGS using the rand() function constraining to -1.6 to 1.6 GFORCE 
    // These ratings are consistent with our road tests of apparently random vibration and Geforce readings not exceeding about 1.6 G's

    float Gforce[150]; // Random Geforce for 3 second window of data
    float velocity[150]; // Hold velocity information
    float displacement[150]; // Hold Displacement information


    float LO = -1.6; // Low GForce limit recorded from 6 road tests at different speeds
    float HI = 1.6; // High GForce limit recorded from 6 road tests at different speeds

    for(int i = 0; i < 150; i++) // 3 Second iwndow of random acceleration data
    {  
            Gforce[i] = LO + (float)rand()/((float)RAND_MAX/(HI-LO)); // Borrowed from Stackexchange : http://stackoverflow.com/questions/686353/c-random-float
            if( i == 0) // Initial values @ first Acceleration
            {
                velocity[i] = Gforce[i] * SAMPLE_RATE * GRAVITYFT_SEC; // Initial velocity
                displacement[i] = velocity[i] * SAMPLE_RATE * INCH_FOOT; // Initial Displacement
            }
            else
            {
                velocity[i] = velocity[i-1] + (Gforce[i] * SAMPLE_RATE * GRAVITYFT_SEC); // Calculate running velocity into buffer
                displacement[i] = displacement[i-1] +(velocity[i] * SAMPLE_RATE * INCH_FOOT); // Calculate running displacement into buffer
            }
            //cout << endl << Gforce[i]; // Debugging
            //cb_push_back(&Acceleration, &Gforce[i]);                   // Push the GeForce into the circular buffer


            Sleep(SAMPLE_RATE*1000); // 20mS delay simulates 50HZ sampling rate Sleep() expects number in mS already so * 1000

    }
    // PRINT RESULTS
    for (int j = 0; j < 150; j++)
            {
                cout << setprecision (3) << Gforce[j] << "\t\t" << velocity[j] << "\t\t" << displacement[j] << endl;
            }

    // READ THE BUFFER





    //cb_free(&Acceleration); // Pervent Memory leaks

    system("PAUSE");
    return EXIT_SUCCESS;
}

Przykładowy przebieg:

    GFORCE          FT/SEC          Inch Displacement Z axis

-0.882          -0.565          -0.136
0.199           -0.437          -0.24
-1.32           -1.29           -0.549
0.928           -0.691          -0.715
0.6             -0.307          -0.788
1.47            0.635           -0.636
0.849           1.18            -0.353
-0.247          1.02            -0.108
1.29            1.85            0.335
0.298           2.04            0.824
-1.04           1.37            1.15
1.1             2.08            1.65
1.52            3.05            2.38
0.078           3.1             3.12
-0.0125         3.09            3.87
1.24            3.88            4.8
0.845           4.42            5.86
0.25            4.58            6.96
0.0463          4.61            8.06
1.37            5.49            9.38
-0.15           5.39            10.7
0.947           6               12.1
1.18            6.75            13.7
-0.791          6.25            15.2
-1.43           5.33            16.5
-1.58           4.32            17.5
1.52            5.29            18.8
-0.208          5.16            20.1
1.36            6.03            21.5
-0.294          5.84            22.9
1.22            6.62            24.5
1.14            7.35            26.3
1.01            8               28.2
0.284           8.18            30.1
1.18            8.93            32.3
-1.43           8.02            34.2
-0.167          7.91            36.1
1.14            8.64            38.2
-1.4            7.74            40
-1.49           6.79            41.7
-0.926          6.2             43.2
-0.575          5.83            44.6
0.978           6.46            46.1
-0.909          5.87            47.5
1.46            6.81            49.2
0.353           7.04            50.8
-1.12           6.32            52.4
-1.12           5.6             53.7
-0.141          5.51            55
0.463           5.8             56.4
-1.1            5.1             57.6
0.591           5.48            59
0.0912          5.54            60.3
-0.47           5.23            61.5
-0.437          4.96            62.7
0.734           5.42            64
-0.343          5.21            65.3
0.836           5.74            66.7
-1.11           5.03            67.9
-0.771          4.54            69
-0.783          4.04            69.9
-0.501          3.72            70.8
-0.569          3.35            71.6
0.765           3.84            72.5
0.568           4.21            73.5
-1.45           3.28            74.3
0.391           3.53            75.2
0.339           3.75            76.1
0.797           4.26            77.1
1.3             5.09            78.3
0.237           5.24            79.6
1.52            6.21            81.1
0.314           6.41            82.6
0.369           6.65            84.2
-0.598          6.26            85.7
-0.905          5.68            87.1
-0.732          5.22            88.3
-1.47           4.27            89.4
0.828           4.8             90.5
0.261           4.97            91.7
0.0473          5               92.9
1.53            5.98            94.3
1.24            6.77            96
-0.0228         6.76            97.6
-0.0453         6.73            99.2
-1.07           6.04            101
-0.345          5.82            102
0.652           6.24            104
1.37            7.12            105
1.15            7.85            107
0.0238          7.87            109
1.43            8.79            111
1.08            9.48            113
1.53            10.5            116
-0.709          10              118
-0.811          9.48            121
-1.06           8.8             123
-1.22           8.02            125
-1.4            7.13            126
0.129           7.21            128
0.199           7.34            130
-0.182          7.22            132
0.135           7.31            133
0.885           7.87            135
0.678           8.31            137
0.922           8.9             139
-1.54           7.91            141
-1.16           7.16            143
-0.632          6.76            145
1.3             7.59            146
-0.67           7.16            148
0.124           7.24            150
-1.19           6.48            151
-0.728          6.01            153
1.22            6.79            154
-1.33           5.94            156
-0.402          5.69            157
-0.532          5.35            159
1.27            6.16            160
0.323           6.37            162
0.428           6.64            163
0.414           6.91            165
-0.614          6.51            166
1.37            7.39            168
0.449           7.68            170
0.55            8.03            172
1.33            8.88            174
-1.2            8.11            176
-0.641          7.7             178
-1.59           6.69            179
1.02            7.34            181
-0.86           6.79            183
-1.55           5.79            184
-0.515          5.46            186
0.352           5.69            187
0.824           6.22            188
1.14            6.94            190
-1.03           6.29            192
-1.13           5.56            193
0.139           5.65            194
0.293           5.84            196
1.08            6.53            197
-1.23           5.75            199
-1.1            5.04            200
-1.17           4.29            201
-0.8            3.78            202
-0.905          3.2             203
-0.0769         3.15            203
-0.323          2.95            204
-0.0186         2.93            205
Press any key to continue . . .

Wygląda na to, że można to rozwiązać przez dość proste filtrowanie. Oto twoje oryginalne dane:

To za dużo, aby zobaczyć, co dzieje się w indywidualnym wydarzeniu na odpowiednim poziomie szczegółowości. Oto tylko dane z drugiego 26 do 28:

Początkowo myślałem, że to filtr dolnoprzepustowy, ale to nie działa, ponieważ nie ma tam sygnału niskiej częstotliwości. Zamiast tego rośnie amplituda sygnału wysokiej częstotliwości. Oto dolnoprzepustowy nałożony na oryginał:

Zauważ, że wynika to z „średniej” sygnału całkiem dobrze, a nie podczas zdarzenia wyboju. Jeśli odejmiemy tę średnią od oryginalnego sygnału, pozostanie nam znacznie więcej odchyleń od tej średniej podczas wydarzenia niż w innym przypadku. Innymi słowy, tak naprawdę chcemy filtra górnoprzepustowego. Zrobimy to, odejmując dolnoprzepustowy od oryginału, ponieważ tak się tu dostaliśmy, ale w systemie produkcyjnym zrobiłbyś to przez jawne filtrowanie górnoprzepustowe. Tak czy inaczej, oto oryginał z filtrem górnoprzepustowym:

Wskazuje to teraz na oczywiste podejście do wykrywania zdarzenia. Podczas imprezy występuje znacznie większa amplituda sygnału niż w innym przypadku. Możemy to wykryć, obliczając RMS i stosując filtrowanie dolnoprzepustowe:

Powiększając wszystkie dane, widzimy:

To wyraźnie identyfikuje pięć zdarzeń w danych, chociaż nie wiem, czy to właśnie te dane mają pokazać. Przyglądając się bliżej wydarzeniom, zauważasz, że każde z nich ma niskie spadki około 1 sekundy przed i po szczytach. Oznacza to, że można zrobić więcej, jeśli próg sygnału RMS, jaki jest teraz, nie jest wystarczający. Na przykład prosty algorytm, który szukał wysokości punktu względem najniższej w ciągu 1 sekundy w obie strony, powinien dodatkowo zmniejszyć szum tła. Innym sposobem na powiedzenie o tym samym jest rozróżnienie tego sygnału w poszukiwaniu wzrostu w ciągu 1 sekundy. Zdarzenie w wyboju zostanie wówczas wykryte przez dublet, co oznacza wysoki pik, po którym następuje niski pik.

Innym sposobem spojrzenia na to jest pasmowe przekazywanie sygnału RMS. Jest już filtrowane dolnoprzepustowo, ale ponieważ szukasz nagłych zdarzeń z silnymi zboczami, odsuwanie niektórych niskich częstotliwości powinno działać również w celu zmniejszenia szumu tła.

Stąd jest wiele sposobów na udoskonalenie sygnału, ale mam nadzieję, że pokazałem, jak uzyskać użyteczny wynik co najmniej pierwszego przejścia.

Dodany:

Byłem ciekawy, jak dobrze będzie wyglądało zanurzenie po obu stronach szczytu, więc spróbowałem. Użyłem filtru nieliniowego, zaczynając od RMS z poprzedniego wykresu. Wartość każdego punktu jest minimalną wartością powyżej najniższego punktu w poprzedniej sekundzie i najniższego punktu w następnej sekundzie. Wynik wygląda całkiem dobrze:

Najniższy z 5 pików jest ponad 3 razy wyższy niż najwyższy szum tła. To oczywiście zakładając, że te 5 nierówności reprezentuje zdarzenia, które chcesz wykryć, a reszta nie.

Dodano w odpowiedzi na komentarze:

Zrobiłem filtry w dziedzinie czasu, więc nie znam bezpośrednio odpowiedzi częstotliwościowej. Dla filtra dolnoprzepustowego splotłem sygnał wejściowy z jądrem filtra COS ^ 2. Jeśli dobrze pamiętam, promień (odległość od środka do krawędzi) jądra wynosi kilka 100 ms. Eksperymentowałem z tą wartością, dopóki fabuła nie wyglądała dobrze. Do filtrowania dolnoprzepustowego RMS użyłem tego samego jądra filtra, ale tym razem o promieniu około sekundy. Nie pamiętam dokładnie. Eksperymentuj, aż uzyskasz dobre wyniki.

Filtr nieliniowy nie wykrył dubletów. Jak powiedziałem, znalazłem różnicę między bieżącym punktem a najniższym ze wszystkich punktów w ciągu 1 sekundy wcześniej, a także różnicę między bieżącym punktem a najniższym ze wszystkich punktów w ciągu 1 sekundy po. Potem wziąłem min z tych dwóch.

Oprogramowanie, którego użyłem, było programem, który zhackowałem w tym celu. Miałem już różne procedury do odczytu i zapisu plików CSV, więc wszystko, co musiałem napisać, to kod filtrujący, co jest bardzo proste. Reszta została wykonana przy użyciu istniejących programów do manipulowania i drukowania plików CSV.

Dziury wykrywające krawędzie mogą powodować problemy. Odpowiedzią jest obwiednia drgań samochodów, ponieważ rzeczywiste drgania widoczne przez czujnik mają znacznie wyższe częstotliwości. Poszedłbym z RMS do prądu stałego, który reaguje na około 15 Hz lub więcej, a niskie przekazywanie rzeczy.

Zamiast szukać filtru domeny częstotliwości lub progu, zalecam próbę znalezienia jądra dla „typowej” dziury i wykonanie z nim bieżącej korelacji. Byłoby to uważane za technikę dopasowywania szablonów i wydawałoby się, że nadaje się na platformę mikrokontrolera.

Szybki przegląd można znaleźć na stronie http://scribblethink.org/Work/nvisionInterface/vi95_lewis.pdf , a być może DOBBS, STEVEN E., NEIL M. SCHMITT i HALUK S. OZEMEK. „Wykrywanie QRS przez dopasowanie szablonu przy użyciu korelacji w czasie rzeczywistym na mikrokomputerze.” Journal of engineering of 9.3 (1984): 197-212.

Jeśli byłeś na grubszej platformie, poleciłbym obrócić falki.

Innym podejściem byłoby obliczenie ruchomej wariancji sygnału, aby sprawdzić, czy dziury naprawdę wystają. Oto funkcja matlaba dla filtra ruchomej wariancji, szerokość N punktów - sprytnie (jeśli sam muszę to powiedzieć) przy użyciu splotu do obliczeń

function y=movingvar(X,N)
% y=movingvar(X,N)
% Calculates N-point moving variance of  Vector X
% Highly recommend that N be odd (no error checking)
% Note: first and last N/2 points will be unreliable.
% Output will be a column vector.


X=X(:);
XSQR=X.*X;
convsig=ones(1,N);
y=(conv(convsig,XSQR)-(conv(convsig,X).^2)/N)/(N-1);

y=y(ceil(N/2):length(X)+floor(N/2));

Początkowo sądzę, że filtr dolnoprzepustowy może być niewłaściwym typem filtra. Dziura jest zasadniczo zdarzeniem o wysokiej częstotliwości - takim jak funkcja kroku lub fala prostokątna. Samo spojrzenie na filtrowane dane 50 Hz sprawia, że ​​myślę, że tracisz informacje o wyboju - wszystko wygląda jak te same skręty bez wyraźnego rozróżnienia dla zdarzenia wyboju. Najpierw użyłbym filtra górnoprzepustowego, a następnie filtra dolnoprzepustowego o znacznie wyższej częstotliwości. Możesz całkowicie uniknąć filtrowania dolnoprzepustowego, jeśli twój akcelerometr jest już filtrowany w dolnoprzepustowym.

Po uzyskaniu filtrowanych danych górnoprzepustowych uważam, że prosty komparator z ustawionym progiem odpowiednio wykryje szczyty w danych przyspieszenia spowodowanych przez dziury i pozwoli je policzyć.