A statistical ltering method for denoising of micro-force measurements

Precise measurement of mechanical forces is crucial to e cient micro-manufacturing. The quality of such measurements depends heavily on the properties of the noise inevitably accompanying every measurement process. In the micro-range the signal-to-noise ratio tends to be very low, and the noise dynamic varies for di erent frequencies. In result, common denoising methods that assume white noise perform poorly in this setting. In this paper a novel, easily implementable denoising method based on a local statistic of the measured data's spectrum is proposed. By testing it on a representative dataset, it is shown that the proposed method is robust and stable. Particularly, it allows for an e cient retrieval of the force signal encountered in micro-milling processes.


Introduction
In recent years, ecient and exible production of small components has become indispensable to high-precision mechatronics. This increased demand in miniaturization had an immense impact on the manufacturing technology and led to the development of many novel and innovative processing techniques for high-precision micro-parts (cf. e.g. [39]).
Micro-milling is a prominent method for rapid and cheap production of micro-components. Together with other material removal processes it has become one of the most promising methods for the production of micro-forming tools as well as micro-components (cf. e.g. [12,39]). This is especially true for the manufacturing of components with high surface quality and with textured or structured surfaces (cf. [3,37]). For end products where the desired functionality is determined by these surface properties micro-milling is a crucial micro-production method. Accordingly, in order to satisfy the desired requirements and properties of the nal product optimal process parameters have to be achieved and maintained within the manufacturing process (cf. [8,30,35,37]).
To determine optimal process parameters ecient models and simulations are needed. The predictive power of such models depends heavily on the modeling of the cutting forces with respect to the process parameters and the desired structures. Furthermore, in order to maintain optimal parameters and thereby in order to avoid unsatisfactory manufacturing results careful monitoring of the process is indispensable. In this case, also, the analysis of micro-forces plays an important role in the detection of aberrations from the desired process parameters.
Altogether, with the focus on surface texture and tool wear it is therefore essential to model and predict micro-forces involved in micro-milling precisely. Consequently, a thorough knowledge of these forces is necessary to guarantee the high quality of the manufactured micro-components.
Usually, the force models in cutting processes are based on the common assumption that the force is proportional to the cross sectional area [2,40]. However, the experimental data indicates that this approach is not well-equipped to fully describe eects occurring in micro-processes. Further, conventional methods and models cannot be scaled eectively to the micro-scale, due to the occurring size eects [38]. Consequently, the development of new force-models for micro-processes is of uttermost importance for manufacturing technology.
The detailed research on turning processes and model development for complex surface generation is presented in [25]. In particular, there, the authors consider a dynamic cutting force model, a regenerative vibration model, a machining system response model and a tool prole model.
Further, a two-degrees-of-freedom model which includes dynamical changes caused by forces acting during the micro-process was presented and discussed for turning processes in [11,27]. This model was successfully further developed and extended to the micro-endball-milling processes in [28]. In particular, the extended, more general force model is based on the tool-workpieceinteraction.
An analytical force model for micro-end-milling including tool run-out and tool wear eects was investigated in a series of papers by Bao et al.,see [46]. In addition, Li et al. [24] proposed a model of a three dimensional cutting force for micro-end-milling operation which includes a new nominal uncut chip thickness algorithm. This model was developed under the consideration of tool run-out and an exact trochoidal trajectory of the tool tip. Similar investigation together with a nite element model for the orthogonal cutting is presented in [1]. Finally, additional investigations of micro-milling forces and their models were carried out in [9,20] and [21].
However, many questions concerning force models for micro-cutting processes are still open. In order to develop new, improved models as well as to qualitatively verify the capability of already existing models to precisely predict the micro-forces encountered in such micro-processes a reliable recovery of these forces from noisy measurements is indispensable.
As mentioned above, unanticipated manufacturing errors of various kinds may occur during micro-processes. In order to correctly and reliably assess the state of such processes monitoring is necessary. However, monitoring methods for miniaturized manufacturing processes are more complex and challenging than their counterparts in the conventional manufacturing. A vast number of methods has been proposed for the monitoring of microprocesses. The most common among them are the measurement of micro-forces; further, the measurement of audio signals generated by the process; and, nally, the recording and analysis of the process by video (cf. [19,26,44]). Of course, all of these methods may be considered separately as well as in combination of each other.
In particular, in the case of the tool monitoring one of the most signicant problems is the tool wear. In addition to the above mentioned methods, a monitoring and tool wear prediction method based on accelerometer measurements was presented by Stavropoulos et al. in [33]. Another, innovative approach to process monitoring for micro-milling operations has been studied by Kuram and Ozcelik in [23].
Finally, in order to determine correct machining conditions chatter and/or vibrations have also to be taken into consideration. To that end and to ensure the maximal eciency in monitoring a reliable measurement of micro-forces is required.
Altogether, it can be summarized that accurate measurement micro-forces in micro-milling processes is quintessential for achieving and maintaining optimal results within the manufacturing process.
However, a reliable quantication of micro-forces within measured signals is a highly challenging task. As every other real-world measurement process, the measurement of process forces in micro-range is inevitably accompanied by noise. Therefore, the quality and subsequently the usefulness of the acquired data depends on the ability to separate and remove the noise part from the measured data.
Forces encountered in micro-milling may have low magnitude as compared to the noise of the measurement. Therefore, a low signal-to-noise ratio may be expected in these measurements. In particular, this tends to be the case for high-speed micro-milling processes. Consequently denoising methods tailored to handle micro-forces in micro-processes are of fundamental interest to practitioners.

Approach
The main aim of this paper is to provide a fast and reliable method for removing the noise parts from the force signals measured in micro-milling processes. The presented method has very low computational complexity, allowing for easy implementation.
Currently, already a plethora of standard denoising methods for time-dependent signals is currently available to practitioners. These range from classical Wienerltering (i.e. Fourier ltering) methods (cf. [41]), more recent Gabor Transform and Short Time Fourier Transform (STFT) based methods (cf. [7,17] and the references therein) up to wavelet based methods (cf. [14,16] and the references therein). Further practices include non-Fourier based, nonlinear transforms like the empirical mode decomposition (cf. [10,18]). Finally, we also mention the recently developed, patch-based methods [13].
These general methods have been successfully applied to the problem of denoising force signals in micromilling. For example results for denoising based on independent component analysis (ICA), Short Time Fourier Transform and wavelets are presented by Zhu et al. in [44]. Another approach based on sparse representation of signals in Short Time Fourier Transform has been studied by Zhu et al. in [45]. Further study and denoising of force signals by means of independend component analysis has been discussed in [34,42,43].
However, a typical, underlying assumption of these standard, general denoising methods is that the noise is white and Gaussian. This means that the noise level is constant for all spectral (Fourier) frequencies. This assumption may not be necessary true for measurement of mechanical forces in the micro-range. Exemplary spectra of the three force components where the white noise assumption is invalid are depicted in Figure 1.
Denoising methods operating in the original domain of the signal are essentially smoothing the signal. Therefore, they are virtually low-pass lters. As such they remove high-frequency features important for the analysis of force models. On the other hand, denoising methods operating in the frequency domain are usually thresholding the spectrum. Referring again to Figure 1 it is observed that in order to remove the noise completely the threshold has to be chosen high. Therefore, a lter based on hard-thresholding will not preserve the information on the higher harmonics. Consequently, it can be stated that standard denoising methods generally fail to produce high-quality results in the setting depicted in Figure 1.
Hence, in order to perform the separation reliably and eciently the characteristic properties of the noise and of the force signal have to be determined. To that end, typical power spectra of the three components of the force depicted in Figure 1 are considered. The following, crucial points can be derived from visual inspection of this gure: As already mentioned above, the noise level is dynamic, i.e. the noise-level present in the data is frequency dependent, and thus, is not white. Consequently, the local noise level varies for dierent frequencies. Further, it is observed that the force parts of the signal consist of peaks or outliers in the power spectrum which are distinctively higher than the surrounding, local noise level. Furthermore, these peaks follow roughly the dynamic properties of the noise, i.e. they are high where the noise is high and are low where the noise is low.
Additionally, the signal exhibits a clear fundamental frequency (which is connected to the rotational frequency f as well as the rotational frequency per cutting edge f k , where k is the number of the cutting edges). This fundamental frequency (and the related peaks in the plots of Figure 1) is accompanied by higher order harmonics (and their related peaks in the plots of Figure 1) at integer multiples of the fundamental frequency. In the data depicted in Figure 1 the fundamental frequency is n = 40000 min −1 ≈ 666.667 sec −1 and therefore higher harmonics are expected at multiples of f = 666.667 Hz (rotational frequency).
Previous literature (cf. [2,36]) and research of the authors (cf. [28]) indicates that preservation of the information about the higher harmonics within the ltering process is crucial to the modeling of the underlying forces. Therefore, in the setting described in this paper, an ecient lter must preserve information of higher harmonics of the fundamental frequency.
As observed above, the essential part of the force signal consists of peaks in the spectrum. Therefore, in order to provide ecient denoising the related lter has to be able to detect peaks for variable levels of noise. Further, the lter has to be fast in order to provide the denoised signal in real time, e.g. for monitoring applications. Finally, in order to be easily implementable the lter should have low mathematical and computational complexity.
We close this section with a rough overview of the proposed method. A detailed description is provided in Section 5. Further, a visualization of the main steps is presented in Figure 2. It is also remarked that the proposed method operates on individual force components separately.
In a rst step the (Fourier) power spectrum is computed by means of the Fast Fourier Transform (FFT) from the measured force signal.
Then, this data is employed to derive a peak indicator. The peak indicator is based on a gliding squared coecient of variation (SCV) of the power spectrum. It can be shown that the value of the SCV in absence of peaks has a value which does not depend on the noise level (i.e. SCV is agnostic with respect to the noise level). Further, where peaks are present the SCV has much higher value. The reasons for this behavior are presented in detail in Section 4.
In the next step the Fourier spectrum is decimated according to the value of SCV. This means that the Fourier coecients for all frequencies for which the SCV is smaller than a certain value (chosen by the user) are set to zero, i.e. the Fourier spectrum is thresholded with respect to the SCV.
Finally, the ltered force signal is derived from the decimated spectrum by Inverse Fast Fourier Transform (IFFT).
The gliding SCV is fast and easy to implement, since it can be reduced to two gliding mean lters (cf. Section 4). Consequently, the main computational steps of the approach, namely FFT, SCV, decimation and IFFT are fast and either already available in standard software packages (FFT, IFFT) or very easy to implement (SCV, decimation). Altogether, the proposed approach fullls all of the above requirements and provides a simple yet powerful denoising method for the measurement of process forces in micro-regime. Results concerning the performance of the lter are presented in Section 6.

Experimental setup
All cutting processes were carried out on a DMG Sauer Ultrasonic 20 linear machine tool. The machine tool design and a high speed spindle (maximum rotational speed 42000 min −1 ) allow to conduct micro-milling operations. The workpiece material for the cutting processes was hardened cold working steel 1.2379 with a hardness of 60 HRC. The steel was produced through powder metallurgy and exhibits a ne grained microstructure without primary carbides, making this appropriate for micro-manufacturing.
A Kistler 9119AA2 MiniDyn piezo-electric multicomponent force transducer was used to measure process forces. It allows for the simultaneous acquisition of three orthogonal forces (F x , F y and F z ) or three orthogonal momentums (MA, MB, MC). The sensitivity of the force transducer is given by the manufacturer as less than 2 mN and the force measurement range with −2.5 kN to 2.5 kN for all measurement directions. The force transducer was mounted to an adapter plate, attached to the machine table. A second adapter plate mounted to the force transducer carried the vise for the clamping of the work piece. The F x force component of the force transducer was aligned in feed direction and the F y force component perpendicular to the feed direction. The normal force of the milling operation coincides with the F z force component of the force transducer. The setup is displayed in Figure 3. The measurement chain comprises a Kistler 5019 multichannel charge amplier as well as a computer equipped with an I/O A/D converter card for data sampling. The sampling rate was 30 kHz for all cutting processes. All considered processes were slot milling processes.

Identication of the cutting processes in the data
A measurement campaign consists of a set of three repeated cutting processes which are carried out for an identical set of parameters. An individual cutting process consists of the measurement of the force for two known depth of cut levels a 1 p and a 2 p with a 1 p > a 2 p . The milling tool is moving through the material with the depth of cut a 1 p and after time t 1 changes to level a 2 p . In sum, all three force components F x , F y , F y are being measured in six individual cutting processes within one campaign. For example, the force measurements in x-direction for one whole campaign are presented in Figure 4 (light gray line).
The cutting processes are separated in time. However, the boundaries of these processes are not synchronized with the force measurement aperture. Therefore, in the rst step, these boundaries have to be recovered from the data measured within the whole campaign.  Because of the high level of noise in the measured data it is not possible to determine the exact boundaries of the milling processes and a smoothing step is carried out. After smoothing the data, individual cutting processes at the dierent depths of cut can be identied. In Figure 4 (black line) the result of smoothing with Daubechies wavelet lter (of order 8) is presented. The boundaries of the individual processes the sections divided by dark gray, vertical lines in Figure 4 are then clearly distinguishable. In the subsequent steps every process is considered separately. Only one force component for one depth of cut are presented exemplarily in Figure 4, as the other two were processed in the same way.
However, it is remarked that the smoothed data is only used in order to identify the boundaries of the individual milling processes. The processed data (black line in Figure 4) is too over-smoothed to be useful in further processing. For this reason, in the subsequent steps the original data restricted to the temporal boundaries of each process is employed again.

Statistical peak detection
The proposed lter for the individual milling processes is developed in two parts. In this section a local peakdetection criterion is discussed. This criterion will then be developed into a full denoising method in the next section.
A pure white noise signal, e.g. a vector of length K dened via → u = [u (1) , . . . , u (K) ], can be interpreted as K samples of a normally distributed random variable with zero mean and variance σ 2 . Then, the vector fft(  u )) are sample vectors of two, normally distributed, stochastically independent random variables with zero mean and variance σ 2 · K 2 . Thus, they can also be considered pure-noise signals. Consequently, the vector 2 ) 2 may be interpreted as a sample vector of a Rayleigh-distributed random variable with the scale parameter σ · K/2. The expected value of that random variable is given via π/2 · σ and the variance by 4−π 2 σ 2 , (cf. [15] p. 173). Hence, the so called squared coecient of variation (SCV) of that random variable, which is the ratio of the variance to the square of the expected value, is independent of the scale parameter of the Rayleigh distribution and has the value 4 π − 1 ≈ 0.27. An estimator for SCV is given by: where mean( and Var( are the standard unbiased estimators for the expected value (mean) and for the variance. Of course, this estimator is by itself not expected to be unbiased. However, due to numerical experiments and the analogue computations for a normal distributed random variable (cf. [32] p. 58.) the error is expected to be negligible.
In other words, for pure-noise signals It should be noted that the assumptions on the distribution properties of → v 1 and → v 2 were only used to derive the peak detector. In order for the lter to be sensible   the only important property is that → w is a sample of a Rayleigh-distributed random variable. As can be seen in Figure 5 this assumption is fullled for high-frequency part of the data (higher than 5 kHz). This frequency was chosen since there the noise level is approximately stationary. Obviously, the spectrum is indeed Rayleighdistributed there. Therefore, it can be stated that the data (locally) fullls the assumptions of the lter.
It should be stressed that the detector described above is scale invariant, i.e. the scale parameter of the Rayleigh-distribution respectively the noise level are not parameters within the detector and therefore have neither to be estimated nor to be known. Further, the lter may be applied to all signals which are (locally) stationary in the power spectrum. Finally, it should also be emphasized that the pointwise absolute-value of the Fourier transform is used as the denition of spectral power of a signal. Another common denition for the spectral power is the square of the absolute-value. Of course in that case the SCV can be employed to construct a scale invariant peak detector, too. However, in our experiments the denition used here achieved better performance.

Description of the lter
As was mentioned before, it is assumed that all three force components F x , F y , F y are measured and thus also have to be ltered. The data for every of the six individual experimental parts of the data (cf. Section 3.1) is denoised separately. The three force components for the measured, noisy data of a single process are denoted by → F x , → F y , → F z an index identifying the process is dropped in order to simplify the notation. As the lter treats every component separately the action of the lter is presented only for the force component ] of the data → P x . Further, the data is extended symmetrically at the boundaries. For example, for N = 2 and m = 1 the extended subset is given by [p (3) x , p (2) x , p (1) x , p (2) x , p The parameter N controls the width of typically detected peaks, while the parameter c described in the last section controls the height of typically detected peaks. Together these two parameters control the quality of the denoising.
The proposed lter is a gliding lter of width K = 2N + 1. In particular, every peak is usually surrounded by a window of width 2N + 1 where the noise is not ltered. In the given dataset this window was helpful for the modeling of the micro-forces. However, if such behavior is not desired additional morphological operations [31] may be applied to narrow or even delete that surrounding window.
Since it is known that → F x is real-valued, the left half of the FFT is complex conjugate of the right half (if the vector is considered to be a row vector). As such, both halves contain the same information.
Therefore, subsequently, the vector → G x could also be cropped to the left half, i.e. it is sucient to only consider entries with index 1 ≤ m ≤ L, where L := M/2 and · denotes the rounding towards innity (also known as ceil-function).

Results
In this section the performance of the lter in the spectral domain and in the time domain based on real data is described. Since in the tested dataset the timestamps were equidistant, FFT was used. The algorithms were implemented in MATLAB and have been tested for several parameter settings.

Performance in the spectral domain
In Figure 6 (top row) the spectrum of a milling process' signal is depicted which mainly exhibits peaks at the rst and second multiple of the rotational frequency n = 40000/min (≈ 666.67 Hz). These peaks are recognizable by visual inspection. The width N of the lter was chosen to be 6 in all three force components while the peak detection parameter c was chosen 10 for the F x force component, 7 for the F y force component and 15 for the F z force component of the force signal. In the bottom row of the gure, it can be seen that the lter successfully recognized and retained these peaks. At this point it is once again stressed that the lter is agnostic with respect to the information on the rotational frequency, i.e. that information was not fed to the lter. In Figure 7 the performance of the lter is shown for signals which exhibit peaks at more than only two multiples of the rotational frequency. In this case the width N was chosen to be 12 and the peak detection parameter c was chosen to be 5 in all three force components. As before, the top row of the gure shows that the lter is at least as good as choosing the peaks by visual inspection. Furthermore, the bottom row shows that also in this case the lter recovers the correct peaks at the multiples of the rotational frequency. It is worth pointing out that characteristic peaks at about 250 Hz are recognizable in all three force components as well. Although the mechanical interpretation of this frequency is currently still an open question. It is clear that these peaks have dierent characteristics than the surrounding noise and, therefore, are correctly recognized as 'not-noise' part of the signal.
Notice that in both cases (Figures 6 and 7) the lter is at least as good as traditional lters. This is especially evident for the force component in x-direction (left row) and z-direction (right row).
The considered processes were slot milling processes with feed velocity in x-direction. As such the F y component of the cutting force is expected to have much less pronounced peaks then the other two components. This is clearly visible in the middle column of Figures 6  and 7. In order to remove the noise caused by the machine table (best visible at frequencies 2.5 kHz-3 kHz) a lter based on spectral thresholding would also remove all peaks belonging to the signal, as these peaks have signicantly smaller magnitude than the surrounding noise. However, the presented lter is clearly able to recover the correct peaks for the signal.
The lter was designed to be scale invariant and therefore to automatically adapt to the level of noise present in the spectrum. This fact is depicted in the subplots in the upper right corner of plots of Figure 7. Firstly, these sub-graphs show that information about the force component is present even at higher multiples of the rotational frequency (in the depicted case at the sixth and seventh multiple). Secondly, it shows that the magnitude of this information is smaller than the magnitude at low multiples (e.g. the rst and second multiple). Finally, this shows that the lter works in the way it was designed and successfully retains the force information while removing the noise from the signal, even if the noise level varies within the signal spectrum.
Summarizing, this shows that the presented lter is often as good and in critical cases superior to the established spectral tresholding methods.

Performance in the time domain and comparison with force model
In order to fully evaluate the performance of the lter, the results in the time domain are considered also. In Figures 8 and 9 the measured and ltered F x (top) and F z (bottom) force components for the same process as in Figure 6 is depicted. The force component in y-direction is omitted, as the process is a slot milling process in the x-direction. Consequently that force component is of less interest as it is expected to be much less variable than the other two components (a view which is conrmed by the much lower magnitudes of the relevant spectral information, as shown in Figure  6).
It is remarked that the unltered signal (gray line in Figure 8) is relatively smooth. Therefore, smoothingbased lters are expected to perform poorly for that signal. However, this also means that a smoothnessbased evaluation of the performance of the lter not is not the right choice. Therefore, instead, the ltered force measurement is compared with an well-established force model.
For conventional cutting processes a number of cutting force models have been described in the literature. In the case of milling processes the acting forces are separated into three components, radial force, tangential force and axial force and are denoted by F r , F t and F a respectively. A common approach relies on fact that the force components are proportional to the cross section area of cut A c and can be computed using the specic cutting force components K rc , K tc and K ac . This results in the forces F r = K rc A c , F t = K tc A c , F a = K ac A c . By applying the coordinate transformation operator the force components F x , F y , F z in the Cartesian tool coordinate system are obtained. This cutting force model was rst introduced by Weck and Teipel [40]. Altintas extend this model by adding a second term which additionally takes into account friction, cf. [2, p. 35-47]. The model of Altintas is applied in the remainder of this section. Accordingly, the model for the (instantaneous) cross section area A c (t), the (instantaneous) chip thickness h(t) and the instantaneous angle of immersion φ(t) are given by where the depth of cut a p , the feed velocity v f , rotational frequency n and the number of cutting edges k are not time varying within a single milling process. However, they may have dierent values for dierent milling processes. The quantity t is the time measured in seconds and φ 0 is the initial angle of immersion. All other quantities have the dimensions given in Figures 6  and 7.
The model of Altintas then gives the following formulas for radial-, tangential-and axial-force components: F r (t) = K rc A c (t) + K re a p (8) F a (t) = K ac A c (t) + K ae a p A visualization of the forces in the x, y-plane is presented in Figure 10. Finally, the forces in the Cartesian tool coordinate system are given via F x (t) = −F t (t) · cos φ(t) − F r (t) · sin φ(t) + c x (10) F y (t) = +F t (t) · sin φ(t) − F r (t) · cos φ(t) + c y (11) F z (t) = +F a (t) + c z (12) where the oset parameters c x , c y and c z are additionally introduced to counter the oset caused by the measurement aperture. It should be noted that there are ten unknown parameters φ 0 , K tc , K rc , K ac , K te , K re , K ae , c x , c y , c z . These are recovered from the ltered data by least square optimization. The results of the optimization are depicted in Figure 9 for the F x and F z force component. The component F y is once again omitted for the reasons mentioned above, i.e. the nature of the milling process at hand.
A very good t between the ltered force component (gray line, top of Figure 9) and the optimized model (black line, top of Figure 9) is achieved for the F x force component. This validates the eectiveness of the lter in that case.
The dierence between the ltered signal (gray line, bottom of Figure 9) and the optimized model (black line, bottom of Figure 9) for the F z direction can be accounted for by the fact that the macro-scale milling Fig. 10 Coordinate system in milling operations: x, y-plane (cf. Subsection 6.2 for a description of the involved symbols). models (like the Altintas model) do only describe the chip thickness h in terms of the feed and feed per tooth (which translates to the parameters v f and k in the above denition of h(t)). Therefore, the rotational frequency is dominant in such models in the axial force component F a and therefore also in the F z force component.
Summarizing, the presented ltering is as good as the model based constraints allow it to be. However, the development of force models for milling processes in the micro-scale is currently an ongoing research eort and the comparison of the lter to these models is an interesting future research endeavor.

Conclusions
A new approach for signal denoising based on a statistical peak detection methods is presented. It can be applied to a wide range of noise models. Moreover, no assumptions from the milling process are needed to model the signal. The only necessary assumption is that the signal parts exhibit a peak-like behavior in the Fourier spectrum. The lter is controlled by two parameters which have to be chosen in advance by the user. The ltering method has succeeded to determine the fundamental frequency as well as its harmonics from the input signal.
The successful application of the presented lter combined with further cutting experiments and theoretical research will allow the development of novel cutting force modeling techniques which account for the higher multiples of the rotational frequency within the signal.