Hi @gvh,

The two orthogonalizations happen in

```
Xext = repmat(tXh,1,nSig) ;
Yext = repelem(tXh,1,nSig) ;
Yext_p = imag(bsxfun(@times, Yext, conj(Xext)./abs(Xext)));
```

This operation is performed for each band, we can let aside the frequency dimension for now. Consider two (column vector) signals `u`

and `v`

with their corresponding analytic signals [`ũ`

, `ṽ`

], and the variable `tXh`

containing those signals, this is to say: `tXh`

= [`ũ`

, `ṽ`

]

Then, `Xext`

is equal to [`ũ`

, `ṽ`

, `ũ`

, `ṽ`

] , and `Yext`

is equal to [`ũ`

, `ũ`

, `ṽ`

, `ṽ`

] .

The operation is performed for each column vector, thus: `Yext_p`

= [**ũ⟂ũ**, **ũ⟂ṽ**, **ṽ⟂ũ**, **ṽ⟂ṽ**]. Lastly, the correlation is computed and reshaped in

```
CorrMat = reshape(diag(corrFcn(abs(Xext),abs(Yext_p))),nSig,nSig) ;
```

In this case, `CorrMat`

is a [2×2] matrix with:

`[corr(|ũ|, |ũ⟂ũ|), corr(|ṽ|, |ũ⟂ṽ|) ; corr(|ũ|, |ṽ⟂u|), corr(|ṽ|, |ṽ⟂ṽ|)]`

The average happens later in the line:

```
A(:,:,t,f) = (abs(CorrMat) + abs(CorrMat'))/2 ;
```

Best,

Raymundo