MCSpreadOption(CallPutFlag, S1, S2, X, T, r, b1, b2, v1, v2, rho , nSimulations)

function MCSpreadOption(CallPutFlag, S1, S2, X, T, r, b1, b2, v1, v2, rho , nSimulations)
%
%   Monte Carlo Simulation Spread Options
%       Key words: finance, option valuation, Monte Carlo simulation
%
%   Author:       Espen Gaarder Haug
%                                               (The author of "The Complet Guide to Option Pricing Formuals," McGraw-Hill 1997
%   E-mail:       [e-mail]
%   hompage:            http://home.sol.no/~eshaug/
%   Environment:  MATLAB 5.3 in Windows 98  
%
%       Application: Monte Carlo Simulation of Spread Option Values:
%
%       The MCSpreadOption function can be used for valuation of spread options
%       Spread options are actively traded on New York Mercantile Exchange (NYMEX) 
%       on the difference/spread between different oil qualities. Spread options
%       are also popular in many other markets.
%
%       The code can easely be modified to also value other option contracts dependent on
%       two correlated assets. Assuming correlated geometric Brownian motions
%
%       Variable Description:
%
%                               CallPutFlag= 'c' gives call option, 'p' gives put option
%                               S1= Asset one
%                               S2= Asset two
%                               X= Strike price
%                               T= time to option maturity in number of years
%                               r= Risk-free-rate
%                               b1= Cost of carry asset one
%                               b2= Cost of carry asset two
%                               v1= Volatility asset one
%                               v2= Volatility asset two
%                               rho= Correlation between assets (asset return)
%                               nSimulations=Number of simulations used for Monte Carlo simulation
%                                            More simulations increase accurancy, typically minimum 10000
%
%                               At maturity the payoff from a call spread option is max(S1-S2-X,0)
%                               At maturity the payoff from a put spread option is max(X-S1+S2,0)
%
%
%           For example :
% 
%           MCSpreadOption('c',122,120,3,0.5,0.1,0,0,0.25,0.2,0.5,10000)
%                               
%
% To learn more about pricing of Spread Options
% and other exotic options see:
% Haug, E. G. 1997: "The Complete Guide To Option Pricing Formulas," McGraw-Hill


        if CallPutFlag == 'c' 
      z = 1; %  Will return call option value
   else
      z = -1;   % Will return put option value
   end
    
    Drift1 = (b1 - v1 ^ 2 / 2) * T; %Drift asset one
    Drift2 = (b2 - v2 ^ 2 / 2) * T; % Drift asset two
    v1Sqrdt = v1 * sqrt(T);
    v2Sqrdt = v2 * sqrt(T);
         sum=0;
    for i = 1: nSimulations,
        Epsilon1 = normrnd(0,1); % Random number from a standard normal distribution
        Epsilon2 = rho * Epsilon1 + normrnd(0,1) * sqrt(1 - rho ^ 2); %
            St1 = S1 * exp(Drift1 + v1Sqrdt * Epsilon1);
            St2 = S2 * exp(Drift2 + v2Sqrdt * Epsilon2);        
            sum = sum + max(z * (St1 - St2 - X), 0); 
     end

    exp(-r * T) * (sum / nSimulations)% Return Monte Carlo estimate of spread option