Possible solutions to avoid a loss of performance due to RF collisions are mainly based on adaptive frequency hopping [5-7] and concern only interferences between a device using static RF channels, such as 802.11 devices and a device using dynamic RF channels allocation such as a Bluetooth device. In such an approach, a Bluetooth device experiencing interferences on a particular channels group will stop using the corresponding frequencies for a certain time.

In the case of Bluetooth piconets, the technique presented above is not suitable as the piconets implement the frequency hopping algorithm independently, without coordination. These interferences are dynamic as each piconet may not know the frequency channel used by the others. Some works have been done to evaluate the performances of the frequency hopping algorithm in Bluetooth device in such a case [8]. Although the results show that the coexistence is possible, these statistical evaluation does not contribute to reduce interferences between piconets.

In this paper, we present a method to ensure a free interference coexistence between Bluetooth interfaces implemented in the same device via a judicious setting of the Bluetooth frequency hopping algorithm. We do not address the static interferences case. The paper deals with the realization of a multiple interfaces Bluetooth access point [9]. To implement multiple Bluetooth interfaces (2 to 5) in a single access point into a single device, interferences between the RF Bluetooth interfaces must be suppressed. Indeed, each interface acts as a Bluetooth master creating several piconets in the same location. As the Bluetooth interfaces are located very close to each other, no propagation channel model nor geographical repartition need to be defined. We consider the worst case where all jammers are located in the same place radiating the same power as the useful signal. We do not seek to exactly define the probability of collision, nor to obtain a statistical model of these interferences but we show how to practically suppress them.

We first present the Bluetooth frequency hopping algorithm. Then we explain why probabilistic models are not adequate to understand interferences between piconets in our case where Bluetooth interfaces must experience no interference in order to create a powerful Bluetooth access point. As we explain in the third part, even if collisions are not very probable, in some specific cases, two piconets can experiment severe collisions.

After a detailed study of the Bluetooth frequency hopping algorithm, we give a solution to obtain a interference-free environment. For clarity's sake, the

method is detailed only in the case of asynchronous single slot packets exchanged in the case of two piconets.

The solution proposed has been implemented for validation purpose in a real prototype system which exhibits very interesting performances. Some basic results obtained with the prototype are presented as an illustration of the method efficiency.

## 2. The Bluetooth Standard

Bluetooth provides ad hoc networking using a low cost, low power, short range radio interface [10]. It uses a master driven time division duplex (TDD) scheme and provides both synchronous and asynchronous channels for voice or data transfers [11]. In order to reduce the effect of interferences with others devices, a slow frequency hopping sequence along with a fast Automatic Repeat reQuest (ARQ) scheme have been planned. The allocated frequency band is divided into 79 RF channels, 1 MHz wide each ranging from 2402 up to 2480 MHz. A Bluetooth communication channel on which RF packets are exchanged is defined by the successive channels used along with slot timing. A set of Bluetooth devices connected together forms a piconet. The slot timing and channel sequence are determined by the master's Bluetooth device address (BD_ADDR) and native clock (CLKN) as can be seen in Figure 1.

The BD_ADDR is a 48 bits address which is composed as follows:

$$BD_ADDR_{47-0} = NAP_{15-0} | UAP_{7-0} | LAP_{23-0} \quad (1)$$

where NAP (non-significant address part) and UAP (upper address part) are fixed (assigned to a company by IEEE) and LAP (lower address part) is freely chosen by the manufacturer.

<figure><img src="image_69.png" /><figcaption>Fig. 1. Typical packets exchange into a Bluetooth piconet between master and slaves.</figcaption></figure>
<figure><img src="image_70.png" /><figcaption>Fig. 2. Selection kernel of the Bluetooth system.</figcaption></figure>

Slaves have to adjust their slot and Bluetooth clock offsets in order to match the master's Bluetooth clock:

$$ \mathrm{CLK}_{\mathrm{slave}} = \mathrm{CLKN}_{\mathrm{master}} + \mathrm{clock\_offset} + \mathrm{slot\_offset} \quad (2) $$

The consecutive RF channels used are calculated in a pseudo-random way to satisfy the Federal Communications Commission (FCC) requirements concerning the 2.4 GHz ISM band. The algorithm, called the Bluetooth selection kernel, uses the master identity (BD_ADDR, fixed for a Bluetooth device) and clock (CLKN, incremented by one every 312.5 µs) as inputs. The slave also uses master identity, and its own clock as given by Equation (2).

The Bluetooth selection kernel is represented in Figure 2, according to Reference [11].

The notations used in this paper are those of Figure 2. Inputs in connection mode are given by:

$$ X = \mathrm{CLK}_{6-2} \quad (3) $$

$$ Y1 = \mathrm{CLK}_{1} \quad (4) $$

$$ Y2 = 32 \times \mathrm{CLK}_{1} \quad (5) $$

$$ A = \mathrm{ULAP}_{27-23} \oplus \mathrm{CLK}_{25-21} \quad (6) $$

$$ B = \mathrm{ULAP}_{22-19} \quad (7) $$

$$ C = \mathrm{ULAP}_{8,6,4,2,0} \oplus \mathrm{CLK}_{20-16} \quad (8) $$

$$ D = \mathrm{ULAP}_{18-10} \oplus \mathrm{CLK}_{15-7} \quad (9) $$

$$ E = \mathrm{ULAP}_{13,11,9,7,5,3,1} \quad (10) $$

$$ F = 16 \times \mathrm{CLK}_{27-7} \mod 79 \quad (11) $$

$$ Z' = (X + A) \mod 32 \quad (12) $$

$$ Z[4] = Z'[4] \quad (13) $$

$$ Z[0:3] = B \oplus Z'[0:3] \quad (14) $$

where CLK represents the input clock, and the ULAP bits are the 32 least significant ones (LSB) of the

BD_ADDR, that is

$$ \mathrm{ULAP} = \mathrm{UAP}_{7-0} | \mathrm{LAP}_{24-0} \quad (15) $$

With this presentation of the basics of the Bluetooth technology, it is possible to describe the method we propose to ensure coexistence and to enhance performance in high interferences situations.

### 3. How a Bluetooth System can Coexist With Other Ones

It is possible to evaluate how Bluetooth piconets can coexist. Each piconet has its own slot timing based on the master's one, and its own frequency sequence based on the master's Bluetooth clock and BD_ADDR. From hereon, it is clear that many piconets in close vicinity to each other will experience mutual interferences. An illustration is given in Figure 3.

In this case, there is a non-zero probability to have overlapping packets in the two piconets. Moreover, the same frequency may be used during a slot. For instance, if $f(k) = g(k)$ the packets exchanged in the two piconets during the slot are lost due to a collision.

As an example, consider two Bluetooth piconets. We assume for simplicity that they are slot aligned together. Since there are 79 channels, the probability that the two piconets use the same channel at the same time is given

<figure><img src="image_71.png" /><figcaption>Fig. 3. Reference timings for two different Bluetooth piconets.</figcaption></figure>
<figure><img src="image_72.png" /><figcaption>Fig. 4. Mobile average computation of number of collision over 40 slots.</figcaption></figure>

<figure><img src="image_73.png" /><figcaption>Fig. 5. Mobile average computation of number of collision over 40 slots. Local view.</figcaption></figure>

[12] by Equation (16):

$$P = \frac{1}{79} \approx \frac{1266}{100} \qquad (16)$$

Statistical simulations were performed to evaluate the residual collisions between piconets and they show a good behavior [3,8]. But the statistical approach masks some particular dramatic situations if we want to implement a multiple Bluetooth access point. For example, problems arise when computing the channels used by two piconets having addresses and clocks ULAP₁ = 0x09808105, CLKN₁ = 0x0, ULAP₂ = 0x0 and CLKN₂ = 0x2731000, respectively. Indeed, by computing the number of collisions using Equations (3) to (11) combined with the selection kernel box of Figure 2, it shows that there is such a high number of collisions that, for a long time, it is not possible to establish a connection. The results of the computation are given in Figures 4 and 5.

This example in a particular case of Bluetooth Addresses and Bluetooth Clocks association shows that in some cases, probabilistic-based models do not really fit the reality of Bluetooth piconets life.

Although, it is possible to find a lot of situations similar to that one, in most cases the number of collision will not be as high as in the above example.

To optimize piconets coexistence, it is necessary to know how often such dramatic situations occur in real conditions. The number of possibilities of frequency hopping sequences for each Bluetooth piconet is so important that it is not realistic to try to calculate all

the possibilities of collisions between two or more piconets during a communication. Since we could not obtain collision statistics for every Bluetooth piconet, we decided to adopt a simulation-based approach to evaluate collisions.

## 4. What About Bluetooth Piconets Collisions in the Real World?

In order to answer this question, we carried out some simulations. The method used is:

1. choose a couple of random ULAP and CLKN;

2. compute the channels used by the two piconets determined by these parameters;

3. calculate the number of collisions between two piconets over the duration of connection simulated.

Two methods were used:

1. **First method:** it consists in calculating the number of collisions between a piconet using parameters (ULAP = 0, CLKN = 0) and a piconet using random parameters (ULAPᵣnd, CLKNᵣnd);

2. **Second method:** it consists in calculating the number of collisions between two random piconets using (ULAPᵣnd,1, CLKNᵣnd,1) and (ULAPᵣnd,2, CLKNᵣnd,2).

Each simulation considers a connection using only one slot packets, and assumes there is no drift between the two piconets clocks. The duration of
<figure><img src="image_74.png" /><figcaption>Fig. 6. Simulation of collisions between two Bluetooth piconets. First method.</figcaption></figure>

each connection is $T_{com} = 2h$. The time slot is $T_{slot} = 625$ µs. So each connection represents $N_{slot} = (T_{com})/(T_{slot}) = 11 520 000$ slots.

The results are given in Figure 6 (first method) and in Figure 7 (second method). To obtain these results, we performed about 660 000 simulations for the results presented in Figure 6 and 560 000 for those in Figure 7.

<figure><img src="image_75.png" /><figcaption>Fig. 7. Simulation of collisions between two Bluetooth piconets. Second method.</figcaption></figure>

Copyright © 2008 John Wiley & Sons, Ltd.

We define the ratio of collisions in per cent as

$$R_c = 100 \cdot \frac{N_c}{N_{Total}} \qquad (17)$$

where $N_c$ is the number of slots with collisions, $N_{Total} = N_{simul} \cdot N_{slot}$ is the total number of simulated slots, and $N_{simul}$ is the number of simulations in each case. We define also the percentage of simulations experimenting the ratio $R_c$ of collisions as

$$R_T = 100 \cdot \frac{N_{R_c}}{N_{simul}}$$

where $N_{R_c}$ is the number of simulations with a ratio of collisions of $R_c$.

The number of simulations performed seems small compared to the great possibilities of couple of frequency hopping sequences selection parameters (ULAP, CLKN) [3]. However, with the number of simulations of connections performed and the duration of each connection, the number of frequency hops simulated is about $500 000 \times N_{slot} = 5.76 \times 10^{12}$ which is significant and large compared with simulations performed in previous work (50 000 in Reference [3]). The results obtained confirm the relevance of our approach as we observe no significant modification of the statistical behavior of the collisions between 5000 and 10 000 simulations. Performing more simulations gives only more clearly defined and regular curves. As we are only interested in a tendency, these results are enough for that study.

As expected, there are statistically only a few per cent of couples of piconets which experience severe interferences showing that the Bluetooth Selection Kernel Algorithm works well. To see details of the repartition of collisions, the y-axis scale has been magnified, and the numbers $R_c(R_T = 0)$ which give the severe interferences percentage are not visible on the curves but mentioned under them.

However, one could argue that in real world perfectly slots aligned with no drift Bluetooth clocks do not exist. To take into account a possible drift between two Bluetooth clocks, we performed new series of simulations.

These simulations use the second method mentioned above. The channels used by two piconets during connection are calculated. Every 12 500 slots, we shifted one piconet by one slot with respect to the other. By doing so, we simulated a drift between the two piconets master's clocks of about 40 ppm. In each simulation, a connection with a $T_c = 2h$ is used for the two piconets. To obtain the curve presented in Figure 8,

Wirel. Commun. Mob. Comput. 2009; 9:673-683

DOI: 10.1002/wcm
<figure><img src="image_76.png" /><figcaption>Fig. 8. Simulation of collisions between two Bluetooth piconets. Simulation of drift.</figcaption></figure>

$N_s = 237 000$ simulations were performed (please pay attention to the permutation of X and Y axes). The same remark concerning the number of simulations and the fast convergence of the results as mentioned above also applies here.

The drift induces no more collisions, but they are more equally distributed over the couples of piconets. The drift stops the current series of collisions normally experimented by a couple of piconets. For instance, if we consider again the case presented in part 3, the series of collisions would normally span over much more than 12 500 slots. The drift in time stops it at 12 500 slots, breaking the series. The frequency channel is then calculated with different clocks for the two piconets. In the mean time, because of the drift, some collisions may occur in cases where they would normally not exist. The time drift averages the number of collisions over all possible piconet couples. When examining the data collected, we see that the average collision rate is about 1.2 percent, which is what we could expect according to the theoretical calculation based on simple probability [10]. But here, the variance of the collision rate is lower than for the simulations without drift.

The results of these simulations are summarized in Table I.

These simulations show that in general collisions between two Bluetooth piconets are not so important and that only a few packets may be lost this way. However, situations where the number of collision is

<table><caption>Table I. Simulation of collisions between two Bluetooth piconets (2 h connection).</caption><thead><tr><th>Simulation</th><th>Number of simulations</th><th>Average collision (percent)</th><th>Collision variance</th></tr></thead><tbody><tr><td>Method 1</td><td>660 948</td><td>1262</td><td>0.7074</td></tr><tr><td>Method 2</td><td>564 003</td><td>1265</td><td>0.6378</td></tr><tr><td>Method 2 with drift</td><td>237 001</td><td>1215</td><td>0.0009</td></tr></tbody></table>

so high that most packets are lost for quite a long time may occur, and must be avoided. In order to reduce collisions in these particular situations, it is necessary to study more deeply the Bluetooth frequency hopping algorithm.

## 5. Properties of the Bluetooth Selection Kernel in Connection Mode

The problem of collisions occurring between two or more close Bluetooth piconets will be considered from different points of view.

In most cases, the structure of the Bluetooth selection kernel allows piconets to coexist without too many collisions as shown previously. However, very critical situations may occur, leading to a poor link quality due to a very high number of collisions. Thus Bluetooth may be used in situations where performance is a very critical parameter, so packet collision which slows down the traffic, must be eliminated. Examples of such situations may include real-time wireless control systems operating in extreme conditions.

Taking this into account, we have devised a new way to improve coexistence by using some specific properties of the Bluetooth selection kernel.

Firstly, we will detail the Bluetooth selection kernel in connection mode, which is the default stationary mode for an active Bluetooth device. This study will confirm that, in most cases, Bluetooth piconets coexistence is insured. However, for some applications it may be critical to be sure that collisions will not happen, even if it concerns a small fraction of the transmitted packets. We will present a method to obtain a zero collision Bluetooth environment.

The structure of the Bluetooth selection kernel in connection mode is defined in Equations (3) to (11) for the configuration described in Figure 2.

We can show the influence of each input on the calculated channel.

Because of the modulo 32, $X$ varies every 625 µs time slots from 0 to 31 and then wraps around. The
inputs $Y_1$ and $Y_2$ permit the selection kernel to take into account two different channels for transmit ($Tx$) and receive ($Rx$) slots.

The $A$ input adds an offset to $X$. That offset has a slow variation in time compared to that of $X$ (every 65.536s). In these conditions, $Z'$ can also take values from 0 to 31.

$B$ performs a XOR on the 4 LSB of $Z'$ which is 5 bits wide. The values of $Z$ can be split into two groups according to the MSB of $Z'$. The most significant bit of $Z'$ is unchanged, thus for values of $Z'$ between 0 and 15, $Z$ also takes values from 0 to 15 but in a non-identical sequence with respect to $Z'$. Similarly, for values of $Z'$ between 16 and 31, $Z$ takes also values from 16 to 31, here again in a non-identical sequence.

$C$ and $D$ are used to initialize the PERM5 block which performs a permutation of all the bits of $Z$, leading to $Z_p$. We can derive that the probability of having $Z = Z_p$ is given by

$$ \mathrm{Prob}(Z = Z_p) = \frac{1}{C_n^p} = \frac{p!(n-p)!}{n!} \quad (18) $$

where $n=5$ is the number of inputs of the PERM5 block, and $p$ is the number of ones in the binary representation of $Z$.

The following step is the addition of $E$ and $F$. It just consists of adding to $Z_p$ a constant, $E$, and an offset $F$ varying each time $X$ wraps around.

Similarly, $D$ takes a new value each time $X$ wraps around. Thus, at such times we have a new permutation PERM5 block and a new final additive offset. It can be observed that when $\mathrm{CLK}_{16}$ toggles, $C$ takes a new value, and that when $\mathrm{CLK}_{21}$ toggles, $A$ takes a new value.

As $X$ changes every two clocks ticks and $Y_1$ changes at each ones, $32 \times 2 \times Y_1 = 64$, different values are calculated by the selection kernel for each period of $X$ and that the sequence is not repetitive over a long time.

Considering the Bluetooth selection kernel, we can now understand how collisions between Bluetooth devices arise, therefore how to improve coexistence between Bluetooth piconets.

## 6. How to Improve Coexistence

Firstly, we will see the required conditions for coexistence of two Bluetooth piconets due to the Bluetooth selection kernel properties. The cases of 1, 3, or 5 slots packets are studied. Then we will extend

results to the coexistence of any number of Bluetooth piconets and show that a limitation to that number exists. Finally, we will draw the consequences of the above discussions for coexistence.

### 6.1. The Case of Two Bluetooth Piconets

We now consider two piconets with their selection kernel boxes with their inputs $X_i$, $Y_{1,i}$, $Y_{2,i}$, $A_i$, $B_i$, $C_i$, $D_i$, $E_i$, and $F_i$, and their outputs $C_{hi}$, where $i \in \{1, 2\}$. The masters of these piconets will be called $M_1$ and $M_2$.

In order to suppress interferences between the two piconets, we must use two different channels for every new slot. One can easily see that there are two realistic and deterministic ways to satisfy this condition: either use the same address for the two devices, and shift their clocks by a judiciously chosen offset, or use the same clock for the two devices, but choose relevant individual addresses.

The first point is not to be considered because it leads to difficult problems while trying to establish a connection. Indeed, a device which tries to connect to one of the two devices described here will see one or the other with the same address and different clocks. Such a situation is confusing and should be avoided. Indeed, one can imagine the following situation, where a slave $\mathcal{S}$ connects to master $M_1$. Problems arise when considering authentication and encryption which use master's BD-ADDRs among other parameters. Suppose that $\mathcal{S}$ connects to $M_1$, requests for authentication [11], and then disconnects. After that, suppose that $\mathcal{S}$ tries again to connect using the same address as previously but connects to $M_2$ instead of $M_1$. This time authentication will fail, and the end-user will not be able to understand why. Even if the PIN (Personal Identification Number) codes are the same for $M_1$ and $M_2$, the connection establishment will fail as there is little chance that the two masters could generate the same random number used in the Bluetooth authentication process. For such reasons, the first solution should be avoided.

Let us consider the second point.

If we adopt a mathematical point of view, $\mathrm{SK}_1$ and $\mathrm{SK}_2$ produce identical outputs if they have identical inputs. Note that this may happen even if the inputs are not identical, since the selection kernel is a surjective function. We do not consider here this last particular case: whereas if it is easy to obtain zero collision for a certain duration, it is difficult to have such a situation for all possible values of the input clock.
By controlling the address ULAP<sub>31-0</sub>, we can control inputs B, C, D, E. We can now discuss the influence of the address on these inputs and on the calculated channel.

For that we will successively consider the different inputs assuming that for each choice, all the other inputs of the selection kernels SK<sub>1</sub> and SK<sub>2</sub> are identical.

1. **B:** as seen previously, Z is made of the MSB of Z' and of the 4 LSB of Z', xor-ed with B. Thus, we can write $Z_i = Z'_i(4) | (Z'_i(3:0) \oplus B_i(3:0))$, $i \in \{0, 1\}$. We see that for two different values $B_1$ and $B_2$, the calculated $Z_1$ and $Z_2$ will be different. To summarize this result: $B_1 \neq B_2 \Rightarrow Z_1 \neq Z_2$, so that, the calculated channels $Ch_1$ and $Ch_2$ are different;

2. **C, D:** these inputs define the effect of the permutation. Since, it is difficult to analytically represent their influence, we analyze this permutation block from a probabilistic point of view. As the result given by Equation (18) is also valid here. One can write: $Prob(Z_{p,1} = Z_{p,2} | Z_1 = Z_2 \text{ and } C_1D_1 \neq C_2D_2) = (C_n^p)^{-1}$, where $n=5$ is the number of inputs of the PERM5 block, and $p$ is the number of ones in the binary representation of Z. The probability is not equal to zero and thus we cannot be sure that the calculated channels will be different even if inputs $C_i$ and/or $D_i$ are different;

3. **E:** the effect of this input is a pure addition. The output channels are equal if and only if $Z_{c,1} = Z_{c,2}$. In order to analyze such a possibility, we calculate $Z_{c,1}$ and $Z_{c,2}$, reminding that SK<sub>1</sub> and SK<sub>2</sub> have all their inputs equal, except for the E input. We have

$$Z_{c,1} = E_1 + F_1 + 32 \times Y_{2,1} + Z_{p,1}$$

$$Z_{c,2} = E_2 + F_2 + 32 \times Y_{2,2} + Z_{p,2}$$

We see that by choosing $E_1 = E_2 \mod 79$, that is to say $E_1 = E_2 + n \times 79$, where $n$ is an integer ($n \in \mathbb{Z}$), then $Z_{c,1} = Z_{c,2}$. Since E is coded with 7 bits, it can only take values from 0 to 127. By combining $E_1 - E_2 = n \times 79$ and $0 \le E_i \le 127$ we can write $-127 \le E_1 - E_2 \le 127$, and thus $-127 \le n \times 79 \le 127$ and $-1, 61 \le n \le 1, 61$. Because $n \in \mathbb{Z}$, this last equation is equivalent to $-1 \le n \le 1$, that is to say $n \in \{-1; 0; 1\}$. We can conclude that we have different calculated channels if and only if $E_1$ and $E_2$ are not equal to each other modulo 79. More precisely, if all the inputs of SK<sub>1</sub> and SK<sub>2</sub> are identical, the calculated channels will be the same if and only if $E_1 = E_2$ or $|E_1 - E_2| = 79$.

<figure><img src="image_77.png" /><figcaption>Fig. 9. Illustration of the use of multiple slot packets in piconets.</figcaption></figure>

To summarize these results: if we use two selection kernels with all identical inputs except B or E we will obtain different output channels. If we use two selection kernels with all identical inputs except C and/or D, then the probability of obtaining the same calculated channels is not equal to zero.

As a conclusion, the way to control the selection kernel in order to obtain different channels is to have different inputs B or E (with respect to the modulo 79 condition), with the same Bluetooth clocks. To have a simple method to suppress piconets interferences, we choose to work on a single parameter (E or B). As E is coded with more bits than B, E gives more possibilities to define adequate addresses for piconets. We decide to work only on the E parameter.

## 6.2. Multi-Slot Packets

We now consider multi-slot packets. In this case, we may face many situations, which we discuss now, only considering the case of two coexisting piconets (Master 1, Slave 1) and (Master 2, Slave 2). The reasoning developed here may be extended with no difficulty to the case of multiple coexisting piconets.

Figure 9 summarizes the various situations we may face when using multiple slot packets in one or both the piconets considered:

1. **Situation A:** A multi-slot packet is exchanged in one piconet whereas only one slot packets are exchanged in the other (dotted packets). In this case, the reception of one slot packets overlapping with the emission of the multi-slot one will be jammed. So, we cannot consider these one slot packets for coexistence issues. However, the one slot packets exchanged during the end of multi-slot emission and following reception are not likely to be jammed and may be considered further (plain packets).
2. **Situation B:** this situation represents a case when the two piconets exchange similar multi-slot packet simultaneously. There is no particular problem with regard to coexistence in this case.

3. **Situation C:** this last situation may be encountered when one piconet's master decide not to send a packet at some slots. We are here in a similar situation as described in Situation A.

According to our previous condition stating that the two piconets considered are using the same values of Bluetooth Clock, in Situation B, we are in the same case as for one slot packets from the point of view of calculated RF channels.

This is not the case in Situations A and C. Here, in Situation C, for example, we have to consider RF channels calculated with different values of Bluetooth clock. So, we have to check that:

$$ Ch(ULAP_1, CLK_1(slot 9)) \neq Ch(ULAP_2, CLK_2(slot 11)) \quad (19) $$

that is to say, given that CLK₁ = CLK₂ for each slot:

$$ Ch(ULAP_1, CLK_1) \neq Ch(ULAP_2, (CLK_1 + 4)) \quad (20) $$

It is not possible to force that condition all the time, even with the requirement derived in the previous paragraph. Nevertheless, we have performed simulations of the Bluetooth frequency hopping algorithm with a proper choice of the Bluetooth addresses for two and more piconets (up to five piconets). These simulations show 3 per cent of collisions for two piconets, 6 per cent for three piconets, about 8 per cent for four piconets, and 10 per cent for five piconets. Without specific choice of the Bluetooth addresses, the percentage of collision is about 21, 25, 28, and 31 per cent, respectively. So our solution improves the coexistence between piconets.

In fact with different slot lengths exchanged, the problem is not the collision in emission but the non-reception by a piconet during the transmission of the others since they are very close to each other. The co-channel and adjacent channel rejections are not enough to avoid jamming in that case. The same difficulties occur when one piconet master must scan its environment without being perturbed by the communication of the other piconets. To solve these problems, a supervisor is implemented jointly with the Bluetooth interfaces to control the emissions. The

algorithm used in the supervisor is similar in its principle to the algorithm proposed by Chiasserini [1,7] which inhibits emission of one Bluetooth interface when it experiments collision in a non-synchronized network. Our supervisor uses the same inhibition process in the case of multi-slot packets or during scanning phases but in a synchronized network. Thanks to the synchronization, the slot length and the scan mode informations provided, the supervisor is able to inhibit transmissions only when necessary without experimenting any collision. In our implementation, the supervisor inhibits the transmission of a slot of one piconet master if the corresponding reception slot is during the emission of a multi-slot of the other piconets. The supervisor also gives the reference clock to the all Bluetooth interfaces in the multiple Bluetooth access point we have developed.

As we want to present the advantage of choosing adequately the parameters (ULAP and CLKN) for each piconet master used in a multiple Bluetooth access point and not the functionality of the supervisor, we will consider here only situations such as Situation B.

So, with that simplification, we are not faced with time delayed channel collision problems and the previous analysis concerning one slot packets is still valid.

## 6.3. Extension to n Bluetooth Piconets

We can easily see that it is possible to extend the above principle to n Bluetooth selection kernels. We only consider the effect of E as mentioned previously at the end of the Subsection 6.1. The above condition concerning the modulo 79 E values is immediately applicable to the case of n selection kernels.

Since there are only 128 possible values for E, we are limited in our coexistence study to at most 128 simultaneous active Bluetooth piconets. We know that in general if $E_i \neq E_j$, with $i \neq j$ we have $Ch_i \neq Ch_j$. As seen previously, we have $Ch_i = Ch_j$ when $E_i = E_j + n \times 79$ with $n \in \{-1; 0; 1\}$. Thus, if we limit our system to at most 79 selection kernels, each having a different $E_i$ and by ensuring that for all $(i, j)$ we have $|E_i - E_j| \neq 79$, there will be not collision at all, even for multi-slot packets.

We will now analyze the consequences of the above observations on the coexistence of Bluetooth piconets.

## 6.4. Consequences for Coexistence

In the previous section, we have shown that by taking identical clocks for the various piconets and by
making a particular choice for the device addresses we eliminate the risk of collision.

However, such a synchronization of the various devices is not always possible. For instance, if we consider a first piconet composed of a headset and of a GSM terminal, and a second piconet composed of an access point and a laptop, there is no easy way to realize the synchronization between the two masters (the GSM terminal and the network access point).

But there are many cases where this is possible. Let us consider for instance two printers close to each other, which can be linked through a Bluetooth connection. In this example, we assume that the printers are configured to act as masters. If this was not originally the case, they could obtain this status through a master-slave role switch procedure. By using the previous results, we see that synchronizing the clocks of the two piconets will prevent collisions, thus providing better overall throughput and decreasing the need for a long connection time.

In order to make the Bluetooth devices no more complex, other independent processes of synchronization can be implemented such as a GPS synchronization or a single distributed clock like in our prototype. In these cases, the clock is given by a specific circuitry (with correction of the propagation time) without introducing significant algorithm complexity.

These procedures may be applied to different devices and particularly to very important ones such as LAN or Dial-Up Networking access points. One can in such a way build a “network” of access points, keeping them synchronized and thus avoiding collisions.

Similarly, if Bluetooth is to be used in situations as described in paragraph 5 (monitoring real-time systems), it is also possible to synchronize many piconet masters together.

## 6.5. A Real Implementation

In the solution proposed to reduce collisions one can note that, from Equation (1), the bits of interest in *B* and *E* are the LAP bits, thus they can be chosen by the manufacturer. So, it is not only a theoretical approach.

The method has been really implemented with no modification of the Bluetooth standard in a prototype including three Bluetooth interfaces in a single access point. The Bluetooth interfaces are 3 cm far from each others and controlled by a supervisor ensuring synchronization of clocks. The prototype has been designed to allow a great number

of Bluetooth connections simultaneously in order to develop wireless LAN access point minimizing the number of located bases and the number of synchronization links between them.

To validate our prototype, some measurements were performed and compared to a single commercial Bluetooth access point [9]. These measurements were made in an office environment including some Wi-Fi interferers. Connections were realized between our Bluetooth multiple access point and other commercial Bluetooth devices such as USB dongles for laptop.

Some typical results obtained with the prototype are:

* connection of more than seven slaves simultaneously is possible;

* when only two slaves are connected, the data throughput is about $55 \text{ kByte} \cdot \text{s}^{-1}$ each, which is 40 percent larger than if the two slaves were connected to a classic Bluetooth access point;

* when seven slaves are connected, the average data throughput is about $17 \text{ kByte} \cdot \text{s}^{-1}$ each, which is 30 percent larger than if the seven slaves were connected to a classic Bluetooth access point;

* when nine slaves are connected, the data throughput is about $13 \text{ kByte} \cdot \text{s}^{-1}$ each.

The time measured to establish a new connection is about 3 s even with six connections already running. All these measurements were performed 200 times to ensure a reproducibility.

In the case of two slaves connected, one might expect a doubling of the data throughput. This is not the case because all the tests were performed by sending random length data packets from the access point to each Bluetooth client. In this case, the situation described in the Subsection 6.2 applies. Master 1 (or 2) can send a multi-slot (3 or 5) packet whereas master 2 (or 1) only sends a 1-slot packet.

Additionally, Bluetooth master devices shall switch periodically into one of the following scan modes:

* inquiry scan mode, which allows others Bluetooth devices to discover the master,

* page scan mode, which allows others Bluetooth devices to connect to the master.

When a master is in scan mode, it cannot send data to its connected slaves, which results in a lower throughput. Moreover, in our implementation, the supervisor allows only one master to be in scan mode at a time and suspends the traffic of others masters.
This strategy has two goals:

* It allows only one master to be discoverable and connectable at a time

* It prevents the scanning master to be jammed by the traffic of adjacent masters.

All these considerations explain the reduction of the data throughput measured.

Nevertheless, the results presented show that as a conclusion, it is possible to build a network of access points, for instance in a public environment (airport, train station, etc....) or in a private environment (small office, home, etc....) in which many Bluetooth devices are synchronized together using our method.

## 7. Conclusion

After a precise study of the Bluetooth selection kernel, we have proposed in this paper a specific use of the Channel Selection algorithm in the Bluetooth system which makes it possible to have piconet coexistence without losing packets. This is achieved by synchronizing Bluetooth clocks and by properly choosing the Bluetooth addresses of the piconets masters.

Assuming a system consisting of many masters fulfilling this condition, we can obtain an interference-free Bluetooth environment. A prototype of a multiple LAN access points was realized providing a better data flow than three access points close to each other.

## References

1. IEEE 802.15 WPAN Task Group 2. IEEE 820.15 tg2 technical presentations. Technical Report, IEEE, 1.1 edn, 2000-2001. http://www.ieee802.org/15/pub/tg2.html

2. Bluetooth Special Interest Group. Technical Report, Bluetooth SIG, Bluetooth SIG.

3. Haartsen JC, Zürbes S. Frequency hop selection in the Bluetooth radio system. IEEE Seventh International Symposium on Spread Spectrum Techniques and Applications, 2002.

4. Howitt I. IEEE 802.11 and Bluetooth coexistence analysis methodology. IEEE Vehicular Technology Conference, 2001; 2: 1114-1118.

5. Bluetooth SIG Coexistence Working Group. AFH Specification. Bluetooth SIG, 0.7 edn, November 2002.

6. Specification of the Bluetooth System Core. Bluetooth SIG, 1.2 edn, January 2003.

7. Rao RR, Chiasserini CF. Coexistence mechanisms for interference mitigation between IEEE802.11 WLANs and Bluetooth. IEEE INFOCOM, 2002.

8. Zürbes S. Considerations on Link and System Throughput of Bluetooth Networks. Personal, Indoor and Mobile Radio Communications, 2000; 2: 1315-1319.

9. Villain A. Etude et réalisation d'un système de points d'accès optimisés en technologie Bluetooth. PhD Thesis, Université Pierre et Marie Curie (Paris 6), 7 July 2004.

10. Haartsen J. The bluetooth radio system. IEEE Personal Communications 2000; 7(1): 28-36.

11. Specification of the Bluetooth System, Core. Bluetooth SIG, 1.1 edn, February 2001.

12. Haartsen JC, Matisson S. Bluetooth—a new low-power radio interface providing short-range connectivity. Proceedings of IEEE 2000; 88(10): 1651-1661.

## Authors' Biographies

<figure><img src="image_78.png" /></figure>

Emmanuel Géron, graduated from the engineering school École Supérieure de Physique et de Chimie Industrielles (E.S.P.C.I) of the City of Paris in 1992 and obtained his Ph.D. in Electronics from the University of Paris 6 in 1997. He is Assistant Professor of the E.S.P.C.I since 1997. He has made many contributions to Telecommunications.

<figure><img src="image_79.png" /></figure>

Alban Villain, graduated from the engineering school École Supérieure de Physique et de Chimie Industrielles (E.S.P.C.I) of the City of Paris in 1999, obtained his Ph.D. in Electronics from the University of Paris 6 in 2004. He has been an engineer in Inventel, a telecommunication research company, from September 1999 to March 2005. He

is an R&D engineer for the Thomson Telecom company since April 2005.
<figure><img src="image_80.png" /></figure>

ISSN 1751-8725

# Fast antenna pattern synthesis using the variational method

G. Collin<fnref n="1" /> E. Géron<fnref n="1" /> J. Lucas<fnref n="1" /> T. Ditchi<fnref n="1" />,<fnref n="2" />,<fnref n="3" />

E-mail: thierry.ditchi@espci.fr

**Abstract:** The radiation pattern is usually a factory set characteristic of telecommunication equipment. By using an antenna array, the radiation pattern can be shaped in order to optimise the covered region. Antenna pattern synthesis consists of determining the phase and the amplitude to apply to each radiating element of the network for obtaining the targeted radiation pattern. A radiation pattern synthesis technique based on the variational method is presented. The goal is expressed in terms of antenna radiation pattern and special attention is given to the computation time. We calculate an optimised mask function from specifications obtained during a surroundings learning step. The resulting radiation pattern is optimised in terms of power waste. Simulated and measured radiation patterns are presented in order to illustrate the performances of the method. The speed of convergence makes this method well suited to systems where radiation patterns must be modified in real time.

## 1 Introduction

The radiation pattern is most of the time a factory set characteristic of a telecommunication base station. As a consequence, the covered region cannot be optimised. Optimising the radiation pattern toward the users reduces the radiated power in useless directions, increases communication range and leads to a greater number of users per base station. Furthermore, minimising the gain of the system toward parasitic sources reduces their influence and consequently improves the quality of the communication.

By using an antenna array in which the phase and the amplitude of each source ('complex weights') can be controlled, the radiation pattern can be shaped nearly at will. The main issue in antenna pattern synthesis consists of determining the phase and the amplitude to apply to each radiating element of the network in order to obtain the targeted radiation pattern. Because of the rocketing development of numerical radio communications, antenna pattern synthesis is being thoroughly investigated.

Mathematical studies of antenna arrays have been published such as in [1] but only numerical processing can

treat complex patterns such as those containing several main lobes for instance. Many techniques have been developed to increase synthesis precision and speed. They use analytical calculus completed by numerical optimisations [2] for instance in the case of non-uniform side lobes [3]. Some optimisation methods use a neuronal network approach [4]. Specific and well-defined field radiation patterns are employed as training configurations for the learning step of the neuronal networks. Other methods optimise the geometrical positioning of the emitting sources in the antenna array leading to a non-uniform final array [5-8]. More complex optimisation methods combine this geometrical approach with the classical complex weights modification, coupled with the implementation of non-identical sources in the network [8]. Among these methods the optimisation of the radiated near field by physically placing zeros near the antenna in order to improve the far field pattern is probably the most original [9]. These geometric methods are not suitable for smart antennas whenever the emitted field pattern has to be modified dynamically according to the electromagnetic environment. A few methods propose to minimise the couplings between sources in the array [10] during the determination of the complex weights in order to optimise

<footnote n="1">¹ESPCI ParisTech, L.E.G., 10 rue Vauquelin, 75005 Paris, France</footnote>

<footnote n="2">²UPMC Univ. Paris 06, UPR 5, LPEM, Paris, France</footnote>

<footnote n="3">³CNRS, UPR 5, LPEM, Paris, France</footnote>
the radiated pattern. All previous methods concern narrowband signals. Calculations are developed considering only the carrier frequency of the signal. However, wideband signals have been studied. A specific method has been proposed to optimise the far-field pattern for such wideband signals [11]. Mangenot and Judasz [12] present a short description of the variational method applied to a two-dimensional antenna array for a satellite. A synthesis example illustrates the performance of the method with a single main lobe and concludes that the variational method is very fast. These papers mostly consider only one main lobe. Our application implies more elaborate specifications such as several main lobes, as well as a null radiated field in some given directions.

In our work, the antenna shape and the distribution of the radiating elements are fixed. The number of radiating elements has been chosen to allow complex radiation patterns while remaining realistic in terms of electronic complexity and antenna size. The synthesis is performed considering two points. The goal is expressed in terms of antenna radiation pattern and the time allocated for the synthesis must be very short to be compliant with the available time of most telecommunication protocols. Among the most widely used methods, the variational method seems to be a good compromise considering the previous points.

In this paper, we present a complete theoretical development of a variational-based method. In this method, a mask function is built from the specifications. The influence of this mask shape on the performance of the method is also discussed. Some synthesis involving elaborated specifications are presented as an illustration. The numerical results obtained are then compared with the radiation pattern of a real antenna array measured in an anechoic chamber. A discussion of algorithm speed and stability concludes.

## 2 Mask function

In an adaptive system, pattern synthesis takes place after a learning stage of the surroundings. In this stage, the system identifies and localises the surrounding electromagnetic sources and sorts them as users or parasitic emissions. It yields a set of triplets ($G_i, \theta_i, \varphi_i$), where $G_i$ is the gain the system must present in each direction ($\theta_i, \varphi_i$) as defined in Fig. 1. This set will be referred to in what follows as the learned specifications. Since a set of phases and amplitudes ($A_n, \psi_n$) satisfying exactly the learned specifications does not exist, an optimisation method is required. The learned specifications must be previously converted in a mask function as presented in Fig. 1. A mask function consists of two functions $G_{max}(\theta, \varphi)$ and $G_{min}(\theta, \varphi)$, respectively, representing the maximal and minimal limits the antenna gain must respect. Fig. 1 presents a typical mask function in a given plane. In this figure, the characteristic parameters of a mask function are presented. A is the amplitude tolerance on the gain in the direction of the users, $A_{2nd}$

<figure><img src="image_81.png" /><figcaption>Figure 1 Typical mask function in the case of one-dimensional synthesis</figcaption></figure>

defines the minimal attenuation in directions where no source is detected, $L_{max}$ is the maximal width of the mask function in the direction of users, $L_{min}$ is the minimal width of the mask function in these directions, $G_{zero}$ is the maximal value of the antenna gain in the direction of the parasitic emissions and $L_{zero}$ is the minimal width of the mask function in these directions.

Before studying the effects of parameters ($L_{max}$, $L_{min}$, A, $A_{2nd}$) on the convergence speed and on the validity of the results, the next section presents the principle of the variational method.

## 3 Variational method

Let us consider the total radiated electrical field at great distance from an antenna network of combining N radiating elements. The polarisation of the electrical field radiated by each element is taken along vector e

$$E_{\text{total}}(\theta, \varphi) = \sum_{n=1}^{N} [a_n e^{j(2\pi/\lambda)OM_n \cdot U} E_n(\theta, \varphi)]e \quad (1)$$

where $a_n = A_n e^{j\psi_n}$ is the complex weights applied to the radiated element n. $OM_n$ is the position of element n relatively to the Cartesian coordinates system (Oxyz) as presented in Fig. 2a. U is the vector pointing in the direction ($\theta, \varphi$) of observation and $E_n(\theta, \varphi)$ is the field radiated by element n.

In order to determine the best complex vector ($a_n$) leading to a radiated field as close as possible to the goal, the variational method minimises the error function J presented in (2)

$$J = \int \int_{(\theta, \varphi) \in D_1} \frac{K(\theta, \varphi) + |K(\theta, \varphi)|}{2} d\theta d\varphi \\ + \int \int_{(\theta, \varphi) \in D_2} |E_{\text{total}}|^2 d\theta d\varphi \quad (2)$$
<figure><img src="image_82.png" /><figcaption>Figure 2 Prototype of antenna</figcaption></figure>

a Coordinates system

b Realised array antenna on double side FR4 substrate (resin epoxy reinforced with a woven fiberglass)

The domain $D_1$ consists of the directions in which the mask function is defined. The domain $D_2$ represents the directions where no mask function is defined. One can note in Fig. 1, that the mask function is defined for all directions. The domain $D_2$ is thus empty in this case.

In (2), K is a function that measures the difference between the modulus of the total radiated field and the boundaries of the mask function

$$ K(\theta, \varphi) = [|E|_{\text{Max}}(\theta, \varphi) - |E|_{\text{total}}(\theta, \varphi)|] \times [|E|_{\text{Min}}(\theta, \varphi) - |E|_{\text{total}}(\theta, \varphi)|] \quad (3) $$

where

$$ |E|_{\text{Max}} = 10^{G_{\text{Max}}/20} \quad \text{and} \quad |E|_{\text{Min}} = 10^{G_{\text{Min}}/20} \quad (4) $$

The error function J presents two terms. The first one depends on the distance in the domain $D_1$ between the total radiated field and the mask function. The second one is proportional to the total radiated power in the domain $D_2$.

One can notice that in the directions where $|E|_{\text{total}}$ is between the two limits $|E|_{\text{Max}}$ and $|E|_{\text{Min}}$ K is negative, hence $K + |K|$ is null. However, when $|E|_{\text{total}}$ is outside this template, K is positive and $K + |K|$ is greater than zero.

In the variational method, a small variation $\delta a_n$ of vector $a_n$ is considered. If a small variation $\delta a_n$ leads to a null variation $\delta J$ of the function J, then the value of vector $(a_n)$ corresponds to a minimum of the J function.

The calculation of $\delta J$ presented in Appendix leads to

$$ \delta J = \Re \left\{ \sum_{n=1}^{N} \delta_{n}^{*} \left[ \iint_{D_1} [-E_{n}^{*} e^{-j(2\pi/\lambda)(OM_{n}U)}] \times E_{\text{total}} \left[ \left( \frac{|E|_{\text{Max}} + |E|_{\text{Min}}}{2|E_{\text{total}}|} - 1 \right)(1+\sigma) + 2 \right] d\theta d\varphi + 2 \iint_{D_1+D_2} [E_{n}^{*} e^{-j(2\pi/\lambda)(OM_{n}U)}] E_{\text{total}} d\theta d\varphi \right] \right\} \quad (5) $$

As $\delta J$ must be null for any variation $\delta a_n$ one must have

$$ \iint_{D_1} [E_n^* e^{-j(2\pi/\lambda)(OM_n U)}] \times E_{\text{total}} \left[ \left( \frac{|E|_{\text{Max}} + |E|_{\text{Min}}}{2|E_{\text{total}}|} - 1 \right) \left( \frac{1+\sigma}{2} \right) + 1 \right] d\theta d\varphi = \iint_{D_1+D_2} [E_n^* e^{-j(2\pi/\lambda)(OM_n U)}] E_{\text{total}} d\theta d\varphi \quad (6) $$

in which $E_{\text{total}}$ can be replaced by $[\sum_{m=1}^{N} a_m e^{j(2\pi/\lambda)(OM_m U)} E_m]$

It yields a non-linear system of equations with N variables $a_n$. N is the number of radiating elements. This system can be written as

$$ b_n = \sum_{m=1}^{N} a_m C_{m,n} \quad (7) $$

where

$$ \begin{aligned} b_n &= \iint_{D_1} [E_n^* e^{-j(2\pi/\lambda)(OM_n U)}] E_{\text{total}} \\ & \times \left[ \left( \frac{|E|_{\text{Max}} + |E|_{\text{Min}}}{2|E_{\text{total}}|} - 1 \right) \left( \frac{1+\sigma}{2} \right) + 1 \right] d\theta d\varphi \\ C_{m,n} &= \iint_{D_1+D_2} e^{j(2\pi/\lambda)(OM_m - OM_n)U} E_n^* E_m d\theta d\varphi \end{aligned} \quad (8) $$

and $a_m$ is the searched complex weight applied to the radiated element m.

The variational method consists in numerically solving (7). The resolution is performed iteratively starting from an initial vector $a_m(0)$ as presented in the flowchart of Fig. 3. The matrix C only depends on the network antenna geometry and is calculated once and for all at the beginning of the algorithm.

Iterations stop upon two conditions: when the maximal number of iterations is reached and when the error is small enough. The maximal number of iteration is given by Maximal_Time/Iteration_Time.
<figure><img src="image_83.png" /><figcaption>Figure 3 Variational method flowchart</figcaption></figure>

‘Maximal\_Time’ is the maximal duration allowed for the calculus. It depends on the telecommunication protocol used. ‘Iteration\_Time’ is the duration of each iteration that depends on the speed of the processor used. As an iteration involves about 10 thousands multiplications, and considering an instruction cycle of a few nanoseconds that is the typical value for a medium performances signal processor, it yields to an iteration time smaller than 0.1 ms. In what follows, we have chosen to use a maximum of 50 iterations, which leads to a calculation time of a few milliseconds.

## 4 Influence of the characteristic parameters of the mask function

Let us now study the influence of characteristic parameters A, Lmax and A2nd on the synthesised radiation pattern and on the speed of convergence. These parameters allow us to define the mask function in all directions from the required antenna gain defined for discrete directions. One must take care to choose parameters that are consistent with the physical capabilities of the antenna network used. A bad choice could prevent the algorithm to find the best solution or could lead to slow convergence.

We consider now a linear antenna network made of eight radiating elements. The pattern synthesis is restricted to the half horizontal space ($\varphi = \pi/2, \theta \in [0; \pi]$, Fig. 2a).

The first parameter A allows the limit functions $G_{Max}$ and $G_{Min}$ of the mask function to be calculated from $G_i^{max}$, the antenna gain values retained in the learning stage in the direction $\theta_i^{max}$ of user i. One has

$$G_{Max}(\theta_i^{max}) = G_i^{max} + A \quad \text{and} \quad G_{Min}(\theta_i^{max}) = G_i^{max} - A$$

As presented in Fig. 1 parameter Lmax accounts for the maximal lobe width in the directions of the users. As the

narrowest half power beam width it is possible to obtain with our array is 15°, we have restricted our study to values of Lmax ranging from 20° up to 35°. The last characteristic parameter A2nd accounts for the minimal attenuation required in all directions {$\theta_{sec}$} other than the directions of the users or the parasitic sources. The side lobes level is defined relative to the maximum of $G_{Max}$ and is given by $G_{Max}(\theta_{sec}) = Max(G_{Max}) - A_{2nd}$.

The influence of some parameters, such as $G_{zero}$ the maximal value of the gain in the direction of parasitic emission, will not be considered here. $G_{zero}$ depends on each parasitic emission power and is chosen in the learning step. The influence of $L_{zero}$, the angular width around the parasitic emissions, and the influence of the minimal angular width $L_{min}$ in the direction of the maximum will not be considered either. Both parameters are fixed at 10°. Lastly, the minimal gain $G_{min}$ in all directions not listed in the learning specification {$\theta_{sec}$} is taken null.

It is necessary to define a criterion in order to quantify the error resulting from the synthesis. We have defined two criteria that only depend on the difference between the synthesised electrical field and the learned mask function. They are calculated over two different angular domains. The firsts one, $\epsilon_{max}$, measures the error on the radiated field in the direction of the main lobes (users directions: $\theta \in \{\theta_i^{max}\}$). The second one, $\epsilon_{sec}$, measures the power radiated in all others directions $\theta \in \{\theta_{sec}\}$. They are defined using normalised fields and presented in (9) and (10)

$$\epsilon_{max} = 100 \sum_{\theta \in \{\theta_i^{max}\}} \frac{||E(\theta)| - |E_{specification}(\theta)||}{|E_{specification}(\theta)|} (\%) \quad (9)$$

$$\epsilon_{sec} = \sum_{\theta \in \{\theta_{sec}\}} |E(\theta)|^2 \quad (10)$$

In order to illustrate the influence of the parameters Lmax, A and A2nd let us use representative sets of learned specifications. The first set corresponds to the case of two users situated at 90° and 130° at different distances. The distances are such that the gain in the direction of the second user must be 1.4 dB smaller than the gain of the first user in order to obtain the same signal level for both users.

### 4.1 Parameter Lmax

As presented in Fig. 4, parameter Lmax affects the side lobes level and the error in the direction of the users. A value of 35° for this parameter allows obtaining a better power value in the direction of the users and a smaller side lobes level than with a value of 20°. Nevertheless, it leads to a greater power waste because the main lobes are wider than necessary. Finally, even if the side lobes level is greater when using a value of 20° for parameter Lmax, the global useless power remains low because the main lobes are better fitted.
<figure><img src="image_84.png" /><figcaption>Figure 4 Influence of parameter $L_{max}$ on the radiated pattern</figcaption></figure>

## 4.2 Parameter A

Parameter $A$ represents the tolerance in the direction of the main lobes. It affects the side lobes as well as the error on the direction of the synthesised pattern main lobes as presented in Fig. 5. Increasing $A$ leads to a bigger error on the main lobes direction meanwhile it reduces the side lobes level.

## 4.3 Parameter $A_{2nd}$

Parameter $A_{2nd}$ represents the difference between the maximum gain and the side lobes level. It affects, of course, the side lobes level, but it also affects the width of the main

<figure><img src="image_85.png" /><figcaption>Figure 5 Influence of parameter $A$ on the radiated pattern</figcaption></figure>

lobes. As a consequence, the error on the gain in the direction of the users increases with $A_{2nd}$. Considering the results presented in Fig. 6, using a value ranging from 18 to 20 dB is a good compromise between power loss and gain value in the direction of users.

## 4.4 Optimal values

Finally, optimal values that have been determined by studying a great number of typical specifications are $L_{Max} = 25^{\circ}$, $A = 10\%$ and $A_{2nd} = 20$ dB.

# 5 Results and measurements

In this section, we firstly discuss the obtained radiation patterns with a real antenna. Secondly, the speed of convergence and stability of the algorithm are specifically discussed. These are indeed some of the key points that make such an algorithm practically usable.

## 5.1 Suitability of obtained radiation patterns

The variational method has been implemented in a system made of an antenna with eight radiating elements and its associated control electronics. Radiating elements are positioned every half wavelength ($\lambda/2$) and are composed of three rectangular patch antennas fed in phase. The realised array antenna is presented in Fig. 2b. The radiation pattern of each element of the array antenna noted $E_n$ in the previous calculation has been measured in an anechoic chamber. It has then been used to calculate the radiation pattern of the whole antenna using (1). The coupling between elements was not taken into account at this point. The optimal parameter values determined in the preceding

<figure><img src="image_86.png" /><figcaption>Figure 6 Influence of parameter $A_{2nd}$ on the radiated pattern</figcaption></figure>
paragraph have been used to build up the mask functions. Figs. 7a and b present two examples of calculated pattern synthesis in the xOy plane (i.e. $\varphi = 0$), in light grey dotted curve, resulting from two different sets of specifications. In Fig. 7a, the users are present at bearings 85° and 115° in the xOy plane ($\varphi = 0$) and at different distances. The specifications set the gain in those respective directions to 0 and -4 dB. In Fig. 7b, users are present at bearings 60°, 90° and 120° with respective gain values of -4, 0 and -2.5 dB. In both cases, the synthesis yields patterns presenting correct gain values in the user's directions, when the radiated power in the other direction is attenuated by 20 dB. In both cases the synthesised patterns have been measured (dark grey curve) and are compared to the previously calculated patterns. As for the calculated ones, the desired gain is obtained, but the side lobes level can be as high as -10 dB at worst. This can be explained by considering the coupling between radiating elements. This coupling modifies indeed the effective amplitude radiated by each element. The result does not invalidate the synthesis method because mutual coupling can be taken into account if a greater precision is required.

## 5.2 Speed of convergence and stability of the method

In order to evaluate the speed of convergence of the method and to verify that the algorithm cannot be trapped in some local minimum of the error function, for many

<figure><img src="image_87.png" /><figcaption>**Figure 7** Examples of pattern synthesis</figcaption></figure>

a Two users at bearings 85° and 115°, in the xOz plane ($\varphi = 0°$)

b Three users at bearings 60°, 90° and 120°, in the xOz plane ($\varphi = 0°$)

specifications we have performed a great number of syntheses with randomly created initial complex weights.

As far as the local minima of the error function $\delta J$ are concerned, we have found that whatever the initial complex weights the synthesised pattern is always nearly the same. In Fig. 8 one can see as an illustration that for given learned specifications, with a thousand sets of randomly generated complex weights, the obtained error $\delta J$ at the end of the calculus is the same. One can notice nevertheless that the convergence is noticeably slower in some of the presented cases. In those cases, the algorithm has encountered a weak slope region that does not correspond anyhow to a local minimal value. Fig. 9 presents the corresponding radiation patterns (in light grey lines) that are here again nearly identical.

If the algorithm is stopped after 50 iterations as proposed in Section 3, the obtained radiation pattern is then different

<figure><img src="image_88.png" /><figcaption>**Figure 8** Error as a function of the number of iterations for 1000 randomly generated complex weights</figcaption></figure>

<figure><img src="image_89.png" /><figcaption>**Figure 9** Radiation patterns obtained after synthesis for 1000 randomly generated complex weights</figcaption></figure>

Light grey line: After 200 iterations. Dark grey dotted line: After 50 iterations, in the case of slow convergence
<figure><img src="image_90.png" /><figcaption>**Figure 10** Error as a function of the number of iterations for 200 randomly generated specifications with initial complex weight amplitudes [1 000 000 0] and phases [0 000 000 0]</figcaption></figure>

from those after 200 iterations. Nevertheless, even in the case of slow convergence, the obtained radiation pattern as presented in dark grey dotted line in Fig. 9 remains sufficiently close to the specifications to be used as is.

Furthermore we have found that when the initial complex weights are set to [1 000 000 0] for the amplitudes and zero for all phases, the previously described phenomena is no longer encountered. This result is illustrated in Fig. 10, which presents the *δJ* obtained whilst this set of initial weights for 200 random specifications and does not present any slow convergence curve.

Finally, oscillations or divergence of the error function with iteration have not been observed even with very elaborate specifications such as multiple main lobes and zeros. The convergence of the method appears quasi monotonous as shown in Figs. 8 and 10.

## 6 Conclusion

In this paper we have presented a radiation pattern synthesis technique based on the variational method. We have shown how to define a mask function from specifications obtained, for instance, from a learning stage of the surroundings. We have shown as well that it is efficient to use optimised characteristic parameters for the mask function in order to obtain radiation patterns fitting the specifications and optimised in terms of power waste. As those parameters depend on the kind of antenna, the determination of the optimal characteristic parameters has to be performed only once for each kind of antenna. Calculated as well as measured radiation patterns are presented in order to illustrate the performances of the method. The speed of convergence is an additional benefit of the presented

method. It makes it a good candidate for use in a system where the radiation pattern must be modified in real time.

## 7 References

[1] SCHELKUNOFF S.A.: ‘A mathematical theory of linear array’, *Bell Syst. Tech. J.*, 1943, **22**, p. 80

[2] SHAVIT R., LEVY S.: ‘A new approach to the Orchard-Elliott pattern synthesis algorithm using LMS and pseudoinverse techniques’, *Microw. Opt. Technol. Lett.*, 2001, **30**, (1), pp. 12–15

[3] ELLIOTT R.S., STERN G.J.: ‘A new technique for shaped beam synthesis of equispaced arrays’, *IEEE Trans. Antennas Propag.*, 1984, **AP-32**, (10), pp. 1129–1133

[4] SHAVIT R., TAIG I.: ‘Comparison study of pattern-synthesis techniques using neural networks’, *Microwa. Opt. Technol. Lett.*, 2004, **42**, (2), pp. 175–179

[5] KHODIER M.M., CHRISTODOULOU G.G.: ‘Linear array geometry synthesis with minimum sidelobe level and null control using particle’, *IEEE Trans. Antennas Propag.*, 2005, **53**, (8), pp. 2674–2679

[6] AKDAGLI A., GUNEY K.: ‘Shaped-beam pattern synthesis of equally and unequally spaced linear antenna arrays using a modified Tabu search algorithm’, *Microw. Opt. Technol. Lett.*, 2003, **36**, (1), pp. 16–20

[7] NG B.P., ER M.H., KOT C.: ‘Linear array geometry synthesis with minimum sidelobe level and null control’, *IEE Proc. Microw. Antennas Propag.*, 1994, **141**, (3), pp. 162–166

[8] OLEN C.A., COMPTON R.T.: ‘A numerical pattern synthesis algorithm for arrays’, *IEEE Trans. Antennas Propag.*, 1990, **38**, (10), pp. 1666–1676

[9] STEYSKAL H.: ‘Synthesis of antenna patterns with imposed near-field nulls’, *Electron. Lett.*, 1994, **30**, (24), pp. 2000–2001

[10] DARWOOD P., FLETCHER P.N., HILTON G.S.: ‘Retrodirective beams and their application to low sidelobe pattern synthesis in small planar phased arrays’, *IEE Proc. Microw. Antennas Propag.*, 1998, **145**, (4), pp. 344–348

[11] GUELLA T.R., DAVIS R.M.: ‘Synthesis of notched patterns for wideband processing’, *IEEE Trans. Antennas Propag.*, 1995, **43**, (12), pp. 1465–1471

[12] MANGENOT C., JUDASZ T.: ‘Power synthesis of shaped beam antenna patterns’. AP-S Digest on Antennas and Propagation Society Int. Symp., June 1989, vol. 1, pp. 420–423
# 8 Appendix

## 8.1 Calculation of the variation $\delta J$ of the error function $J$ as defined in (2)

Let us calculate $\delta J$ by calculating $J + \delta J$. For clarity sake, $K(\theta, \varphi)$ and $E(\theta, \varphi)$ will, respectively, denote $K$ and $E$. More generally the angle dependency will be omitted in following equations.

$J + \delta J$ can be written in the following way

$$J + \delta J = \iint_{(\theta, \varphi) \in D1} \frac{K + \delta K + |K + \delta K|}{2} d\theta d\varphi + \iint_{(\theta, \varphi) \in D2} |E_{\text{total}} + \delta E_{\text{total}}|^2 d\theta d\varphi \quad (11)$$

Furthermore, one has

$$K + \delta K = [|E|_{\text{Max}} - |E_{\text{total}} + \delta E_{\text{total}}|] \times [|E|_{\text{Min}} - |E_{\text{total}} + \delta E_{\text{total}}|] \quad (12)$$

which can be developed in

$$K + \delta K = |E|_{\text{Max}}|E|_{\text{Min}} - |E|_{\text{Max}}|E_{\text{total}}| + \delta E_{\text{total}}| - |E|_{\text{Min}}|E_{\text{total}}| + \delta E_{\text{total}}| + |E|_{\text{total}} + \delta E_{\text{total}}|^2 \quad (13)$$

Noting $A^*$ the conjugate of $A$, one has

$$|E_{\text{total}} + \delta E_{\text{total}}|^2 = (E_{\text{total}} + \delta E_{\text{total}})(E_{\text{total}} + \delta E_{\text{total}})^* \\ = |E_{\text{total}}|^2 + |\delta E_{\text{total}}|^2 + 2\Re(e(E_{\text{total}}\delta E_{\text{total}}^*) \quad (14)$$

The conditions $|\delta K| \ll |K|$ and $|\delta E_{\text{total}}| \ll |E_{\text{total}}|$ must be realised to ensure the variational method convergence.

So one can write at first order

$$|E_{\text{total}} + \delta E_{\text{total}}|^2 \approx |E_{\text{total}}|^2 + 2\Re(e(E_{\text{total}}\delta E_{\text{total}}^*) \quad (15)$$

thus

$$|E_{\text{total}} + \delta E_{\text{total}}| \approx (|E_{\text{total}}|^2 + 2\Re(e(E_{\text{total}}\delta E_{\text{total}}^*))^{1/2} \\ \approx |E_{\text{total}}| \left( 1 + \frac{2}{|E_{\text{total}}|^2} \Re(e(E_{\text{total}}\delta E_{\text{total}}^*) \right)^{1/2} \quad (16)$$

which yields once again at first order

$$|E_{\text{total}} + \delta E_{\text{total}}| \approx |E_{\text{total}}| \left( 1 + \frac{1}{|E_{\text{total}}|^2} \Re(e(E_{\text{total}}\delta E_{\text{total}}^*) \right) \\ \approx |E_{\text{total}}| + \frac{1}{|E_{\text{total}}|} \Re(e(E_{\text{total}}\delta E_{\text{total}}^*) \quad (17)$$

Finally

$$K + \delta K \approx (|E|_{\text{Max}} - |E_{\text{total}}|)(|E|_{\text{Min}} - |E_{\text{total}}|) - \frac{1}{|E_{\text{total}}|} \Re(e(E_{\text{total}}\delta E_{\text{total}}^*) \\ \times (|E|_{\text{Max}} + |E|_{\text{Min}} - 2|E_{\text{total}}|) \quad (18)$$

It is then possible to calculate $|K + \delta K|$ considering that

$$|K + \delta K| = \left| K \left( 1 + \frac{\delta K}{K} \right) \right| = |K| \left( 1 + \frac{\delta K}{K} \right) \\ = |K| + \sigma \delta K \quad (19)$$

where $\sigma = |K|/K = \text{sgn}(K)$.

We have

$$|K + \delta K| \approx |(|E|_{\text{Max}} - |E_{\text{total}}|)(|E|_{\text{Min}} - |E_{\text{total}}|)| \\ - \frac{\sigma}{|E_{\text{total}}|} \Re(e(E_{\text{total}}\delta E_{\text{total}}^*) \\ \times (|E|_{\text{Max}} + |E|_{\text{Min}} - 2|E_{\text{total}}|)) \quad (20)$$

Using (16), (18) and (20), one can write $J + \delta J$ in the following way

$$J + \delta J = \frac{1}{2} \iint_{D_1} K + |K| d\theta d\varphi + \iint_{D_2} |E_{\text{total}}|^2 d\theta d\varphi \\ + \iint_{D_2} 2\Re(e(E_{\text{total}}\delta E_{\text{total}}^*)) d\theta d\varphi \\ - \iint_{D_1} \frac{1 + \sigma}{2|E_{\text{total}}|} \Re(e(E_{\text{total}}\delta E_{\text{total}}^*)) \\ \times (|E|_{\text{Max}} + |E|_{\text{Min}} - 2|E_{\text{total}}|) d\theta d\varphi \quad (21)$$

By identifying with the expression of $J$ given by (5) one can deduce $\delta J$

$$\delta J = \iint_{D_2} 2 \Re(e(E_{\text{total}} \delta E_{\text{total}}^*)) d\theta d\varphi \\ - \iint_{D_1} \frac{1 + \sigma}{2|E_{\text{total}}|} \Re(e(E_{\text{total}} \delta E_{\text{total}}^*)) \\ \times (|E|_{\text{Max}} + |E|_{\text{Min}} - 2|E_{\text{total}}|) d\theta d\varphi \quad (22)$$
which can be rewritten as

$$ \delta J = \Re e \left[ \int \int_{D_1+D_2} 2 E_{\text{total}} \delta E_{\text{total}}^* d\theta d\varphi \right] \\ - \Re e \left[ \int \int_{D_1} E_{\text{total}} \delta E_{\text{total}}^* \left[ \left( \frac{|E|_{\text{Max}} + |E|_{\text{Min}}}{2|E|_{\text{total}}} - 1 \right) \times (1+\sigma) + 2 \right] d\theta d\varphi \right] \quad (23) $$

As $\delta E_{\text{total}}^* = \sum_{n=1}^{N} \delta a_n^* e^{-j(2\pi/\lambda)(OM_n U)} E_n^*$ it yields

$$ \delta J = \Re e \left\{ \sum_{n=1}^{N} \delta a_n^* \left[ \int \int_{D_1} -[E_n^* e^{-j(2\pi/\lambda)(OM_n U)}] \times E_{\text{total}} \left[ \left( \frac{|E|_{\text{Max}} + |E|_{\text{Min}}}{2|E|_{\text{total}}} - 1 \right) (1+\sigma) + 2 \right] d\theta d\varphi \right. \\ \left. + 2 \int \int_{D_1+D_2} [E_n^* e^{-j(2\pi/\lambda)(OM_n U)}] E_{\text{total}} d\theta d\varphi \right\} \quad (24) $$
<figure><img src="image_91.png" /></figure>

ISSN 1751-8725

# Electronically controlled asymmetric microstrip line coupler underpinned by an hybrid right- / left-handed line

E. Géron<fnref n="1" /> T. Ditchi<fnref n="2" /> J. Lucas<fnref n="1" /> S. Holé<fnref n="2" />

E-mail: emmanuel.geron@espci.fr

**Abstract:** This work presents a microstrip coupler based on an hybrid right-/left-handed (RH/LH) line whose coupling ratio is electronically controlled. After some considerations on the LH line theory, the asymmetric RH/LH coupler is presented focusing on the parameters yielding the coupling ratio and the frequency range. The dependence of the coupling ratio on the number of cells constituting the coupler is highlighted. From these results, a prototype of a coupler with an electronically switched ratio has been developed. The variable coupling ratio is obtained by switching ON or OFF some cells in the structure of the coupler. The isolation of the OFF mode of the electronic switch implemented is carefully studied. Measurements and calculus are presented in order to estimate the isolation factor required to obtain a significant variation of the coupling ratio as well as to illustrate the original behaviour of such couplers. A functional electronically controlled coupler with one switch implemented is finally characterised and its performances are discussed.

## 1 Introduction

Since the 1990s, the feasibility of materials presenting simultaneously a negative permeability and permittivity has been demonstrated in built up macro-structures named metamaterials. A negative phase velocity for the electromagnetic waves is observed in such structures. Consequently, the wave vector **k**, the electrical field vector **E** and the magnetic field vector **H** constitute an indirect trihedron. Such metamaterials are called left-handed materials (LH) in opposition to the standard right-handed propagation (RH) [1]. The singularity of this negative phase velocity has triggered a great interest and a lot of works have been carried out lately to develop these metamaterials in one-, two- and three-dimensional (1D, 2D and 3D) structures.

Whereas 2D and 3D structures are interesting for optical and more generally for applications where a radiated phenomenon is involved [2, 3], 1D structures such as hybrid transmission lines are potentiality relevant for microwave circuit applications [4, 5]. For instance, symmetrical coupling structures using two LH coupled lines [6] or asymmetrical structures using one LH line and one RH line coupled together [7] present interesting characteristics to develop new circuits for telecommunications.

In a great number of applications both a large coupling bandwidth and a high coupling ratio are required. At high frequencies these electromagnetic couplers are often designed by closing up two RH lines. The frequency range

of such couplers is limited since it is directly dependent on the length of the two lines in regard. In this case, the coupling ratio between the input port and the coupled port can only be adjusted via the gap between the two lines [8] which does not allow variable coupling ratio.

This work presents a microstrip coupler based on an hybrid RH/LH line whose coupling ratio is electronically controlled. After some considerations on the LH line theory, the asymmetric RH/LH coupler is presented focusing on the parameters yielding the coupling ratio and the frequency range. We highlight the dependence of the coupling ratio on the number of cells constituting the coupler. From these results, a prototype of a coupler with an electronically switched ratio has been developed. The variable coupling ratio is obtained by switching ON or OFF some cells in the structure of the coupler. The isolation of the OFF mode of the electronic switch implemented is carefully studied. Measurements and calculus are presented in order to estimate the isolation factor required to obtain a significant variation of the coupling ratio as well as to illustrate the original behaviour of such couplers. A functional electronically controlled coupler with one switch implemented is finally characterised and its performances are discussed.

## 2 Theoretical approach of the LH/RH coupler

A transmission line in which the role of the capacitor and of the inductor have been swapped compared to the elementary model of a conventional RH transmission line behaves as a

<footnote n="1">¹Laboratoire d'Electricité Général, ESPCI-ParisTech, Paris, France</footnote>

<footnote n="2">²Laboratoire de Physique et d'Etude des Matériaux, UPMC Univ Paris 06, ESPCI-ParisTech, CNRS UMR8213 10, rue Vauquelin 75005, Paris, France</footnote>
LH line. As these LH lines do not exist naturally, they must be artificially realised using lumped capacitors $C_L$ connected in series with lumped inductors $L_L$ shorting each connection to the ground. The realisation of such a line requires to use short RH lines to connect the lumped components. At high frequencies the propagation through these RH lines cannot be neglected [9] which makes it an hybrid line.

The transmission line used in our prototype is an hybrid RH/LH line. It consist in a sequence of the elementary cells presented Fig. 1 connected in series. Each cell is composed of two capacitors $2C_L = 2.4$ pF serially connected, separated by a self $L_L = 5.1$ nH connected in parallel to the ground. The two RH line sections used to connect these components are $d_R/2 = 2.35$ mm long. The distributed components values of the RH lines are $C_R = 82$ pF m⁻¹ and $L_R = 410$ nH m⁻¹.

Assuming that the wavelength in the line is largely greater than the cell length, and replacing the line segments with their equivalent capacitances and inductances per length unit, one can determine the dispersion relation without losses [9]

$$ \cos(\theta) = 1 - \left( (L_R d_R)\omega - \frac{1}{C_L \omega} \right) \left( (C_R d_R)\omega - \frac{1}{L_L \omega} \right) \quad (1) $$

In this relation $\theta$ is the phase shift between each cell of the hybrid line during propagation for an infinitely long line.

With the condition $-1 < \cos(\theta) < 1$, one obtains the three following characteristic frequencies

$$ \begin{aligned} f_{c1} &= \frac{1}{4\pi\sqrt{L_L C_L}} \\ f_{c2} &= \frac{1}{2\pi\sqrt{C_R L_L d_R}} \\ f_{c3} &= \frac{1}{2\pi\sqrt{L_R C_L d_R}} \end{aligned} \quad (2) $$

Observing that $\theta$ is real between $f_{c1}$ and $f_{c2}$, and over $f_{c3}$ and purely imaginary in the other cases, propagation occurs between $f_{c1}$ and $f_{c2}$ and over $f_{c3}$ as well. The waves are evanescent otherwise. As the initial phase shift per cell is $-\pi$ before $f_{c1}$, the line is LH between $f_{c1}$ and $f_{c2}$. The behaviour of the line considering the losses and the propagation in the RH line segments can be obtained by simulation. It exhibits the typical transmission response which is presented Fig. 2. The result presented in this figure was obtained as all other simulations presented in this paper using the ADS circuit simulator software from Agilent Technologies. Simulation allows to tune the characteristic impedance of the line segment $Z_0 = \sqrt{L_R/C_R}$ to 71 Ω

<figure><img src="image_92.png" /><figcaption>Fig. 1 Hybrid RH/LH transmission line elementary cell</figcaption></figure>

<figure><img src="image_93.png" /><figcaption>Fig. 2 Simulated $|S_{21}|$ of the hybrid RH/LH transmission line</figcaption></figure>

which yields an overall characteristic impedance of about 50 Ω for the hybrid line.

The targeted applications are telecommunication systems in the 2.4 GHz ISM band. Consequently the central frequency of the coupling range was chosen at $f_{coupling} = 2.45$ GHz. The circuit board has been realised with an ordinary FR4 substrate ($\epsilon_r = 4.3$, $h = 0.8$ mm, copper thickness 35 µm, $\tan(\delta) = 0.02$) where the guided wavelength of a RH microstrip line at the central frequency is $\lambda_{line} = (1/\sqrt{L_R C_R})/f = 70.34$ mm

The obtained RH/LH structure can be considered as distributed line for a length $d_R$ much smaller than $\lambda$ the wavelength in the structure. We have $d_R < (\lambda/14)$ at the central frequency ($f_{coupling}$) which makes it possible to consider the RH/LH line as an artificial homogeneous material (i.e. a metamaterial). It is consequently relevant to consider propagation phenomena in the hybrid line.

As expected from the theory, the RH behaviour is separated from the LH behaviour by a forbidden band-gap between $f_{c3}$ and $f_{c2}$ where the waves are evanescent. With the values selected for the components of the line, the three characteristic frequencies are $f_{c1} = 1$, $f_{c3} = 3.34$ and $f_{c2} = 3.63$ GHz. Hence the central working frequency is in the LH band.

The coupler presented in this work takes advantage of the LH properties of the hybrid line presented above. Practically, it is realised as in a classical coupler, by taking the hybrid line close to a standard RH line. The study of the coupling mechanism [10, 11] shows that it only occurs via the RH line segment of the hybrid line and not through the lumped components, the coupler can then be modelled as a succession of elementary cells as for the hybrid line by adding two short coupler as presented in Fig. 3. In this

<figure><img src="image_94.png" /><figcaption>Fig. 3 Elementary hybrid coupler cell</figcaption></figure>
figure, one can notice that the length of the lumped capacitors has been taken into account by adding two small line segments of the same length whereas the total length of the two short couplers corresponds to that of the lines in regard. The width of each lines as well as the gap used are also indicated in this figure.

At this point, we have verified that the realised hybrid line still exhibits an LH behaviour even in the coupler structure. Fig. 4 presents the phase of the transmission coefficient $S_{32}$ of the LH line in the hybrid coupler, divided by the number of cells (blue curve), measured using an Agilent 8722ES network analyser. The measured phase in the forbidden band [0; 1 GHz] is no relevant because of the extremely small signal. This measurement has been backed up by calculating the $S_{32}$ argument from its amplitude with the technique based on Kramers–Kronig relations presented in [12] (dotted green curve). One can then verify that the line is LH up to about 3.5 GHz where the phase shift turns from a negative value to a positive one. For frequencies above 3.5 GHz, the line becomes RH.

Whereas the functioning of the classical coupler is well understood, the one of the LH coupler is still being investigated. 3D electromagnetic simulations of this kind of coupler seem to show that the coupling effect is the result of a vortex-like interface mode between the hybrid RH/LH and the classical RH line [13]. Anyhow, the coupling effect between the RH line and the hybrid line is noticeable when the hybrid line exhibits its LH characteristic [10]. For a given gap value, the hybrid coupler presents a greater power transfer between the two lines than for classical structures which allows to use larger gap values [6]. In addition as the coupling effect is controlled by the LH behaviour of the hybrid line, its frequency bandwidth is set by the characteristics frequencies controlled by the lumped capacitors and inductors of the elementary cell of the line. Finally the bandwidth of the coupler does not depend on the number of elementary cells of the structure but only on their length. In this work we take advantage of this property to build up a coupler whose coupling ratio is controlled by the number of cells involved without perturbing the bandwidth. The next section focuses on the relation between the number of cells and the coupling ratio.

<figure><img src="image_95.png" /><figcaption>**Fig. 4** Dispersion in the hybrid line of the asymmetric coupler</figcaption></figure>

<figure><img src="image_96.png" /><figcaption>**Fig. 5** Layout of the hybrid coupler composed of 10 cells</figcaption></figure>

### 3 Coupling ratio as a function of the number of cells

In order to study the coupling ratio as function of the number of cells, a simple hybrid RH/LH microstrip line coupler composed of 10 cells whose footprint is presented in Fig. 5 was realised. In [13] the input port is located on the LH line. In the coupler studied in this work, the input and transmit ports (ports 1 and 4) have been chosen on the RH microstrip line whereas the coupled (port 2) and the isolated (port 3) are on the LH/RH line as in [14]. With such a choice it is possible to connect or disconnect some of the last cells while maintaining the transmitted power way. All the parameters of the coupler: value of the components and dimensions can be found in Fig. 3.

In order to verify that the coupling ratio can be controlled by varying the number *N* of cells connected, simulations and measurements were performed for various values of *N*.

The simulations were performed with the structure of the coupler presented Fig. 5. Assuming that the coupling effect which depends on LH behaviour is localised in the small RH couplers, it is relevant to use a microwave circuits simulator like advanced design system (ADS) from Agilent Technologies instead of a 3D electromagnetic simulator. The results are as accurate but obtained with a shorter computational time. A good agreement between measurements and simulations validates this choice. Simulation results are presented Fig. 6.

In a complementary way, measurements were performed with 1 up to 10 cells by populating them with the lumped components one by one on the prototype whose footprint is presented in Fig. 5. The five cells measurement, for instance, was performed with all the capacitors and inductors of the five first cells soldered leaving the remaining component footprints empty. Fig. 7 presents the coupling ratio and the transmission ratio of the coupler against frequency. Measurements were performed using an HP8722ES network analyser. The measured bandwidth of the coupler is nearly as large as 500 MHz centered at 2.45 GHz. Such a large bandwidth is adequate even for ultra-wide band telecommunication applications.

Fig. 8 presents on the same graph, at the central frequency, the increase of the coupling ratio with the number of populated cells and its counterpart decrease of the transmission ratio. The graph presents an asymptotic behaviour for more than nine cells populated. The obtained coupling ratio goes from –15 (1 cell) up to –3 dB (10 cells).

By comparing Figs. 6 and 7, one can deduce that the hybrid coupler presents some power losses resulting from propagation along the 10 cells. This power losses are owing to the FR4 substrate which presents dielectric losses in the frequency band considered. This explains as well the asymptotic behaviour of Fig. 8, indeed after nine cells the power that has not been coupled has been simply lost and

<footnote>IET Microw. Antennas Propag., pp. 1–6
doi:10.1049/iet-map.2011.0426</footnote>
<figure><img src="image_97.png" /><figcaption>**Fig. 6** Coupling and transmission ratio simulated for 1–10 connected cells</figcaption></figure>

<figure><img src="image_98.png" /><figcaption>**Fig. 7** Measured coupling and transmission ratio for 1–10 populated cells</figcaption></figure>

can no longer be coupled whatever the LH behaviour of the hybrid line.

## 4 Electronically switched ON/OFF coupling ratio

To introduce some electronic switches in the structure to modify the number of coupled cells thus controlling the coupling ratio, a switch has to be chosen. The LH behaviour of the RH/LH line and the characteristic frequencies $f_c2$ and $f_c3$ are very sensitive to the length of the RH lines section. The lines required to integrate both the switch and its polarisation and control circuit must be as short as possible. New simulations must then be performed to revalidate the structure.

Various switches can be used, such as micro electromechanical systems (MEMS) switches or PIN diode switches. The isolation provided by the electronic switch must be high enough to well control the coupling ratio in the OFF mode. Simultaneously the losses in the switch in the ON mode must be as low as possible to transmit the power correctly. Recent MEMS switches exhibit promising performances with an isolation better than –25 dB, with very low insertion losses (about –0.25 dB) and long life cycles. Unfortunately most commercial MEMS, such as the ones from Radant MEMS company for example, require a high control voltage (±90 V) which makes them more difficult to implement than PIN diodes. Recent researches [15] show that lower control voltage will be available soon. Moreover the current commercial MEMS are provided as die which require bonding making them again less simple to use than diodes. Finally to quickly realise a prototype we have used a PIN diode as a switch.

The PIN switch has been inserted in the hybrid RH/LH line between the second and the third cell of the coupler as

presented in Fig. 9. When the switch is ON, all cells are connected and participate to the coupling process from port 1 to port 2. When the switch is OFF, only two cells are involved in the coupling.

<figure><img src="image_99.png" /><figcaption>**Fig. 8** Coupling and transmit ratio at 2.45 GHz as a function of the number of component populated cells</figcaption></figure>

<figure><img src="image_100.png" /><figcaption>**Fig. 9** One switch asymmetric controlled coupler</figcaption></figure>
A new prototype has been designed for the implementation of a PIN Diode switch with its polarisation circuit (Fig. 9). Simulations were performed to define the polarisation circuit and a PIN diode with a very small package was chosen in order to minimise the parasitic elements.

The polarisation circuit includes a shock self of 10 nF to connect the PIN diode anode to the ground and a resistor of 1.2 kΩ to connect its anode to the polarisation voltage. A shock self has been preferred to a λ/4 decoupling circuit to minimise the size of the whole polarisation circuit. According to the simulations, with this optimised design the C_L and L_L components and the length of the right line section of the hybrid line can remain unchanged without modifying the characteristics (coupling band and ratio) of the coupler.

According to the previous measurements the activation of the switch should take |S<sub>21</sub>| from -10 up to -3 dB. However the serial resistance of the PIN diode R<sub>OFF</sub> in OFF mode has a great influence on the behaviour of the coupler. Fig. 10 obtained by testing various resistors instead of PIN diode on the prototype, presents the measured coupling ratio as a function of the resistor value. For a resistance greater than 510 Ω, the coupling ratio is constant. The coupling effect is not no longer owing to the resistor value but to parasitic electromagnetic coupling. Consequently, this value of resistor is a good estimation of the minimal R<sub>OFF</sub> value required for the PIN diode.

The isolation S<sub>21</sub> provided by an impedance Z serially inserted in a loaded line of characteristic impedance Z<sub>0</sub> is given by (3).

$$S_{21} = \frac{2Z}{Z + 2Z_0} \qquad (3)$$

Assuming that the characteristic impedance of the hybrid line is approximatively 50 Ω, a resistor of 510 Ω in series leads to an isolation coefficient S<sub>21</sub> = -15, 6 dB. Consequently the PIN diode must present such an isolation at 2.45 GHz. The great majority of commercial RF PIN diodes exhibits an isolation around -10 dB above 2 GHz. The diode selected for the prototype is the NXP BAP50-02 which presents an isolation around -17 dB at 2.45 GHz.

With the layout of Fig. 9 and the diode mentioned previously, the prototype which can be seen in Fig. 11 reaches the performances presented in Fig. 12. The electronic

<figure><img src="image_101.png" /><figcaption>**Fig. 10** Measured coupling ratio as a function of the R<sub>OFF</sub> serial resistance</figcaption></figure>

<figure><img src="image_102.png" /><figcaption>**Fig. 11** Picture of the controlled asymmetric coupler</figcaption></figure>

<figure><img src="image_103.png" /><figcaption>**Fig. 12** Coupling ratio in the ON/OFF modes of the RF PIN diode switch</figcaption></figure>

switch controls the coupling ratio S<sub>21</sub> from -4 dB (diode ON) to -8 dB (diode OFF) which is in accordance with the previous measurements obtained by varying the number of cells (Fig. 7) with a resistor connected instead of the diode (Fig. 9). The slightly lower than expected coupling ratio in ON mode is due to the losses of the additional RH line sections. In the OFF mode the resulting coupling ratio is 1 dB greater than the expected one because of parasitic coupling via by the diode polarisation circuit.

The isolation S<sub>31</sub> has not been discussed so far. It is an important parameter. It has been measured in all presented situations in the paper and it is always better than -10 dB.

## 5 Conclusions

Using an hybrid LH/RH line is a very efficient way to realise couplers for telecommunication applications taking advantage of its LH behaviour. When well dimensioned, the coupling ratio of such structure can be modified by driving the number of cells of the coupler serially connected as it has been demonstrated in this work. A coupler with an electronically controlled coupling ratio has been realised and tested. Its performances are those expected from the initial analyse of its fixed coupling ratio counterpart. The realisation of such a coupler is highly dependent on the ability to maintain the LH behaviour of the LH/RH line once coupled and once the electronic
switch has been introduced. As discussed in this work, the isolation parameter of the switch implemented is a critical parameter. Most RF PIN Diode for instance do not fit this application. Nevertheless, by choosing the diode accordingly to the discussion and analyse presented in this paper we have realised a functional prototype.

The great overall length of such coupler that could be invoked as a limiting factor for its industrial use. Nevertheless for typical industrial application where the maximal coupling ratio can be less than -3 dB, a smaller number of cells is necessary which reduces drastically the length of the coupler. The length of such couplers is therefore not the blocking factor to their industrial use.

## 6 References

1. Veselago, V.G.: 'The electrodynamics of substances with simultaneously negative values of epsilon and mu', *Sov. Phys. Usp.*, 1968, **10**, (4), pp. 509–514

2. Coulombe, M., Nguyen, H.V., Caloz, C.: 'Substrate integrated artificial dielectric (SIAD) structure for miniaturized microstrip circuits', *IEEE Antennas Wirel. Propag. Lett.*, 2007, **6**, pp. 575–579

3. Pendry, J.B., Holden, A.J., Robbins, D.J.: 'Low frequency plasmons in thin wire structures', *J. Phys. Condens. Matter*, 1998, **10**, (4785), pp. 4785–4809

4. Zhang, D., Liu, F., Zhang, Y.: 'Experimental investigation of a power divider based on microstrip and metamaterials with L-C lumped-elements', *Microw. J.*, 2005, **48**, (11), pp. 114–117

5. Scher, A.D., Rodenbeck, C.T.: 'Compact ring resonators using negative-refractive-index microstrip line', *Microw. Opt. Technol. Lett.*, 2005, **45**, (4), pp. 294–295

6. Caloz, C., Sanada, A., Liu, L.: 'A broadband left-handed (LH) coupled-line backward coupler with arbitrary coupling level', *IEEE MTT-S Digest, Strony*, 2003, pp. 317–320

7. Caloz, C., Itoh, T.: 'A novel mixed conventional microstrip and composite right/left-handed backward-wave directional coupler with broadband and tight coupling characteristics', *IEEE Microw. Wirel. Compon. Lett.*, 2004, **14**, (1), pp. 31–33

8. Tsui, J.B.: 'Microwave receivers and related components' (Peninsula Publishing, 1985)

9. Lai, A., Caloz, C., Itoh, T.: 'Composite right/left-handed transmission line metamaterial', *IEEE Microw. Mag.*, 2004, **5**, (3), pp. 34–50

10. Wang, Y., Zhang, Y., He, L.: 'Discussion on coupling mechanism of asymmetric CRLH/RH TL coupler', *J. Zhejiang Univ. Sci. A*, 2006, **7**, (1), pp. 95–98

11. Wang, Y., Zhang, Y.: 'Direct observation of negative phase velocity and positive group velocity in time domain for composite right/left-handed transmission lines', *J. Appl. Phys.*, 2006, **100**, (113503), pp. 1–4

12. Lucas, J., Géron, E., Ditchi, T., Holé, S.: 'Practical use of the Kramers-Kronig relation at microwave frequencies, Application to photonic like lines and left handed materials', *PIERS 2001, Marrakesh, 22 Marzec 2011*, pp. 1003–1009

13. Wang, Y., Zhang, Y., He, L.: 'Time-domain study of vortexlike interface mode in metamaterials', *Appl. Phys. Lett.*, 2007, **91**, (221907), pp. 1–3

14. Engheta, N., Ziolkowski, R.W.: 'Metamaterials physics and engineering explorations' (Wiley-Interscience, 2006), pp. 198–204

15. Hoc-Chul, L., Jae-Hyoung, P., Hyo-Jin, N., Jon-Uk, B.: 'Design, fabrication and RF performances of two different types of piezoelectrically actuated Ohmic MEMS switches', *J. Micromech. Microeng.*, 2005, **15**, pp. 2098–2104
## 4.5 A fast Fourier transform implementation of the Kramers-Kronig relations: Application to anomalous and left handed propagation

AIP ADVANCES 2, 032144 (2012)

# A fast Fourier transform implementation of the Kramers-Kronig relations: Application to anomalous and left handed propagation

Jérôme Lucas,ª Emmanuel Géron, Thierry Ditchi, and Stéphane Holé

Laboratoire de Physique et d'Étude des Matériaux, UPMC Univ Paris 06, ESPCI-ParisTech,
CNRS UMR8213, France

(Received 22 May 2012; accepted 10 August 2012; published online 17 August 2012)

In this paper we quickly derive the Kramers-Kronig relations from simple causal-
ity considerations and propose a simple way to implement them using the Fast Fourier
Transform. This work focuses on how to make these relations a usable tool even when
their conditions of validity are not strictly respected. In this respect we emphasize
on their application to the constant low level loss approximation at microwave fre-
quencies. The method presented is demonstrated on various typical cases of fancy
propagation: high velocity, negative phase velocity and evanescent waves. Copyright
2012 Author(s). This article is distributed under a Creative Commons Attribution 3.0
Unported License. [http://dx.doi.org/10.1063/1.4747813]

I. INTRODUCTION

Exotic propagation phenomena called anomalous dispersion at optic frequencies have been reported long ago. Those phenomena often imply superluminal group velocity such as in photonic crystals. Though Brillouin and Sommerfield¹ early stated the conditions required for the group velocity to be the velocity of the energy, there are still numerous papers presenting and discussing superluminal group velocities,²⁻⁴ More recently left handed material lines or structures introduced by Velesago⁵ and Pendry⁶ came into focus. They are typically media where some anomalous dispersion occurs: there are bands where the phase velocity is negative and forbidden bands where the waves are evanescent. The thread binding those phenomena is the relation which links the phase shift to the amplitude of the signal. When the phase shift is due to propagation one has to consider the dispersion relation along with the amplitude. At lower frequencies when no propagation is involved it is referred to as the frequency response. The Kramers Kronigs relations⁷⁻⁹ are an interesting theoretical tool to meet those considerations. They translate to the frequency domain the simple temporal paradigm that the effect cannot occur before the cause. By doing so they relate the phase shift to the amplitude.

The Kramers-Kronig relations are commonly used in optics,¹⁰ acoustic¹¹ or at terahertz¹² frequencies to verify the consistency of performed measurements or to extrapolate for instance reflection coefficient measurements at frequencies where they are difficult to perform. Nevertheless the methods employed are not directly presented or restricted to a given kind of application. For instance, Kendall¹⁰ presents a method without explicitly giving the algorithm used for the calculation and applies it to relate the attenuation coefficient to the phase velocity of acoustical waves. In the same way, left handed propagation and evanescent waves are not considered.

On the basis of previous work,¹³ we present how to apply the Kramers-Kronig relations even when they cannot be straightforwardly applied because of limited bandwidth or because of Fourier transform convergence problem for instance. Its applications to S parameter measurements or to data extrapolation are discussed. The method is then practically applied to some interesting typical cases

<footnote>ªElectronic address: jerome.lucas@espci.fr</footnote>
such as photonic crystals and right and left handed media, where the waves can be superluminal, evanescent or left handed.

## II. THEORETICAL BACKGROUND

### A. Kramers-Kronig relations

Whatever the introductory phenomenon presented (photonic crystal, left handed material and so on) they occur in real thus causal media. Consequently the media response $r(t)$ to any excitation beginning at time $t = 0$ verifies:

$$r(t) = u(t)r(t) \qquad (1)$$

where $u(t)$ denotes the unity step function.

By taking the Fourier transform of this relation, and considering Cauchy principal value, one obtains:

$$R(\omega) = -\frac{2j}{\omega} * R(\omega) = j\mathcal{H}_T(R(\omega)) \qquad (2)$$

In this Equation, $j = \sqrt{-1}$, $\omega$ is the circular frequency, $R(\omega)$ is the Fourier transform of $r(t)$, $\mathcal{H}_T$ the Hilbert transform and $*$ denotes the convolution product. Note that the minus sign in this relation may change accordingly to the definition used for the Fourier Transform. In this paper we use the sign that is relevant with the FFT algorithm. By presenting separately the imaginary part $R_i$ and the real part $R_r$ of $R(\omega)$, one has:

$$\begin{aligned} R_i(\omega) &= \mathcal{H}_T(R_r(\omega)) \\ R_r(\omega) &= -\mathcal{H}_T(R_i(\omega)) \end{aligned} \qquad (3)$$

These relations are known as the Kramers-Kronig (K-K) relations. In a linear and time invariant media, Equation (1) also applies to the impulsive response $h(t)$, whose Fourier transform $H(\omega)$ is the frequency response. Consequently a frequency response $H(\omega)$ respects the pair of relations (3).

In a very general way, a relation such as (4), valid over all frequencies between two physical quantities as the electrical field $\vec{E}(\omega)$ and the displacement field $\vec{D}(\omega)$ where $\varepsilon(\omega)$ is the complex dielectric permittivity, describes a linear system.

$$\vec{D}(\omega) = \varepsilon(\omega) \vec{E}(\omega) \qquad (4)$$

Indeed with the inverse Fourier transform, it comes $\vec{D}(t) = \varepsilon(t) * \vec{E}(t)$, Permittivity $\varepsilon(t)$ can be interpreted as an impulsive response and consequently $\varepsilon(\omega)$ verifies the K-K relations (2) and (3).

### B. Practical computation of the Hilbert Transform

Let us consider a function $f(x)$ of variable $x$ whose Fourier transform is $\hat{f}(k)$ where $k$ is the corresponding variable in the Fourier space. Since the Hilbert transform of $f(x)$ is the convolution product with $-2/x$, it can be calculated as follows:

$$\begin{aligned} \mathcal{H}(f(x)) &= -\frac{2}{x} * f(x) \\ &= \mathcal{F}^{-1}\left(\mathcal{F}\left(-\frac{2}{x} * f(x)\right)\right) \\ &= \mathcal{F}^{-1}\left(j \operatorname{sgn}(k) \hat{f}(k)\right) \end{aligned} \qquad (5)$$
where $\mathcal{F}$ denotes the Fourier transform, $\mathcal{F}^{-1}$ the inverse Fourier transform and where the sgn($k$) function is +1 when $k \ge 0$ and -1 when $k < 0$.

The numerical Hilbert transform HT corresponding to Equation (5) can be computed for a given sampled signal $S(x_i)$ using the FFT algorithm:

$$HT(S) = IFFT(j \times SGN(k_i) \times FFT(S)) \qquad (6)$$

In this Equation SGN($k_i$) emulates the signum function, It is 1 for the samples corresponding to positive values of $k_i$ in the Fourier space and it is -1 for the others.

It is also possible with some precautions to use the direct convolution algorithm to compute the Hilbert transform. It does not yield better results than with the FFT algorithm and it is much slower.

## 1. Making things work

The K-K relation of Equation (2) implies an integration from $-\infty$ to $\infty$. It is indeed not possible numerically and, worse, no measurement can be performed over such a wide band. The FFT algorithm circumvent this problem by calculating the numerical Fourier transform of the periodized sampled signal S. This technique leads to “spectral leakage” when the original signal cannot be periodized with matching junctions.¹⁴ In our case we compute the Fourier transform $R(\omega)$ of some real signal $r(t)$ as defined in section II A. Consequently, $R(\omega)$ is Hermitian which can be written as in Equation (7), where $R^*$ is the conjugate of $R$.

$$R(-\omega) = R^*(\omega) \qquad (7)$$

This condition is equivalent to the even symmetry of the real part of $H(\omega)$ and to the odd symmetry of its imaginary part. Consequently the response $R(\omega)$ must be symmetrized respecting Equation (7) before applying the algorithm (6) to take Equation (2) into account where Cauchy’s principal value implies symmetry over the integration variable $\omega$.

As it will be described later in section III, the “spectral leakage” is dealt with by taming up the response $R(\omega)$ when required.

## C. Computing a phase shift from an amplitude measurement and vice-versa

This is the main practical application of the K-K relations. Indeed $R(\omega) = G(\omega)e^{j\theta(\omega)}$ is often determined by measuring its amplitude $G(\omega)$ and its argument $\theta(\omega)$. The K-K relation relates $G(\omega) \cos(\theta(\omega))$ to $G(\omega) \sin(\theta(\omega))$. But it is interesting from an experimental point of view to obtain $G(\omega)$ from the measurement of $\theta(\omega)$ and vice-versa. It is all the more interesting when one of the measurements is difficult or impossible to perform such as $\theta(\omega)$ when using a spectrum analyzer.

In order to separate $G(\omega)$ and $\theta(\omega)$, it is usual to consider the logarithm of $R(\omega)$: $\ln(R(\omega)) = \ln(G) + j\theta(\omega)$. It would be interesting if the following pair of relations could be applied:

$$\begin{cases} \ln(G(\omega)) = -\mathcal{H}_T(\theta(\omega)) \\ \theta(\omega) = \mathcal{H}_T(\ln(G(\omega))) \end{cases} \qquad (8)$$

This is unfortunately not the general case. As an example, let us consider a system consisting in a piece of ideal coaxial cable of length $L$. For the TEM mode, the group and phase velocities are equal and constant. Therefore the frequency response is: $H(\omega) = \exp(-j \frac{\omega}{v_{\phi}} L)$. Function $H(\omega)$ is the frequency response of a causal system since the transit time in such a system is $\frac{L}{v_{\phi}} \ge 0$.

Considering $\ln(H(\omega))$, one obtains:

$$\begin{cases} \ln(G(\omega)) = 0 \\ \theta(\omega) = -\frac{\omega}{v_{\phi}} L \end{cases} \qquad (9)$$

Therefore $\mathcal{H}(\ln(G(\omega))) = 0$ which is obviously different from $\theta(\omega)$. Consequently the Relations (8) are not respected. In this case $H(\omega)$ corresponds to a causal system and the K-K relations apply when they do not apply to $\ln(H(\omega))$.
One can find in the literature methods using Cauchy's integral theorem to demonstrate the K-K pairs of relations between ln (G(ω)) and θ(ω) of Equation (8) by adding some convergence condition to H(ω). For instance Wooten⁸ assuming that:

$$|H(z)| < |b \times z^{-s}| \qquad (10)$$

where b and s are real strictly positive constants, and z the analytical extension to the complex plane of ω with ℜ(z) = ω and ℑ(z) ≤ 0, proves that:

$$\theta(\omega) = -\frac{2\omega}{\pi} \mathcal{P} \int_{0}^{+\infty} \frac{\ln(G(u))}{\omega^2 - u^2} du \qquad (11)$$

In Equation (11), P stands for Cauchy's principal value.

Considering that G(ω) is an even function, it can be rewritten as: θ(ω) = ℌ$_T$ (ln(G(ω))). Calculation details can be found in appendix. Finally by taking the Hilbert transform of this relation one obtains the other Equation of the pair of relations (8).

The condition presented in Equation (10) is often invoked to justify the application of the K-K relations to the logarithm of the frequency response because it is a realistic physical condition for the response of a physical system. Nevertheless it seems to be a sufficient condition and not a necessary one. Let us consider a simple first order high pass filter (HPF) with a cutoff frequency of ωc:

$$H_{HPF}(\omega) = \frac{1}{1 - j \frac{\omega c}{\omega}} \qquad (12)$$

Nevertheless it corresponds to what is measured when the frequency remains sufficiently low. This system is causal because it respects the causality condition of Equation (1), but it does not respect the condition of Equation (10). However one obtains with a very good accuracy ln (G<sub>HPF</sub>(ω)) from θ<sub>HPF</sub>(ω) and vice-versa by simply applying Relations (8) using the algorithm based on the Fourier transform. That example shows that the condition (10) is a sufficient but not necessary one.

### III. APPLICATIONS

In this section we demonstrate how to successfully apply K-K relations to various examples chosen to cover most cases that can be encountered. These examples show how to deal with algorithm discrepancy. They emphasize on the case of simple propagation which in any case has to be dealt particularly. Finally most kind of wave propagation such as in photonic crystal like lines, evanescent waves, left handed waves are illustrated through application to measurements.

#### A. Second order low pass Filter:

Let us consider as an example a plain second order low pass filter of central frequency ω₀ and quality factor Q, whose frequency response is:

$$H_{LPF2}(\omega) = \frac{1}{(1 + j \frac{\omega}{Q\omega_0} - \frac{\omega^2}{\omega_0^2})} \qquad (13)$$

It is a causal system, and it respects the condition (10). However applying straightforwardly the K-K relation pair (8), does not yield the expected result as presented in Figure 1. In this figure one can compare ln (GLPF2(ω)) (1(a) dashed black) to - ℌ (θ<sub>LPF2</sub>(ω)) (1(a) solid red) and θ<sub>LPF2</sub>(ω) (1(b) dashed black) to ℌ (ln (G<sub>LPF2</sub>(ω))) (1(b) solid red) and be convinced of the difference. As condition (10) is respected, the discrepancy can only be a consequence of the algorithm used, which is to say in our case some “frequency leakage”.

In order to solve this problem one can use a gain compensation technique.

Let H₁(ω) be the frequency response of a given physical quantity. According to the case, whether the real part of ln (G₁(ω)) = ℜ<sub>c</sub> (ln (H₁(ω))) or the imaginary part of θ₁(ω) = ℑ<sub>m</sub> (ln (H₁(ω))) are known or can be measured between 0 rad/s and ω<sub>max</sub>. Let H₂(ω) be a frequency response such that the pair of Equations (8) applies.
<figure><img src="image_104.png" /><figcaption>FIG. 1. Application of the compensation gain method to $H_{LPF2}(\omega)$.</figcaption></figure>

By calculating the product $H(\omega) = H_1(\omega) \times H_2(\omega)$, one has $\theta(\omega) = \theta_1(\omega) + \theta_2(\omega)$ and $\ln(G(\omega)) = \ln(G_1(\omega)) + (G_2(\omega))$. To compute $\theta_1(\omega)$ from $\ln(G_1(\omega))$ for instance, we choose $H_2(\omega)$ such that: $\ln(G_1(\omega)) + \ln(G_2(\omega))$ is asymptotic to zero when $\omega \to \omega_{max}$. With this condition one reduces the “frequency leakage” problem when computing the FFT. It is simply sampled at $1/\omega_{max}$ rate. Consequently one can expect to compute the Hilbert transform with a minimal error. The phase $\theta_1(\omega)$ can then be retrieved from the computation of $\theta(\omega) = H_T(\ln(G_1(\omega)) + \ln(G_2(\omega)))$ by simply computing: $\theta_1(\omega) = \theta(\omega) - \theta_2(\omega)$.

This has been applied successfully using $H_{HPF}^2(\omega)$ as a compensation gain. The result is presented in Figure 1: “K-K Comp” plots using blue o markers. One can see it behaves as expected. Let us emphasize again that it is possible to use $H_{HPF}^2(\omega)$ as a compensating frequency response because $H_{HPF}(\omega)$ respects Equation (8). It is not possible to use any built up function of a complex variable for whose Hilbert transform may yield anything.

Of course there are many ways to respect the condition: “$\ln(G_1(\omega)) = \ln(G_2(\omega))$” is asymptotic to zero when $\omega \to \omega_{max}$” even when simply using $H_{HPF}(\omega)^2$ as a compensation gain. Indeed the cutoff frequency of $H_{HPF}(\omega)$ has to be chosen at best. Noticing that when there is “frequency leakage”, some imaginary part appears when computing the Hilbert transform, the compensation optimization can be performed by minimizing $\int_m (H_T(\ln(G(\omega))))$ as the Hilbert transform of a real function must be real. This is what has been performed in the example presented in Figure 1.

## B. Plain propagation in lossy FR4 media

In the case of plain propagation, such as the propagation of TEM waves in the piece of cable previously presented, the condition (10) does not apply. We have verified in section II C that the pair of relations (8) does not apply neither. Nevertheless, as propagation cannot be neglected when the frequency gets sufficiently high, it is necessary to handle the case.

### 1. A quiet general case: the Debye model

It is interesting at this point to shift from an ideal coaxial wire model to a more realistic model where the wire presents dielectric losses. Those losses result from the response time of the dielectric used in the wire and of its conductivity $\sigma$.

With the Maxwell’s equations the propagation constant $\gamma$ for a TEM wave and the associated frequency response $H(\omega)$ can be obtained:¹⁵

$$ \gamma^2 = -\omega^2\varepsilon\mu_0 + j(\mu_0\sigma\omega) \quad (14) $$

$$ H(\omega) = \exp(-\gamma L) $$

where $\mu_0$ is the permeability of vacuum.
<figure><img src="image_105.png" /><figcaption>FIG. 2. Plane wave propagation in FR4 substrate considering the Debye model compared to the constant loss angle approximation, Δθ is the phase difference between models</figcaption></figure>

With for instance the Debye model,¹⁶ along with the numerical values from,¹⁷ one obtains (dashed black) the dispersion plot of Figure 2(a) and the attenuation plot of Figure 2(b). Let us introduce at this point the classical approximation where the losses are sufficiently small so that the transfer function can be approximated as follows: $H(\omega) = G(\omega) \times e^{-j\frac{\omega L}{v_{\Phi}}} = e^{-\frac{1}{2}\frac{\omega}{c}\sqrt{\epsilon_r}\tan(\delta)L} \times e^{-j\frac{\omega L}{v_{\Phi}}}$, In this Equation c, is the speed of light in vacuum, $\epsilon_r$ the relative dielectric constant, L the length of the wire and $\delta$ the loss angle considered as constant. This yields:

$$ \frac{1}{\Omega} = -\frac{\omega}{2c} \sqrt{\epsilon_r \tan(\delta)L} \\ \ln(H(\omega)) = -j \frac{L}{v_{\Phi}} \omega - \frac{1}{\Omega} \qquad (15) $$

Let us note that this approximation is not physical. Indeed, considering the pair of Equations (8), the phase velocity cannot remain constant when the amplitude varies. The phase shift and the amplitude corresponding to this approximation has also been plotted in Figure 2 (green solid). One can verify that the approximation is very good insofar $\theta(\omega)$ is concerned (both curves seem superimposed), and correct for $\log(G(\omega))$ up to around 15 GHz.

At this point, it is interesting to rewrite Equation (14) splitting the imaginary part $\epsilon_i(\omega)$ and real part $\epsilon_r(\omega)$ of $\epsilon(\omega)$: $\gamma^2 = -\omega^2\epsilon_r\mu_0 + j(\mu_0\sigma\omega - \omega^2\epsilon_i\mu_0)$. Without losses only the real part of this Equation remains which corresponds to the phase variation $\frac{L}{v_{\Phi}}\omega$ of the small and constant $\delta$ approximation: $\frac{1}{v_{\Phi}} = \sqrt{\lim_{\omega\to\infty} \epsilon_r(\omega)\mu_0}$. Finally by considering the difference between the constant loss angle model and the Debye model one obtain the blue $\Delta\theta$ curve plotted in Figure 2. We have verified that when applying the K-K relations to $\Delta\theta$, the attenuation corresponding to the Debye model (Figure 2(b) Black dashed) is obtained. This is in fact very general and is demonstrated in the next section.

## 2. Making K-K relations work with the small and constant loss angle approximation.

Many times the measurement bandwidth is restricted enough so that the constant loss angle approximation applies, and one obtains the aspect of the plots presented in Figure (3) for the dispersion and the attenuation. In this figure dispersion curves with and without some so called “anomalous dispersion” are presented.

With that approximation, it is realistic to consider media where the two following limits exist:

$$ \left\{ \begin{array}{l} \lim_{\omega \to \pm\infty} -\frac{\theta(\omega)}{\omega} = \Theta \\ \\ \lim_{\omega \to \pm\infty} \frac{\ln(G(\omega))}{|\omega|} = -\frac{1}{\Omega} \end{array} \right. \qquad (16) $$
<figure><img src="image_106.png" /><figcaption>FIG. 3. Dispersion and attenuation with the low and constant loss angle approximation.</figcaption></figure>

In plain words, it means that both $G(\omega)$ and $\theta(\omega)$ are asymptotic to straight lines whose slopes are respectively $-\frac{1}{\Omega}$ and $\Theta$. Using those two limits, one obtains deviations to standard propagation that are limited in frequency:

$$ \left\{ \begin{array}{l} \lim_{\omega \to \pm\infty} (-\theta(\omega) - \omega\Theta) = 0 \\ \lim_{\omega \to \pm\infty} \left( \ln(G(\omega)) + \frac{|\omega|}{\Omega} \right) = 0 \end{array} \right. \qquad (17) $$

Using the pair of conditions (17), via the alternate classical demonstration of the K-K relations$^{8,18}$ which consists in integrating in the complex plane using Cauchy integral theorem, one can obtain the following new relations:

$$ \ln(G(\omega) + \frac{|\omega|}{\Omega}) = -\mathcal{H}_T(-\theta(\omega) - \omega\Theta) \qquad (18) $$

$$ -\theta(\omega) - \omega\Theta = \mathcal{H}_T \left( \ln(G(\omega)) + \frac{|\omega|}{\Omega} \right) \qquad (19) $$

The absolute value in Equations (16) and (19) turns $\ln(G(\omega))$ into an even function when symmetrizing over negative frequencies.

In order to use these two Equations one must have some knowledge of the system studied. The limits $\Theta$ or $\frac{1}{\Omega}$ must be known. Note that $\Theta = \lim_{\omega \to +\infty}(-\frac{\theta(\omega)}{\omega})$, is sometime called the front velocity,$^{19}$

In what follows these results are applied to various cases chosen to cover some of the most typical situations that can be encountered when working at microwave frequencies.

### C. Photonic Crystal at microwave frequencies

A photonic crystal as the one presented in the Figure 4 has been realized. This structure is made of microstrip lines whose characteristic impedance is alternately $Z_0 = 50 \Omega$ and $Z_1 = 100 \Omega$. The substrate used is standard FR4 substrate 0.8 mm thick. The length of each segment has been chosen to be equal to a quarter of the wave length for each characteristic impedance at 2 GHz. This way at 2 GHz, the wave reflected at the end of each $Z_0$ line segment interferes destructively with the wave entering the segment thus creating the photonic effect. This is often referred to as the Bragg condition. It can be written for this line $\omega = \omega_0 + 2n\omega_0$. Where $\omega_0$ is the frequency for which the length of each segment is $\frac{2}{n}$ (2 GHz) and $n$ is a positive integer. At these frequencies, the amplitude of the waves is greatly diminished which leads to frequency bands where the transmission coefficient is very low. Those structures have known some renown recently because the group velocity is superluminal in these bands.

In order to demonstrate the application and efficiency of the method we have performed the measurement of $\ln(G) = \ln(|S_{21}|)$ using spectrum analyzer (SPA) and a generator. The same measurement has also been performed using a vector analyzer (VNA) along with the measurement
<figure><img src="image_107.png" /><figcaption>FIG. 4. Photonic crystal built up from microstrip lines.</figcaption></figure>

<figure><img src="image_108.png" /><figcaption>FIG. 5. $|S_{21}|$ measured using a Vectorial Network Analyzer (VNA) and a Spectrum Analyzer (SPA).</figcaption></figure>

of $\theta$. The results are presented Figure (5). The results match up to 12 GHz. Beyond this frequency a difference is visible. This difference can be accounted for the calibration procedure (simple two ports SOLT) used for the VNA which is certainly not the most adapted for such a large bandwidth. The measurement from the SPA were compensated using a THRU measurement from the generator.

One can verify in Figure 5 that $\log(|S_{21}|)$ behaves asymptotically. By mean of asymptotic analysis a “rectified” gain as defined by the left member of Equation (18) can be obtained: $G(\omega) + \frac{|\omega|}{\Omega}$. This gain has been filled in with zeros for frequencies below the first effectively measured frequency by the analyzer in order to extrapolate the signal down to the null frequency. It has also been symmetrized over the negative frequencies to take into account the considerations presented in section II C. Finally the argument of $S_{21}$ has been calculated using Equation (19). Figure 6 presents both the calculated and the measured argument of $S_{21}$. Of course to compare both results one needs to know the asymptotic behavior of $\theta$ as well. To obtain a figure with an interesting magnification factor the arguments are presented without the asymptotic behavior.

The results are identical as long as the amplitudes measured by the VNA and the SPA are identical which validates the method.

### D. Left/right handed hybrid media presenting evanescent waves

A left handed media is a media where the phase velocity $V_\varphi = \frac{\omega}{k}$ is negative over a frequency band.⁵ A simple way to build up a left handed line is presented in Figure 7 between Port 2 and Port 3. The lumped capacitors and inductors yield the left handed properties of the line at low frequencies. The transmission lines segments used to connect the lumped components turn the line into a standard right handed line by adding their distributed capacitance C and inductance L when the frequency gets higher. At even higher frequency, when the physical length of the line gets small compared the wavelength, propagation has to be considered. A thorough study of that kind of line can be found in.²⁰ This kind of line can be used to realize couplers²¹ with special properties by taking advantage of the left handed behavior. These couplers are made as presented in Figure 7 of a left handed lines coupled to a standard microstrip line. In this figure for instance the coupled ports are ports 1 and 2 whereas port 4 is the direct port.
<figure><img src="image_109.png" /><figcaption>FIG. 6. Phase measured using a VNA compared to the phase calculated from a spectrum analyzer measurement using K-K.</figcaption></figure>

<figure><img src="image_110.png" /><figcaption>FIG. 7. Hybrid left handed line coupled to a standard microstrip line.</figcaption></figure>

At low frequencies, when the effect of line segments and coupling are negligible, the capacitors and inductors induce a phase shift that can be calculated using the Kirchhoff's Laws:

$$ \cos(\theta) = 1 - \frac{1}{LC\omega^2} \qquad (20) $$

where $\theta$ is the phase shift per cell. At frequencies below $\frac{1}{2\pi\sqrt{2LC}}$, the waves in the structure are evanescent and no phase shift occurs. Consequently the phase shift remains fixed at $-\pi$. Furthermore being in a forbidden band implies a very high attenuation at consequently a very low signal to noise ration. In Figure 8(a), the measurement of $|S_{32}|$ has been plotted. One can see indeed that below about 1 Ghz only noise is measured.

When the frequency gets higher, both the distributed capacitance and inductance of the line and the coupling effect, introduce an attenuated frequency band (centered about 1.8 Ghz) which has consequences on the dispersion. At even higher frequencies the propagation occurring in the line segments becomes preponderant, and $|S_{32}|$ remains constant. This line is interesting regarding the application of K-K relations. Indeed, the waves are successively evanescent, left-handed and right handed. Further more the measurement starts with noise which is an additional interesting issue.

In Figure 8(b) the solid blue curve represents the measured argument of $S_{32}$. As the signal is lost in noise at low frequency, the argument of $S_{32}$ is not significant below 1 GHz. When the signal is great enough to be measured above 1 Ghz as the initial phase value has been lost, the argument is biased and only the phase variations are correct.

In order to be able to apply K-K relations with success and to obtain the black dashed plot of Figure 8(b) it is first necessary to deal with the data below 1 GHz which are not significant. We use here the theoretical value of $|S_{32}|$ obtained from Equation (20) which has been plotted (black dashed) in Figure 8(a) to extrapolate the data below 1 Ghz. The black curve of Figure 8(b) is then obtained by successively computing the Hilbert transform and adding the phase due to the propagation knowing the phase velocity in the line segments. The result obtained is very good and has allowed us to
<figure><img src="image_111.png" /><figcaption>FIG. 8. Application of the K-K relations to measurements obtained on an hybrid left/right handed line (per cell data).</figcaption></figure>

verify the left handed behavior of this line between 1 and 3 Ghz when the argument is negative and increasing with the frequency.

## IV. CONCLUSION

Though Kramers-Kronig relations are an interesting tool, it simply relates the real part to the imaginary part of the response of a physical system. Practically, it is often required to be able to extrapolate the phase shift from the attenuation and vice-versa. In a general way, it is not possible to apply the Kramers-Kronig relations to perform such calculation. Nevertheless, by means of some assumptions on the modulus of the frequency response, it is possible to apply the Kramers-Kronig relation on many physical systems. Furthermore we have shown that the case of plain propagation has to be dealt with in a particular way.

In this work, we have firstly recalled how the Kramers-Kronig relations can be obtained from causality considerations. We have proposed and thoroughly tested a simple and practical method to compute the Hilbert transform that subtend the Kramers-Kronig relations. This method takes into account algorithm discrepancy and the corresponding required manipulation of the data.

We have applied successfully the method to various special cases chosen because of their representatives behaviors. We have verified that method can be applied whether there is propagation or not. In the case of propagation the waves can be evanescent, right or left handed, or be the steady state resulting from interferences as in the case of photonic crystals.

The precision reached for data extrapolation is very good. It greatly depends nevertheless on the knowledge of the asymptotical behavior of the dispersion curve or of the attenuation curve of the media at high frequencies.

## APPENDIX

### 1. Converting from K-K integral form over positive frequencies to its Hilbert transform form

$$ \theta(\omega) = -\frac{2\omega}{\pi} \int_{0}^{+\infty} \frac{\ln(G(u))}{\omega^2 - u^2} du \\ = -\frac{1}{\pi} \int_{0}^{+\infty} \left[ \frac{(\omega - u)\ln(G(u))}{\omega^2 - u^2} + \frac{(\omega + u)\ln(G(u))}{\omega^2 - u^2} \right] du \\ = -\frac{1}{\pi} \left[ \int_{0}^{+\infty} \frac{\ln(G(u))}{(\omega + u)} du + \int_{0}^{+\infty} \frac{\ln(G(u))}{(\omega - u)} du \right] $$
Changing the integration variable from $u$ to $-u$ in the left integral with $G(-u) = G(u)$:

$$ \theta(\omega) = -\frac{1}{\pi} \left[ \int_{-\infty}^{0} \frac{\ln(G(u))}{(\omega - u)} du + \int_{0}^{+\infty} \frac{\ln(G(u))}{(\omega - u)} du \right] $$

$$ \theta(\omega) = \frac{1}{\pi} \int_{-\infty}^{\infty} \frac{\ln(G(u))}{(u - \omega)} du $$

$$ \theta(\omega) = \mathcal{H}(\ln(G(\omega))) $$

<footnote n="1">¹Brillouin, *Wave propagation and group velocity*, Pure and applied physics, third edition, (Academic press inc., 111 Fifth avenue New York, 1960).</footnote>

<footnote n="2">²A. Haché and L. Poirier, “Anomalous dispersion and superluminal group velocity in a coaxial photonic crystal: Theory and experiment,” *Physical Review E* **65**, 036608 (2002).</footnote>

<footnote n="3">³CHEN Kai, WU Ling-An, and SHIH Yan-Hua, “Is ‘Superluminal’ Light Propagation Possible in Dispersive Media?,” *Chin Phys Lett* **21**(5), 770 (2004).</footnote>

<footnote n="4">⁴S. P. Ojha, K. B. Thapa, and S. K. Singh, “Superluminal propagation in plasma photonic band gap materials,” *Optics* **119**, 81–85 (2008).</footnote>

<footnote n="5">⁵V. G. Veselago, “The Electrodynamics of substances with simultaneous negative values of $\epsilon$ and $\mu$,” *Soviet Physics Uspekhi* **10**(4), 509–514 (1968).</footnote>

<footnote n="6">⁶J. B. Pendry, A. J. Holden, D. J. Robbins, and W. J. Stewart, “Magnetism from conductors and enhanced nonlinear phenomena,” *IEEE Transactions on Microwave Theory and Techniques* **47**(11), 2075–2084 (1999).</footnote>

<footnote n="7">⁷C. F. Bohren, “What did Kramers and Kronig do and how did they do it?,” *Eur. J. Phys* **31**, 573–577 (2010). doi: 10.1088/0143-0807/31/3/014.</footnote>

<footnote n="8">⁸Frederic Wooten, *Optical properties of solids*, (Academic Press Inc, 111 fifth avenue, New York, 1972).</footnote>

<footnote n="9">⁹J. S. Toll, “Causality and the dispersion relation: Logical foundations,” *Physical review* **104**(6), 1760–1770 (1956).</footnote>

<footnote n="10">¹⁰R. Kendall, J. Mobbey, and J. G. Miller, “Causality-Imposed (Kramers-Kronig) Relationships Between Attenuation and Dispersion,” *IEEE trans on Ultrasonic, Ferroelectrics, and Frequency Control* **52**(2), 822–833 (2005).</footnote>

<footnote n="11">¹¹S. Temkin, “Attenuation and dispersion of sound in bubbly fluids via the Kramers-Kronig relations,” *J Fuid Mech* **211**, 61–72 (1990).</footnote>

<footnote n="12">¹²H. Tuonomen, E. Gomov, J. A. Zeitler, J. Aaltonen, and K.-E. Peiponen, “Using modified Kramers-Kronig relations to test transmission spectra of porous media in THzTDS,” *Optics Letters* **35**(5), 631–633 (2010).</footnote>

<footnote n="13">¹³J. Lucas, E. Géron, T. Ditchi, and S. Holé, Practical use of the Kramers-Kronig relation at microwave frequencies, Application to photonic like lines and left handed materials, *PIERS*, (2011).</footnote>

<footnote n="14">¹⁴Ronald N. Bracewell, *The Fourier transform and its applications*, third edition, (Mac Graw-Hill international editions, 111 Fifth avenue New York, 2000).</footnote>

<footnote n="15">¹⁵Jackson J. D., *Classical Electrodynamics*, sixth edition, 1967.</footnote>

<footnote n="16">¹⁶Ramo S., Whinnery J. R., and Van Duzer T., *Fields and waves in communication electronics*, (John Wiley and Sons, 1984).</footnote>

<footnote n="17">¹⁷M. Y. Kolendisteva, K. N. Rosanov, A. Orlandi, and J. L. Droniak, “Extraction of Lorentzian and Debye of dielectric and magnetic dispersive materials for FDTD modeling,” *Journal of Electrical Engineering* **53**(9/S), 97–100 (2002).</footnote>

<footnote n="18">¹⁸S. J. Orphanidis, *Electromagnetic Waves and Antennas*, ECE Department Rutgers University, http://www.ece.rutgers.edu/orfanidi/ewa/, in line edition, 2008.</footnote>

<footnote n="19">¹⁹Papoulis A., *The Fourier integral and its applications*, 1962.</footnote>

<footnote n="20">²⁰Christophe Caloz and Tatsuo Itoh, *Electromagnetic Metamaterials: Transmission Line Theory and microwave applications/The Engineering Approach*, (Wiley Interscience, 2006).</footnote>

<footnote n="21">²¹E. Géron, T. Ditchi, J. Lucas, and S. Holé, Asymmetric Microstrip Right/Left-Handed line coupler with variable coupling ratio, *PIERS*, (2011).</footnote>