Measurement of distorted exponential signal components using maximum likelihood estimation

Measurement of distorted exponential signal components using maximum likelihood estimation

Measurement 58 (2014) 503–510 Contents lists available at ScienceDirect Measurement journal homepage: www.elsevier.com/locate/measurement Measureme...

1MB Sizes 1 Downloads 77 Views

Measurement 58 (2014) 503–510

Contents lists available at ScienceDirect

Measurement journal homepage: www.elsevier.com/locate/measurement

Measurement of distorted exponential signal components using maximum likelihood estimation Linus Michaeli a,⇑, Ján Šaliga a, Jozef Liptak a, Marek Godla a, István Kollár b a b

Dept. of Electronics and Multimedia Communications, FEI Technical University of Košice, Slovak Republic Budapest University of Technology and Economics, Dept. of Measurement, Budapest, Hungary

a r t i c l e

i n f o

Article history: Received 12 December 2013 Received in revised form 9 September 2014 Accepted 11 September 2014 Available online 19 September 2014 Keywords: Maximum likelihood method Exponential testing stimulus signal Dynamic ADC testing Signal decomposition

a b s t r a c t Exponential signal is a suitable stimulus for dynamic ADC testing because of the simplicity of the generating RC circuit. A potential distortion source of the ideal exponential shape is the dielectric absorption of the capacitor, whose effect can be represented by additional superimposed exponential components with longer time constant and smaller peak value. Measurement of the distortion of the exponential signal by using a reference waveform recorder with known nonlinearity is the initial step in the calibration of an ADC testing stand with exponential stimulus, along with the assessment of its uncertainty. Lack of the orthogonality of stimulus signal components makes classical analysis methods difficult to apply. This paper presents a method for measurement of multiexponential signal components as an example of the more general task of signal decomposition where signal components are non-orthogonal. The proper optimization procedure based on the ML method will be presented, which usually reaches the global minimum of the cost function. Effectiveness will be shown by simulation, and by application to measurement of a multiexponential signal acquired by a reference waveform recorder with known error parameters. Ó 2014 Elsevier Ltd. All rights reserved.

1. Introduction Exponential input stimulus signal is attractive for dynamic ADC tests because of the simplicity of the generating RC circuit [1,2]. The estimation of the testing error requires measurement of the distortion of the stimulus signal. Despite the simplicity of the RC circuit there is a distortion in the exponential signal: distortion components are caused by dielectric imperfections of the capacitors in the discharging circuit. The exponential output voltage is connected without any interfacing as dynamic stimulus signal to the input of ADC under test. Distorting elements ⇑ Corresponding author. E-mail addresses: [email protected] (L. Michaeli), jan.saliga@ tuke.sk (J. Šaliga), [email protected] (J. Liptak), [email protected] (M. Godla), [email protected] (I. Kollár). http://dx.doi.org/10.1016/j.measurement.2014.09.024 0263-2241/Ó 2014 Elsevier Ltd. All rights reserved.

in the circuit model of a capacitor, which represent the dielectric absorption, are additional serial RiCi (for i = 2,3. . .) branches connected in parallel to the primary capacitor C1 (Fig. 1). Each parasitic RiCi branch generates a superimposed exponential component with longer time constant and smaller peak value than the main capacitor, and also slightly changes the main time constant [3]. In order to determine the testing uncertainty, the distorting exponential components must be measured beforehand by a reference waveform recorder (RWR) with known nonlinearity. For the determination of these parameters from the data record a convenient method for estimation of exponential components is needed. Since the exponential functions are not orthogonal, the task to determine parameters of this function is not easy. For this specific problem there are known some determination methods which work well for good Signal-to-Noise Ratio

504

L. Michaeli et al. / Measurement 58 (2014) 503–510

C1

R2

Ri

C2

Ci

R1

xS(t)

U0

Fig. 1. Circuit model of the exponential signal generator.

(SNR). The most known is Prony´s method [4,10]. In general the number of the superimposed exponential components is unknown. An additional superimposed component is the noise from various interfering sources. Here the maximum likelihood method is suggested as an alternative to Prony’s analytical method which is the straightforward method of identification of the signal components of the stimulus signal acquired by the RWR. Smaller sensitivity to the superimposed noise and the possibility to eliminate the effect of the nonlinearity of the recorder are the main advantages of maximum likelihood method in comparison to Prony’s method. The ML method is characterized by larger computational complexity but with the advantage of lower sensitivity to clipping of the signal and to the shape of the quantization noise [6,7,8]. As a parameter estimation method, ML is asymptotically optimal. It also allows to involve the a priori known nonlinearity of RWR in the calculations. A further advantage is the possibility to avoid the convergence into the local minima by appropriate setting of the noise parameters. The computational complexity is related to the multidimensional optimization of the parameters (Ai, Bi). Anyway, the aim is to achieve best matching of the considered input signal to the true code intervals corresponding to the recorded digital shape. The best matching is assessed by the maximum likelihood procedure, taking into account also the true transfer function of the RWR. The speed of convergence of the parameter estimation also depends on the setting of the initial values at the beginning of the optimization procedure. Determination of the multiexponential signal components by the ML method will be studied in Section 2. The efficiency of the proposed ML method under various conditions will be studied by simulation, and will be compared to Prony’s method in Section 3. Results of the experimental verification and assessment of the estimation uncertainty will be shown in Section 4. A brief summary of the achieved results and further development of the ML method in general signal identification will be provided in Section 5.

2. Mathematical model As mentioned in [2], dynamic ADC testing using exponential stimulus signal is dominated by the time constant and by the peak value of the basic exponential component. The parameters of basic components are determined by the precise circuit model which generates exponential stimulus, taking into account dielectric absorption. Distorting elements in the circuit model of a capacitor,

which represent the dielectric absorption, are the additional serial RiCi (for i = 2, 3. . .) circuits, connected in parallel to the primary capacitor C1 (Fig. 1). Each parasitic RiCi circuit generates a superimposed exponential component with longer time constant and smaller peak value than the main capacitor, and also slightly changes the main time constant [3]. The output signal xS(t) for only one absorption circuit R2,C2 in parallel with the primary capacitor C1 and discharging resistor R1 can be determined analytically (1) by formula

xs ðtÞ ¼ A1 eB1 t þ A2 eB2 t

ð1Þ

The values of time constants and the peak values of two exponential components are:

B1;2 ¼

1 ð1  aÞ; T

sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2C 1 R1 C 2 R2 4C 1 R1 C 2 R2 ; a¼ 1 T¼ C 1 ðR1 þ R2 Þ þ C 1 R1 ½C 1 ðR1 þ R2 Þ þ C 1 R1 2   U0 1 A1 ¼  B2 ; A2 ¼ U 0  A1 B1  B2 C 1 R1 ð2Þ The final exponential signal with distorting superimposed components must be known in order to assess the ADC testing stand uncertainty. The method proposed in the paper [4] utilizes Prony’s analytical method for the signal acquired by RWR for equidistant sampling. The measurement precision was increased by the measurement of the time instants for known values of DC threshold voltages and by the application of Prony’s modified method [9,10]. Drawback of both analytical methods is their high sensitivity on the superimposed noise and their inability to estimate the DC parameter. 2.1. Modeling of the input samples The exponential signal for ADC testing distorted by the additional exponential components with bigger time constants is represented by

xS ðtÞ ¼ A1 eB1 t þ

L X Ai eBi t þ C þ nðtÞ

ð3Þ

i¼2

where Ai << A1, Bi << B1, i = 2, 3,. . . While parameters A1 and B1 represent the basic exponential signal determined by the primary values R1, C1, the parameters Ai and Bi (i = 2,3,. . .) represent the distorting exponential components. Number L is the number of all exponential components assumed in stimulating signal (3). The constant C in (3) describes the offset of the whole exponential signal. The distorted multiexponential signal is corrupted by additional (mainly thermal) noise of the analog components and by interferences from external sources (n(t)). The superimposed input noise n(t) is assumed to have Gaussian distribution with zero mean and variance rn. The sampled input signal xs(j) with sampling period ss is represented by (4) before its quantization in the RWR.

xs ð jÞ ¼

L X Ai  eBi jsS þ C þ nð jsS Þ i¼1

ð4Þ

L. Michaeli et al. / Measurement 58 (2014) 503–510

2.2. Derivation of the likelihood function in the digital domain The quantization levels Tid(k) in the ideal case are equidistant: Tid(k) = Q(k + 0.5), where Q is the code bin width, and k = 0, 1, 2,. . ..,2N1 Real quantization levels T(k) are corrupted by the integral nonlinearity INL(k) of the RWR: T(k) = Tid(k) + INL(k)Q. Let us suppose that the shape of the channel edge for level T(k) is of Gaussian distribution with the variance r. The transient code levels T[k] and T[k1] are neighboring quantization levels for code k at the RWR output. The probabilistic profile of the code bin k with actual transient code level is Pðk; xS ; rÞ ¼ F ðT½k  1; xS ; rÞ  F ðT½k; xS ; rÞ, where F(z, l, r) is the normal cumulative distribution function with mean l and standard deviation r. Let the j-th sample from the RWR output be d(j), and let the signal at the input of RWR be xs(j) described by (4). The best estimation of the input signal xs(j) on the basis of the output d(j) is the one which best matches the corresponding digital samples d(j). Mathematical expression for the matching probability for any time j is expressed by one of the probabilities in (5).

Pðdð jÞ ¼ 0jxS ðjÞ; rÞ ¼ 1  F ðT½0; xS ð jÞ; rÞ Pðdð jÞjxS ð jÞ; rÞ ¼ F ðT½dð jÞ  1; xS ð jÞ; rÞ  F ðT½dð jÞ; xS ð jÞ; rÞ      h i  P dð jÞ ¼ 2N  2 jxS ð jÞ; r ¼ F T 2N  2 ; xS ð jÞ; r ð5Þ The probability of the input sample xS(j) matching with the channel profile for the corresponding sample d( j) is determined by one from the equation (5). Overrun of RWR Full Scale Range at lower and upper scale limit is described by the first and third rows in (5), respectively. The joint probability Pfinal equals the product of the individual probabilities for any of J samples d( j) (independent events):

Pfinal ¼

J Y

Pðdð jÞjxS ð jÞ; rÞ

ð6Þ

j¼0

2.3. Maximization of the likelihood function The optimization task is to find constants (Ai, Bi, C) in xs( j) (4) and r in (5) for which the joint probability Pfinal (6) achieves its maximum. In order simplify the task, the minimum of the negative log-likelihood cost function CF(p) = log(Pfinal) is sought, where p represents the set of parameters Ai, Bi, C, r. J X CFðpÞ ¼  log Pfinal ¼  log Pðkð jÞ; xð jÞÞ

ð7Þ

j¼0

There are several optimization methods, e.g., the gradient based method [5], improved in [8] to minimize this. The optimization method chosen for following experiments was the differential evolution strategy [11,12,13]. This algorithm represents optimization strategy suitable for the objective functions nondifferentiable, noisy with many local minima. It is slower than gradient-based method but is resistant on the stochasticity of the cost

505

function. Artificial adaptation of the standard deviation r in the probabilistic description of the RWR code profile allows to control the optimization procedure with the goal to avoid the convergence to local minima. Increasing r smoothes the cost function. That is why the optimization procedure starts with the larger value r (smoothed CF) and its value is decreased along with the optimization process.

3. Implementation of the ML method for signal decomposition in simulation The robustness of the proposed ML method against the number of samples and against superimposed noise was studied by the simulation tool developed in the LabVIEW environment (Fig. 2). The multiexponential input signal xs(jss) with adjustable signal components (Ai, Bi, C) was generated by a software generator (SW generator). The signal xs(jss) was acquired by the RWR model with adjustable INL(k). The nonlinearity of RWR INL(k) was modeled using nonlinearity obtained by the standardized testing procedure from DAQ board by NI PCI 6024 card (Fig. 4). The resolution of NI PCI 6024 is 12 bit and its input voltage range is FSR = ±10 V. The equivalent quantization step is Q = 4.88 mV. The applied sampling frequency was 100 kS/ s. These RWR parameters: sampling frequency, FSR, and resolution were used invariably in all simulations. The real transient code levels T(k) were used in the RWR model and in the channel profiles (5). The difference between known parameters (Ai, Bi, C) of the SW generator and estimated ones (eAi, eBi, eC) allowed to study the impact of RWR parameters (resolution, acquisition frequency and the number of samples) on the ML optimization results. Three exponential components were simulated at the output of the discharging circuit. The simulations (Fig. 2) were used first because in this case values of the signal parameters are known exactly, and accuracy of their estimation by the proposed method can be observed. The results for two simulated signals for noiseless signals and signals (4) with superimposed Gaussian noise n(t) with standard deviation h are shown in Table 1. The differences between estimated and true parameters for noiseless signals and noisy signals (4) are small. However, simulations showed high sensitivity of the optimization process on the initial conditions. Increasing the number of the exponential functions (4) requires good choice of the initial conditions. The suitable way how to estimate the initial values Ai, Bi, is utilization of the modified Prony’s method [4] for the first estimate of exponential components from a reduced amount of recorded samples, in spite of its sensitivity to the superimposed noise. The final optimization is performed by the minimization of the ML cost function CF(p) (7). The simulations were performed for different values of standard deviations r of the superimposed noise. The estimation errors of parameters (Ai, Bi, C) were evaluated by relative error dA = (eAiAi)/Ai, dB = (eBiBi)/Bi. Achieved averaged results of 10 repetitive simulations for different values of noise are shown in Fig. 3. Simulations confirmed

506

L. Michaeli et al. / Measurement 58 (2014) 503–510

SW generator of multiexponential signal

Optimization procedure

Ai, Bi,C

Generation of best estimated multiexponential

xS(t) RWR parameters: Resolution, INL(k) Sampling freq. Number of samples

e

RWR model

e

Ai,eBi,eC

x(j)

d(j)

Calculation of CF

CF

Fig. 2. Block scheme of the simulation environment.

0.5

Estimation error of parameters A and C

Table 1 Results of estimation by ML method for two signals with chosen parameters with and without noise for 5000 samples taken with sampling period ss = 10 (ls).

0.3

A1

0.2

A2

A

[%]

0.4

C 0

0.2

Standard deviation

0.4

0.6

0.8

1

of superimposed noise [LSB]

Estimation error of B parameters

2.5

[%]

2 1.5

B

B1 1

B2

0.5

B3

0

0

0.2

Standard deviation

0.4

0.6

0.8

1,2

0,6

INL (k)

Estimated values at condition: h=0

h = 0.8 LSB

Signal 1 A1 (V) A2 (V) A3 (V) B1 (ms1) B2 (ms1) B3 (ms1) C (V)

5 1 0.5 5 0.50 0.05 0.1

5.00 0.97 0.52 4.99 0.51 0.06 0.09

5.01 1.06 0.39 4.99 0.41 0.20 0.14

Signal 2 A1 (V) A2 (V) A3 (V) B1 (ms1) B2 (ms1) B3 (ms1) C (V)

9.8 5.3 4.8 0.89 0.57 0.37 0

9.93 5.46 4.67 0.90 0.58 0.366 0.01

9.34 5.66 5.07 0.89 0.55 0.362 0.03

1

of superimposed noise [LSB]

Fig. 3. Estimation errors for data acquired by the 12 bit RWR for different standard deviations of superimposed Gaussian noise. The parameters of signal used in simulation are equal to parameters of signal 1 in Table 1.

0

The increase of the time window with a constant sampling interval improves the estimation precision as long as the signal does not die out in the window. For the next study, the nonlinearity of RWR was achieved by the standardized ADC tests and is shown in Fig. 4. Let us consider modeled parameters (Ai, Bi, C) of the simulated signal (4) for L exponential components. The estimated parameters are labeled as (eAi, eBi, eC). Further, the error of the parameter measurement e is determined by the formula

1 e¼ 2L þ 1

-0,6

-1,2

Simulated value of parameter

A3

0.1 0

Signal parameter

sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi XL e Ai  Ai 2 XL e Bi  Bi 2 e C  C 2 þ þ 1 1 C Ai Bi ð8Þ

0

1024

2048

3072

4096

code bin k Fig. 4. Integral nonlinearity of the utilized 12 bit RWR.

the assumption that errors increase with increasing noise level. The precision of the measurement of signal components is influenced by the many factors like time window, sampling frequency and RWR resolution and nonlinearity.

Above mentioned error definition was chosen because of main goal of the studied method to estimate signal parameters without focusing on the approximation of the signal shape. Nevertheless, the results of different approximations of the experimental data are presented in Section 4. The next step of simulation study was the comparison between the proposed ML estimation method and Prony’s method. Fig. 5 shows errors in the estimation signal parameters by both methods with and without superimposed noise to the input signal xS(t). The estimation errors e are quite similar if there is no noise (Fig. 5.a). In that case

507

0,25

Error of parameter measurement [%]

Error of parameter measurement [%]

L. Michaeli et al. / Measurement 58 (2014) 503–510

(a) ML method

0,20 Pronys' method

0,15

0,10

0,05

0,00

5

10

15

20

25

0,50 0,45

(b) ML method

0,40

Pronys' method

0,35 0,30 0,25 0,20 0,15 0,10 0,05 0,00

30

5

10

Time window [ms]

15

20

25

30

Time window [ms]

Fig. 5. The error of parameter estimates (a) without Gaussian noise at the input and (b) with superimposed Gaussian noise with standard deviation r = 0.8 Q. Simulated number of bits N = 12.

12

1st iteration

2nd iteration

3rd iteration

10 8 A1

6

A2 A3

4 2 0 20

5

2

1

0.5 20

5

2

1

0.5 20

5

2

1

0.5

[Q] 1.2

1st iteration

2nd iteration

3rd iteration

1 0.8 0.6

B1 B2

0.4

B3

0.2 0 20

5

2

1

0.5 20

5

2

1

0.5 20

5

2

1

0.5

[Q] Fig. 6. Estimation of multiexponential signal with Ai = (7 (V), 3, (V) 1 (V)), and Bi = (1 (ms1), 0.3 (ms1), 0.1 (ms1)) for adapting values of the standard deviation r of the channel profile (6).

ML method has no advantage when compared to Prony’s method. Noise superimposed to the input signal shows superiority of the ML method (Fig. 5.b). It is caused by the fact that results of Prony’s calculus are sensitive to the quality of samples of signal xS(t). On the other hand the ML method requires good initial conditions. Prony’s method could be used as method to find initial conditions for the successive optimization by ML method which converged under each

noise condition to the optimal solution with highest accuracy. ML method is used in the next step for the improvement of measurement accuracy. The drawback of the procedure described above is the computational complexity of Prony´s method, used to obtain a good initial condition. Therefore the authors propose another optimization procedure which starts from a simply estimated peak value and time constant of the fun-

508

L. Michaeli et al. / Measurement 58 (2014) 503–510

damental exponential component. These were estimated from two time samples d(j1) and d(j2) by Eq. (9):

dðj1 Þ; dðj2 Þ ) B1 ¼

  1 dðj1 Þ ; ln ðj1  j2 Þ dðj2 Þ

A1 ¼

dðj1 Þ e eB1 J1

The next iteration step is looking for additional parameters (A3, B3, C) of signal xS(j) using fixed parameters achieved from previous iterations (A1, B1, A2, B2). If the number of exponentials is L = 3 even the constant value C is optimized in that iteration. Here again the optimization of CF(p) (7) is performed with decreasing r. We have observed that the parameters converged to the global minimum of the CF(p) after each step. Fig. 6 shows results achieved in three consecutive iterations for variable values of r. A fundamental problem in the measurement of multiexponential signal components from the time record is that the number of components in Eq. (3) is not known. In the previous paragraphs an iterative method for gradually increasing the number of components in the stimulus signal (3) was proposed and evaluated. The experiments confirm the assumption that with step by step iterative estimation of components of the input signal xS(t), the parameter values achieved in the previous iteration remain constant as long as the number of exponentials in the estimated signal is less or equal to the number Lm of exponentials in the input signal. Exceeding the number of estimated components over the number Lm of the input signal is accompanied by changes of the already estimated values (Ai, Bi) from the previous iteration steps. The proposed strategy was assessed on a simulated input signal with two exponential components (Lm = 2) with A1 = 11.85, A2 = 1.6, B1 = 1, B2 = 0.3, C = 0. As shown in Fig. 7, the parameters of estimated signal achieved from

ð9Þ

The procedure reaches the optimum iteratively. Each iteration employs a fixed number of exponential components. The first iteration uses only one exponential. Each next step adds one additional component. Within one iteration the standard deviation r in (5) is decreased step by step, what results in the channel profile sharpening. The first iteration by the ML method is aimed to optimize the parameters (A1, B1) of signal (4), starting from its rough estimates (9) and zero condition for higher components (Ai, Bi, C) for index i = 2,. . . L. The optimization of the cost function CF(p) starts with large value of the standard deviation r in the channel profile (5). The next step in the first iteration performs optimization of CF(p) (7) with decreased standard deviation r of channel profile. It represents the approach when the approximated signal shape xS(j) matches the measured record d(j) more precisely with sharper channel profile. The second iteration step starts from the final values (A1, B1) achieved at the end of previous step. The optimization of CF(p) is performed similarly to the preceding iteration for decreasing value of r. At the end of the second iteration the optimal parameters (A1, B1, A2, B2) are determined. 16

1st iteration

2nd iteration

3rd iteration

14 12 10 8

A1

6

A2

4

A3

2 0 -2

20 5

2

1 0.5 0.2 20 5

2

1 0.5 0.2 20 5

2

1 0.5 0.2

[Q] 1.2

1st iteration

2nd iteration

3rd iteration

1 0.8 0.6

B1 B2

0.4

B3

0.2 0 20 5

2

1 0.5 0.2 20 5

2

1 0.5 0.2 20 5

2

1 0.5 0.2

[Q] Fig. 7. Multiexponential signal estimation with Lm = 2 exponentials of input signal with A1 = 11.85, A2 = 1.6, B1 = 1, B2 = 0.3, C = 0 for and leakage among signal parameters of the over-determined estimated signal. Parameter r represents standard deviation of the channel profile (6).

509

L. Michaeli et al. / Measurement 58 (2014) 503–510

7

3,33 3,32

6

ML

x s [V]

3,31

x s [V]

5

3,3

4

3,29

3

3,28

2

2020

2030

2040

2050

2060

Sample number

RWR

1 0

3,27 2010

ML

RWR

0

5000

10000

Sample number Fig. 8. Measured output voltage xS(t) (V) and estimated by ML method.

previous iteration remain constant and are converging to the true values for iterations with L = 1 and L = 2. If the number of iteration L is higher than Lm, estimated values of exponential components are changing from those achieved in the previous iteration. This effect can be considered as leakage among exponential components of the over-determined estimated signal.

Table 2 Measurement of the signal parameters of the circuit in Fig. 1. calculated analytically (1), (2), and achieved by ML estimation. Parameters

Calculated analytically

Measured by ML method

A1 A2 B1 B2 C

7.06 (V) 1.94 (V) 4.71 (ms1) 0.488 (ms1) 0 (V)

6.96 (V) 1.97 (V) 4.78 (ms1) 0.488 (ms1) 0.022 (V)

4. Experimental results The conclusions drawn from the analytical study and from simulations were verified by experimental measurement and using results from the SPICE model of a real circuit in Fig. 1. In order to suppress unpredictable distortion effects, polypropylene capacitors by WIMA with extremely low absorption were used. Capacitor C1 = 22 nF represented the main capacity, the R2C2 circuit modeled the absorption effects. The values of capacitors C2 = 18 nF and resistors R2 = 58 kX ensured that the residual absorption effects of the WIMA capacitors were hidden. (Fig. 1). The discharging resistor was R1 = 12 kX. The analytical expression (1), (2) for basic and second exponential function in the signal xS(t) was compared with the signal from simulator SPICE. The simulator served as a proof of conformity of the real circuit parameters with those assumed in the analytical expression (1). Finally, the measured transient output from circuit simulator was compared with the output d(j) acquired by the RWR from the real circuit. The errors caused by transient effects of the switch were excluded from processing by taking samples d(j) from the data record only after the transient time of the switch. The length of the time record was 5000 samples, as in the simulation. This length covers the most important part of processed signal and excludes the fade-out signal tail from processing. Fig. 8 shows the measured output signal and measured signal with two exponential components. Measurement of the signal parameters from the output of the tested circuit are shown in Table 2. Minimization of the cost function (6) was performed using the differential evolution algorithm with adaptive improvement of r. The optimization procedure started with r = 50Q and was finished with r = 2Q. Table 2 com-

Fig. 9. Errors between measured and estimated function by Prony’s method, ML and LS estimation.

pares the signal parameters calculated analytically (1), (2) and achieved by the ML method. The experiment also proved that the measurement of the DC parameter C is improving proportionally with lengthening of the almost fade-out part of the data record. But this lengthening does not improve the accuracy of the parameter measurement. Finally, the quality of the estimation was also approved by the comparison of the recorded signal achieved by experimental measurement and by the least squared (LS) polynomial approximation of 10th order, Prony’s method and ML method. The quality was evaluated by difference between measured and its estimations by mentioned methods (Fig. 9).

510

L. Michaeli et al. / Measurement 58 (2014) 503–510

The ready-made routine for LS polynomial approximation was taken from LabVIEW signal processing toolbox and applied on the same record. The different orders of polynomial approximation were examined. The polynomial approximation of 10th order could be supposed as optimal balance between computation complexity and achieved precision. Prony’s method utilized 12 subpolynomials as it was used for a similar signal in [10]. The ML method provides estimation with lowest error and without oscillation between measured output signal and its estimation. 5. Conclusions The well known commonly used method for estimation of multiexponential signal is Prony’s method. The paper presents new alternative estimation method based on ML optimization. The paper was aimed on the comparison of both methods using simulation of taking into account a real RWR and by measurement performed on the circuit with known distorting elements. The performed simulations and experimental results confirmed the assumption that ML method is less sensitive on the superimposed noise than Prony’s method. The ML method is optimal even in presence of quantization noise, although this is not normally distributed. Such a noise is always present in the acquired quantized multiexponential signal. Moreover ML method allows inclusion of the a priori known RWR nonlinearity into estimation process. The proposed simulation model of the RWR and signal generator allowed to study sensitivity of both methods at different conditions in the process of the parameter estimations. Simulation indicates the existence of optimal time window of the recorder signal, whose length and number of samples depend on the measured signal components. Another advantage is the possibility of multiexponential signal measurement without a priori knowledge about the number of its components. The simulated results were verified by a real discharging circuit. The circuit was designed with known parasitic components. The output signal was compared with analytically calculated form as well as with the results from the circuit simulator SPICE. The proposed ML method is also suitable to identify parasitic components representing dielectric absorption in the capacitors. The primary capacitor which is in parallel to the parasitic components represented by RiCi, masks their influence on the input impedance. This fact deteriorates classical impedance measurement methods by

decreasing their accuracy. The measurement uncertainty could be improved by the implementation of signal component measurement using signal decomposition by the ML method. Acknowledgments The work is a part of the project supported by the Science Grant Agency of the Slovak Republic (No. 1/0281/ 14) and the project supported by the Slovak Research and Development Agency under the contract No. APVV0333-11. Moreover this work is the result of the project Center of Information and Communication Technologies for Knowledge Systems (Project Number: 26220120020) supported by the Research & Development Operational Program funded by the ERDF. References [1] J. Šaliga, L. Michaeli, R. Holcer, Noise sensitivity of the exponential histogram ADC test, Measurement 39 (3) (2006) 238–244, http:// dx.doi.org/10.1016/j.measurement.2005.11.004. [2] J. Šaliga, L. Michaeli, M. Sakmar, J. Buša, Processing of bidirectional exponential stimulus in ADC testing, Measurement 43 (8) (2010) 1061–1068, http://dx.doi.org/10.1016/j.measurement.2010.03.007. [3] K. Kundert, Modeling dielectric absorption in capacitors (2008), . [4] D.L. Carni, et al., Measurement of the exponential signal distortion, in: IEEE Proceedings of I2MTC 2012: May 13–16, Graz, Austria. – Graz, pp. 2191–2195. doi: 10.1109/I2MTC.2012.6229424. [5] L. Balogh, I. Kollár, A. Sárhegyi, Maximum likelihood estimation of ADC parameters, in: IEEE Proceedings of IMTC 2010. Austin, United States of America, 03/05/2010-06/05/2010. pp. 24–29. doi: 10.1109/ IMTC.2010.5488286. [6] Y. Gendai, The maximum-likelihood noise magnitude estimation in ADC linearity measurement, in: IEEE Trans. on IM, vol. 59, no.7, July 2010, pp. 1746–1754. doi: 10.1109/TIM.2009.2028218. [7] I. Kollar, J.J. Blair, Improved determination of the best fitting sine wave in ADC testing, IEEE Trans. on IM 54 (5) (2005) 1978–1983, http://dx.doi.org/10.1109/TIM.2005.855082. [8] J. Šaliga, I. Kollár, L. Michaeli, J. Buša, J. Lipták, T. Virosztek, A comparison of least squares and maximum likelihood based sine fittings in ADC testing, Measurement 46 (2013) 4362–4368, http:// dx.doi.org/10.1016/j.measurement.2013.05.004. [9] K. Holmstrom, J. Petersson, A review of the parameter estimation problem of fitting positive exponential sums to empirical data, Appl. Math. Comput. 126 (1) (2002) 31–61, http://dx.doi.org/10.1016/ S0096-3003(00)00138-7. [10] L. Coluccio, A. Eisinberg, G. Fedele, A Prony-like method for nonuniform sampling, Signal Process. 87 (10) (2007) 2484–2490, http:// dx.doi.org/10.1016/j.sigpro.2007.03.011. [11] R.M. Storn, J.A. Lampinen, Differential evolution: a practical approach to global optimization, Springer Verlag, 2005, ISBN 9783-540-20950-8. [12] U.K. Chakraborty (Ed.), Advances in Differential Evolution, Springer Verlag, 2008, ISBN 978-3-540-68827-3. [13] R.M. Storn, Differential evolution homepage: April 2014.