You are here
Time of Arrival (ToA)
In the previous chapter, we've established that blade vibration causes the tips to move relative to the shaft's position. A deflected tip causes a shift in the proximity probe's waveform. The size of the shift is proportional to the deflection.
So, how do we calculate the size of the shift?
To determine this, we need the exact instant a blade moved past a probe. This instant is referred to as the Time of Arrival (ToA). Each ToA is the result of applying a trigger criterion to each pulse. Let's consider four trigger criteria.
Chapter outcomes
Understand that we need to implement a trigger criterion to determine the ToA of each pulse.
Understand the constant-threshold trigger method's principle.
Understand the detrimental effects of noise on the trigger criterion.
Understand how to implement hysteresis to prevent multiple trigger events.
Trigger criteria
Several ToA trigger criteria have been discussed in literature 1 2. Four of the most common ones are shown in Figure 1 below.
Here's the principle behind each criterion:
- The constant threshold trigger method: A ToA is triggered when the signal crosses a predefined threshold.
- The maximum slope method: A ToA is registered when the pulse's rate of change, or slope, is at its maximum.
- The constant-fraction crossing method: Similar to method A), but the threshold level is calculated as a percentage drop from the pulse's maximum voltage level.
- The maximum voltage method: The ToA occurs when the pulse is at its maximum.
We start our journey to convert raw BTT data into blade frequency, amplitude and phase by selecting one of these methods. Based on its simplicity, accuracy 2 and widespread adoption, we choose the constant threshold method.
Constant threshold trigger method
This method continuously compares a probe's signal to a voltage threshold. Each time the threshold is crossed, a ToA is recorded.
An idealized pulse generated by a blade as it passes a probe is shown in Figure 2. A trigger threshold of 0.4 V has been set. The moment the voltage signal crosses this level, at approximately 40 \(\mu\)s, the ToA is registered.
Before we get into the code, let's get you operational with Python and this tutorial's supplementary material.
Bladesight package
I have released a pip
installable package, called bladesight
, you can use to access this tutorial's datasets and functions.
The Pypi package can be found here: https://pypi.org/project/bladesight/. The Github repo can be found here: https://github.com/Bladesight/Bladesight.
Follow along with the worksheet
This BTT tutorial is code-centric.
The concepts are explained with reproducible code. I've collected the code examples for each chapter into its own Jupyter worksheet. The worksheets can be accessed at this Github repo:
https://github.com/Bladesight/bladesight-worksheets
Please resist the urge to skip the worksheets. The biggest stumbling block to learning BTT (and many subjects for that matter) is not the theory. It's to understand how the theory is implemented. Getting into the code will pay off in the long run.
You can execute the worksheets in two environments:
- Google Colab, and
- a local installation
Google Colab
Google Colab allows you to run Jupyter notebooks in the cloud. You can open the Google Colab notebook for this chapter by clicking on the link below:
You need a Google account to use Colab. If you don't have one, you can create one for free here.
Local Python installation
At some stage you'll want to set up your local environment to do your own development. I'll guide you through the process. If you can already set up a local Python environment, skip this section.
Python version
Here is an excellent (and entertaining) video on how to install Python on Windows. Please note: the bladesight
package only works on Python 3.9, 3.10 and 3.11. It does not yet work on versions 3.12 and 3.13. If you do not have a compatible version of Python installed, you can download it here.
Virtual environments
Virtual environments are an excellent way to isolate different Python projects from one another. I highly recommend setting one up.
On Windows, you can use these steps to set up a virtual environment:
-
Go to the directory where you want to create the virtual environment. In my case it is
"C:\Users\Dawie Diamond\intro_to_btt\"
. -
Open a command prompt in the directory. In Windows 11, this can be achieved in one of two ways:
-
Type
cmd
in the address bar of the directory:Press enter after typing
cmd
.A command prompt should open.
-
You can also right click in the directory and select "Open in Terminal":
-
-
Type the following command into the command prompt to create a virtual environment called
venv
: -
Run the below command to activate the virtual environment:
-
Execute the below instruction to install the bladesight package:
-
Install Jupyter:
-
We make extensive use of
plotly
for plotting: -
Now you can launch Jupyter:
Tip
One of the first people to test drive this tutorial told me he almost punched a hole through his screen when trying to set up the virtual environment.
I understand this. Setting up a virtual environment is often the most difficult part of starting with a new project.
Please reach out to me at dawie.diamond@bladesight.com and I'd be happy to help you get set up.
If you're on a different platform, here's an excellent video on how to set up virtual environments. The video covers different kinds of operating systems. I've created links for each operating system below:
Vectorized implementation of the constant threshold trigger method
Here's an artificial signal that contains three pulses:
Let's extract the ToAs.
How the code is presented
These are the tutorial's first code examples. I therefore repeatedly display the complete code. Different lines are highlighted and discussed each time.
Step 1: Load the probe signal
- This line loads the table into memory. It returns a
Pandas DataFrame
. We will be making extensive use of Pandas DataFrames throughout this tutorial. Its documentation can be found here: https://pandas.pydata.org/docs/
On Line 1, we import the bladesight
package. In addition to hosting the functions developed in this tutorial, the package makes it simple to download and open the datasets.
On lines 3 and 5, we download this chapter's dataset and read the three_generated_pulses
table. The first ten rows of the table are shown below:
time | data |
---|---|
0 | -0.000594587 |
0.0333444 | 0.00106139 |
0.0666889 | 2.60792e-05 |
0.100033 | -0.00232804 |
0.133378 | -0.00056633 |
0.166722 | 0.00040294 |
0.200067 | -0.000591657 |
0.233411 | -0.000962895 |
0.266756 | 0.000452657 |
0.3001 | -8.50727e-05 |
Units
Column | Units |
---|---|
time |
seconds |
data |
Volt |
The DataFrame has 2 columns: time
and data
. The time
column contains each data value's timestamp.
Step 2: Set the threshold direction and value
In Line 7, we specify the direction of the trigger. If TRIGGER_ON_RISING_EDGE
is True
, we trigger when the signal crosses the threshold on the rising edge. If TRIGGER_ON_RISING_EDGE
is False
, we trigger when the signal crosses the threshold on the falling edge.
In Line 8 we set the threshold level to 0.4 Volt.
Step 3: Determine when the signal has crossed the threshold
- We use the method
.astype(int)
at the end of this line because, by default, comparison operators such as>=
and<=
result in boolean values. We need an integer column for the steps that follow.
In lines 10 - 17, we determine when the signal is "over" the threshold level. The definition of "over" depends on the direction of the trigger. If we trigger on a rising edge, the signal is "over" when it is larger than the threshold. If we trigger on a falling edge, the signal is "over" when it is smaller than the threshold.
The variable sr_threshold_over
contains an array of ones and zeros that indicate whether the signal is above or below the threshold. This variable is shown on top of the original signal in Figure 4 below:
Step 4: Determine when the indicator has changed
- The
bfill
method is used to backward fill any missing values in a pandas Series. We do this because the first value ofdiff_sr_threshold
isNaN
. Why is itNaN
? Because the first value, by definition, has no prior value to be compared with. We usebfill
here to force the first value equal to the second value.
Since we are interested in the exact instant the signal crosses the threshold, we calculate the change in sr_threshold_over
. In Line 19, the diff
method calculates the consecutive differences between adjacent measurements. In Line 23, the >
operator is used to find where the indicator has changed from zero to one.
The result of this operation is shown in Figure 5 below:
Step 5: Select the ToAs
Line 25 is the culmination of our algorithm. We select only the time stamps where our over/under indicator has changed from zero to one. These time stamps are our ToAs.
These are the values in sr_toas
:
How can we set the ideal threshold level?
I don't have a one-size-fits-all answer. Here are some general guidelines I use:
-
The threshold level must intersect the pulse where its slope is steep. Do not set the threshold near a turning point. You can read about it in my paper 2.
-
The blades on your rotor have different lengths. The pulse shape depends on how far the tip is from the probe. Some pulses are noticeably shorter than others. The threshold level must be set to enable the maximum number of blades to be "seen".
-
Many BTT probes produce a first-derivative shaped signal. In other words, each pulse:
-
starts at zero
-
increases in the positive direction
-
rapidly dips through zero again
-
extends into the negative direction
-
returns back to zero
It is tempting to set the threshold level at zero, because you want to capture the ToA when the signal dips through zero in the middle. However, this is a bad idea. You are almost guaranteed of spurious ToAs because of noise around zero before and after the pulse settles.
If you do want to set the threshold level at zero in the middle of the pulse, you should implement hysteresis. We discuss hysteresis later in this chapter.
-
Sequential implementation of the constant threshold trigger method
The vectorized implementation shown above is great... but we can do better. We now implement a sequential version of the constant threshold trigger method.
The sequential implementation is the preferred one for two reasons:
- Contrary to how many people (...I'm referring to myself here) have been trained to think, it is faster than the vectorized implementation.
- It is simpler to understand.
How can a sequential, for-loop based approach be faster than the vectorized method? Python is frequently regarded as a slow language. This is false. It is only slow if you use it like nobody ever intended for it to be used.
Many approaches exist to enhance its speed. One can often achieve performance levels comparable to compiled languages.
We use Numba. Numba is a powerful Just-In-Time (JIT) compiler for Python, compiling portions of your code into machine code at runtime. This often leads to blistering fast functions. Numba is a dependency of the bladesight package, so its already installed if you've followed the installation instructions at the start of this chapter.
To use Numba, we import the njit
decorator from the package, and wrap our function in it.
Simple example
The simplest implementation of a sequential algorithm involves a for
loop to monitor when the signal passes a constant threshold.
- The array of time values. Must be a Numpy array.
- The array of voltage values. Must be a Numpy array.
- The threshold value. Must be a float.
- The estimated number of ToAs in this signal. Defaults to None. This number is used to pre-allocate the array of ToAs. If this number is not provided, the array will be pre-allocated as the same dimension as arr_t and arr_s. You should specify this value for large signals.
- Whether to trigger ToAs on the rising or falling edge. Defaults to True. If True, the ToA is triggered on the rising edge.
- This type annotation
-> np.ndarray
indicates the return value of the function is a Numpy array. - We pre-allocate the array of ToAs. This is a performance optimization. If we don't pre-allocate the array, the function will have to resize the array each time a ToA is found. It's easy to estimate the expected number of ToAs in the signal, and should be done for any real signal.
- Initialize the number of ToAs found in this signal.
i_toa
will increase by one each time a ToA has been found. - The sequential approach compares each sample to the previous sample. Here, we initialize the
prev_sample
value to the first value in the array. - We loop through the entire signal and compare each sample to its predecessor.
- Here, we check if the threshold has been crossed. If so, we store the ToA in the
arr_toa
array. This if-clause is executed for a rising edge trigger. - If the threshold has been crossed, we store the ToA in the
arr_toa
array. - Increment the
i_toa
value by one. - The check for a falling edge trigger.
- We're done with this sample. This command prepares us for the comparison in the next loop pass.
- Return only the ToAs that have been registered. The rest of the array is filled with -1 values.
Type annotations in Python
Throughout this tutorial, you'll notice in many of the functions, the function arguments are accompanied by type annotations. In the first argument above, arr_t : np.ndarray
, arr_t
is the variable name and : np.ndarray
is the type annotation . These type annotations specify the expected data types of the arguments.
It's important to note Python itself does not enforce these type annotations. They are primarily utilized by type checkers and serve as a helpful guide for users of the function.
We often need to import annotations from the typing
library. In the example above, we need to import the Optional
type with the statement from typing import Optional
.
If we pass our three pulses signal through this function, we get the same result as the vectorized implementation:
>>> toas = seq_simple_threshold_crossing(
df_proximity_probe['time'].values, #(1)!
df_proximity_probe['data'].values,
threshold=0.4
)
>>> print(toas)
[27.44248083 52.45081694 77.45915305]
- Note how we use the
.values
attribute of the Pandas object. This returns a Numpy array, which is required by the function.
Interpolate on voltage
Our signals are sampled by a data acquisition system (DAQ). Suppose the DAQ's sampling rate is 100 kHz. This means one sample is acquired every 10 \(\mu\)s. What happens if the threshold is crossed between two samples? The methods derived up to this point will miss the exact ToA. The ToA will only be triggered at the first available sample after the threshold has been crossed.
If we have the raw analogue signal stored on disk, it's wise to capitalize on the continuous nature of the signal.
To improve the accuracy of ToA determination, we employ interpolation between the two nearest samples when the threshold is crossed.
The below function adds interpolation to the trigger algorithm. It is a modified version of the previous function. The only difference is the interpolation step, marked with a comment (+)
:
- This part performs linear interpolation between the two samples either side of the threshold. It is the only difference between this function and the previous one.
Why use linear interpolation?
Linear interpolation is the simplest form of interpolation. It is also the fastest. In my experience, the accuracy gain from a higher order interpolation method is marginal. It's not worth the additional code complexity.
Let's use this function on the same signal:
>>> toas = seq_threshold_crossing_interp(
df_proximity_probe['time'].values,
df_proximity_probe['data'].values,
threshold=0.4
)
>>> print(toas)
[27.42940329 52.44126697 77.45548277]
These values are slightly different than the previous ones ( the previous ones are [27.44248083 52.45081694 77.45915305]
). The new values are more accurate than the old ones.
Hysteresis
The signals we've used so far contain minimal noise. In the absence of noise, the ToAs are easily extracted. However, in real-world scenarios, signals contain non-negligible noise components.
Sometimes, excessive noise is consistently present throughout the signal. Other times it appears in short, sporadic bursts.
To illustrate this point vividly, I've generated a noisy signal that contains three pulses, shown in Figure 6 below:
Zoom
If there is a "Reset zoom" button on the bottom of the figure, you can zoom and pan on the plot. Drag across the screen to zoom , hold ctrl
and drag across the screen to pan.
If you zoom into each pulse around 0.4 V, you'll observe the signal crosses the threshold multiple times. Let's apply our trigger algorithm to this signal. We first need to load the noisy signal dataframe:
Having applied our algorithm to the noisy signal, these are the ToAs:
[27.20280579, 27.26724046, 27.45521975, 27.54516822, 32.4110573, 32.48523812,
32.59630977, 32.73320484, 52.10882236, 52.46204252, 52.59380688, 57.42079388,
57.64309888, 57.69638523, 57.78086048, 57.87833845, 77.46723661, 77.66015932,
78.09842428, 82.3623517, 82.50989437, 82.74746473]
Our algorithm incorrectly returns many ToAs around the 0.4 V level.
We can implement hysteresis to solve the issue. Hysteresis introduces a second threshold and a state variable to our algorithm. The algorithm therefore uses a lower and an upper threshold. The state can only change based on certain rules. The rules are:
- If the state is low and the current sample is above the upper threshold, the state becomes high. A ToA is triggered when the state moves from low to high.
- If the state is high and the current sample falls below the lower threshold, the state becomes low.
- Otherwise, the state is maintained.
Figure 7 below illustrates the state variable and ToA trigger events for a (absurdly) noisy signal.
The below algorithm implements hysteresis to determine the ToAs on the rising edge:
- The height of the hysteresis. It has the same units as the signal. This value is used to calculate the lower threshold.
- The trigger state is a boolean value indicating whether the trigger is currently high or low. We initialize it to True if the first sample is above the lower threshold, and False otherwise.
- If the trigger state is True, we check if the current sample is below the lower threshold. If it is, we set the trigger state to False.
- If the trigger state is False, we check if the current sample is above the upper threshold. If it is, we set the trigger state to True and calculate the ToA using interpolation.
We can use this function to extract the ToAs from the noisy signal:
>>> toas = seq_threshold_crossing_hysteresis_rising(
df_proximity_probe_noisy['time'].values,
df_proximity_probe_noisy['data'].values,
threshold=0.4,
hysteresis_height=0.2
)
>>> print(taos)
[27.20280579 52.10882236 77.46723661]
How to select hysteresis height
and n_est
-
The hysteresis height should be at least twice as large as the noise. A sufficient estimate can usually be made by observing the largest noise spikes.
-
The
n_est
value should be larger than the number of pulses in the signal. I use this calculation:
The *2
at the end of the equation doubles the pre-allocated array length. We can therefore be sure there's enough space for all the ToAs.
Performance of sequential vs vectorized implementations
We've now implemented vectorized and sequential algorithms to extract the ToAs. Which one is faster? Let's use both approaches on a real signal.
The dataset used here contains a signal acquired from an eddy current probe. It is stored in the table table/aluminium_blisk_1200_rpm
. The signal was acquired for a five blade rotor running at 1200 RPM. The acquisition was performed at a sampling rate of 2 MHz. The signal is shown in Figure 8 below:
We can use the Jupyter %%timeit
magic method to time a piece of code. Here's the result:
Vectorized | Sequential |
---|---|
301 ms | 18 ms |
The sequential approach is almost 20 times faster than the vectorized method. The sequential algorithm also contains hysteresis and interpolation that the vectorized function does not.
Real Time Hardware
This chapter assumes you can store large analogue signals on your disk. If your rotor is small and turns slowly, this method is suitable.
When your blades are large or your rotor spins fast, you'll need to determine the ToAs in real-time. This typically requires custom hardware. Here are some providers of BTT hardware:
Conclusion
Congratulations, you've completed the second chapter of this tutorial. You're approximately 30% of the way through. You've now learned how ToA extraction works.
In my experience, getting your Python version set up is a massive challenge. From here, we can churn out code to our heart's content.
Lets push on...
Chapter outcomes
Understand that we need to implement a trigger criterion to determine the ToA of each pulse.
Understand the constant-threshold trigger method's principle.
Understand the detrimental effects of noise on the trigger criterion.
Understand how to implement hysteresis to prevent multiple trigger events.
Acknowledgements
Thanks to Justin Smith and Alex Brocco for reviewing this chapter and providing feedback.
Coding exercises
Here are two coding exercises to solidify your understanding of the material covered in this chapter. Please attempt them before you expand the solution bar.
Problem 1: Automatic range detection ![:yawning_face: 🥱](https://cdn.jsdelivr.net/gh/jdecked/twemoji@15.0.3/assets/svg/1f971.svg)
So far, the threshold level has been set manually from eyeballing the signal. This is often fine, but we can do better...
Write a new function, called determine_threshold_level
, that receives the signal and the % of the range to use as the threshold level. The function should return the threshold level.
Reveal answer (Please try it yourself before peeking)
Problem 2: Hysteresis on the falling edge ![:neutral_face: 😐](https://cdn.jsdelivr.net/gh/jdecked/twemoji@15.0.3/assets/svg/1f610.svg)
We've showcased how a trigger criterion with hysteresis works on the rising edge. Implement a trigger criterion that uses hysteresis on the falling edge.
Write a new function, called seq_threshold_crossing_hysteresis_falling
, that extracts the ToAs on the falling edge using hysteresis.
Reveal answer (Please try it yourself before peeking)
-
Zimmer, A.K., 2008. Investigation of the impact of turbine blade geometry on near-field microwave blade tip time of arrival measurements (PhD thesis). ↩
-
Diamond, D., Heyns, P.S., Oberholster, A., 2021. Constant speed tip deflection determination using the instantaneous phase of blade tip timing data. Mechanical Systems and Signal Processing 150, 107151. ↩↩↩