function cost = GSM01_cost(XX);

global rrs_obs aw bbw;
% This version of the GSM01 model uses constants from Maritorena et al.  (2002)

chl = XX(1); ag = XX(2); bbp = XX(3);

g1 = 0.0949; g2 = 0.0794; S = 0.0206; S_bbp = 1.0337;

% wave = [412 443 490 510 555]'; % seawifs
wave = [412 443 488 531 551]'; % modis
%bbw = [0.00333 0.00239 0.00155 0.00130 0.000925]';	%%%%% Smith/Baker1981
%aw =  [0.00469 0.00721 0.01500 0.0325  0.05920]';    %%% Pope/Fry1997
aph = [0.00665 0.05582 0.02055 0.01910 0.01015]';   %%% Maritorena et al (2002)


atotal = aw + aph.*chl + ag.*exp(-S*(wave-443));
bbtotal = bbw + bbp .* ((443*ones(size(wave)) ./ wave).^S_bbp);

x = bbtotal ./ (atotal + bbtotal);
    
rrs = g1 .*x + g2 .* x .* x;

cost = sum((rrs-rrs_obs).^2);
