Fixing a Half-Value Result in Fourier Call Pricing

# Heston model with underlying BS dynamics always gives 1/2 of the right value, what am I doing wrong?


# Heston model with underlying BS dynamics always gives 1/2 of the right value, what am I doing wrong?












Just as an exercise I'm trying to follow this paper: https://arxiv.org/ftp/arxiv/papers/1502/1502.02963.pdf

In the section 2.2 it calculates the value of a Call using the characteristic function of the underlying under Black-Scholes dynamics.

I'm using Python, and strangely my result seems to be exactly 1/2 of the value of the call option, because both the Pi1 and Pi2 I get are also 1/2 of what is supposed to be.

For instance I'm calculating the value of a call option with parameters:

- S = 100

- K = 10

- t = 1

- vol = 0.2

- r = 0

which should give me a Pi1 (delta) equal to 1 (since it's deep ITM) but instead I get a Pi1 = 0.49999999 and a call value of 44.99999

Could you help to to point out where might I have an error? I really hope is not a silly mistake, I've been trying to figure it out for the last 2-3 hours but I cannot see anything wrong. Here is the code I'm using:

```
import numpy as np
from scipy import stats, special, integrate

def char_func(x,s,vol,t=1,r=0):
    mean = np.log(s) + (r - 0.5*vol*vol) * t
    var = vol*vol*t
    w = np.exp(1j*x*mean - x*x*vol*vol*0.5)
    return w.real

def call_value(s,k,vol,t=1,r=0):
    def integrand(x,s,k,vol,t=1,r=0):
        I = np.exp(-1j*x*np.log(k)) * char_func(x-1j,s,vol,t,r)/(1j*x*char_func(-1j,s,vol,t,r))
        return I.real

    def integrand2(x,s,k,vol,t=1,r=0):
        I = np.exp(-1j*x*np.log(k)) * char_func(x,s,vol,t,r)/(1j*x)
        return I.real

    int1 = integrate.quad(integrand,0,np.inf,args=(s,k,vol,t,r))
    int1 = 0.5 + int1[0]/np.pi

    int2 = integrate.quad(integrand2,0,np.inf,args=(s,k,vol,t,r)))
    int2 = 0.5 + int2[0]/np.pi

    return s*int1 - np.exp(-r*t)*k*int2

print(call_value(100,10,0.2))
```

Thanks

## Answer by Kevin (score 2, accepted)

https://quant.stackexchange.com/a/61453

You're implementing the formula $$C=S_0e^{-qT}\Pi_1-Ke^{-rT}\Pi_2$$ where \begin{align*} \Pi_1&=\frac{1}{2}+\frac{1}{\pi}\int_0^\infty \text{Re}\left(\frac{e^{-iu\ln(K)}\varphi_{\ln(S_T)}(u-i)}{iu\varphi_{\ln(S_T)}(-i)}\right)\text{d}u,\\ \Pi_2&=\frac{1}{2}+\frac{1}{\pi}\int_0^\infty \text{Re}\left(\frac{e^{-iu\ln(K)}\varphi_{\ln(S_T)}(u)}{iu}\right)\text{d}u. \end{align*}

A few words about your code. You should not output the real part of the characteristic function, but the complex value. That is, you need `return w` and not `return w.real`. Also, you define `var` but you don't use it.

I get 89.99999999990168 as option value if I run the below code. This makes sense, the option is so deep ITM, its value is essentially $S_0e^{-qT}-Ke^{-rT}\approx90$

```
import numpy as np
from scipy import stats, special, integrate

def char_func(x,s,vol,t=1,r=0):
    mean = np.log(s) + (r - 0.5*vol*vol) * t
    var = vol*vol*t
    w = np.exp(1j*x*mean - x*x*var*0.5)
    return w

def call_value(s,k,vol,t=1,r=0):
    def integrand(x,s,k,vol,t=1,r=0):
        I = np.exp(-1j*x*np.log(k)) * char_func(x-1j,s,vol,t,r)/(1j*x*char_func(-1j,s,vol,t,r))
        return I.real

    def integrand2(x,s,k,vol,t=1,r=0):
        I = np.exp(-1j*x*np.log(k)) * char_func(x,s,vol,t,r)/(1j*x)
        return I.real

    int1 = integrate.quad(integrand,0,np.inf,args=(s,k,vol,t,r))
    int1 = 0.5 + int1[0]/np.pi

    int2 = integrate.quad(integrand2,0,np.inf,args=(s,k,vol,t,r))
    int2 = 0.5 + int2[0]/np.pi
    
    return s*int1 - np.exp(-r*t)*k*int2

print(call_value(100,10,0.2))
```

Note: The option pricing formua you're using is very simple. It was one of the original Fourier methods discovered. Try instead something like \begin{align} c(S_0;K,T) &= c^\mathrm{BS}(S_0;K,T) +\frac{\sqrt{K}e^{-rT}}{\pi}\int_0^\infty \frac{\text{Re}\left( e^{-iku} \left(\varphi_T^\mathrm{BS}(u-0.5i)-\varphi_T(u-0.5i) \right)\right)}{{u^2+0.25}}\mathrm{d}u, \end{align} where $c^\mathrm{BS}$ is the Black Scholes option price and $\varphi_T^\text{BS}$ the corresponding log-normal characteristic function. That formula uses a control variate and integrates along a different contour.