function [aph] = get_ciotti_aph(wl,Sf)

%
% Jeremy Werdell, NASA Goddard Space Flight Center, July 2013
% from Ciotti and Bricaud 2006, Limnology and Oceanography Methods
%

if exist('Sf') == 0; Sf = 0.5; end
	
pico = [1.7439,1.8264,1.9128,1.9992,2.0895,2.1799,2.2702,2.3684,2.4666,2.5687,2.6669,2.7612,2.8437,2.9183,2.9890,3.0479,3.1029,3.1500,3.1854,3.2089,3.2247,3.2325,3.2286,3.2168,3.1932,3.1540,3.1029,3.0361,2.9576,2.8712,2.7848,2.6944,2.5137,2.4273,2.3488,2.2781,2.2486,2.2192,2.1720,2.1328,2.1013,2.0660,2.0267,1.9835,1.9285,1.8657,1.7989,1.7203,1.6339,1.5357,1.4336,1.3276,1.2176,1.1076,1.0016,0.8994,0.8013,0.7109,0.6284,0.5538,0.4870,0.4320,0.3782,0.3307,0.2875,0.2486,0.2137,0.1842,0.1599,0.1402,0.1233,0.1080,0.0935,0.0789,0.0656,0.0530,0.0424,0.0344,0.0290,0.0260,0.0258,0.0268,0.0304,0.0320,0.0331,0.0347,0.0355,0.0363,0.0382,0.0401,0.0416,0.0428,0.0432,0.0432,0.0432,0.0424,0.0416,0.0408,0.0408,0.0424,0.0452,0.0503,0.0562,0.0628,0.0695,0.0758,0.0821,0.0880,0.0939,0.1002,0.1060,0.1123,0.1178,0.1229,0.1261,0.1280,0.1288,0.1296,0.1308,0.1331,0.1371,0.1422,0.1493,0.1591,0.1728,0.1909,0.2137,0.2416,0.2757,0.3178,0.3692,0.4281,0.5499,0.6009,0.6324,0.6402,0.6324,0.6245,0.5892,0.5342,0.4674,0.3967,0.3276,0.2635,0.2078,0.1618,0.1249,0.0958,0.0746,0.0601,0.0503];

% 2002 values pico = [1.682,1.734,1.8,1.89,1.978,2.057,2.162,2.269,2.327,2.398,2.457,2.533,2.614,2.663,2.749,2.804,2.84,2.915,2.947,2.978,3.014,3.032,3.011,2.965,2.937,2.888,2.816,2.783,2.706,2.655,2.59,2.526,2.455,2.402,2.331,2.281,2.205,2.136,2.063,2.049,1.998,1.93,1.918,1.897,1.867,1.812,1.776,1.701,1.648,1.522,1.439,1.373,1.27,1.162,1.04,0.961,0.886,0.794,0.734,0.665,0.617,0.544,0.522,0.486,0.448,0.391,0.375,0.336,0.305,0.292,0.288,0.261,0.245,0.214,0.194,0.187,0.138,0.137,0.111,0.094,0.095,0.07,0.053,0.076,0.064,0.043,0.05,0.051,0.065,0.067,0.084,0.111,0.072,0.073,0.073,0.099,0.07,0.095,0.085,0.09,0.086,0.068,0.078,0.069,0.09,0.096,0.094,0.084,0.105,0.128,0.119,0.126,0.138,0.146,0.135,0.175,0.189,0.176,0.203,0.19,0.19,0.191,0.174,0.197,0.176,0.168,0.16,0.217,0.244,0.286,0.381,0.437,0.52,0.66,0.716,0.824,0.846,0.816,0.891,0.869,0.812,0.741,0.605,0.496,0.372,0.278,0.215,0.113,0.075,0.047,0.009]

micro = [1.574,1.584,1.6,1.617,1.633,1.654,1.669,1.674,1.684,1.697,1.708,1.71,1.716,1.737,1.763,1.793,1.812,1.827,1.83,1.834,1.824,1.8,1.771,1.741,1.712,1.685,1.667,1.65,1.641,1.631,1.631,1.623,1.616,1.606,1.592,1.568,1.542,1.509,1.481,1.459,1.437,1.415,1.399,1.387,1.377,1.367,1.349,1.338,1.319,1.301,1.271,1.242,1.222,1.196,1.169,1.141,1.118,1.096,1.075,1.057,1.035,1.013,0.992,0.977,0.959,0.944,0.927,0.909,0.888,0.868,0.847,0.826,0.806,0.785,0.764,0.737,0.711,0.682,0.653,0.626,0.604,0.58,0.555,0.535,0.514,0.501,0.487,0.478,0.475,0.468,0.464,0.459,0.452,0.452,0.449,0.443,0.433,0.424,0.416,0.406,0.401,0.4,0.403,0.408,0.416,0.429,0.443,0.458,0.473,0.487,0.495,0.499,0.504,0.514,0.521,0.525,0.532,0.535,0.534,0.535,0.532,0.528,0.526,0.528,0.538,0.549,0.574,0.605,0.655,0.72,0.798,0.889,0.979,1.068,1.147,1.207,1.243,1.249,1.227,1.174,1.096,1.004,0.893,0.767,0.635,0.516,0.409,0.323,0.253,0.2,0.158]; 
% compute the basis vectors
pico = pico * 0.023 / 0.891;
micro = micro * 0.0086 / 1.249; 

w0 = 400:2:700;

a0 = (Sf * pico) + ((1. - Sf) * micro);

aph = interp1(w0,a0,wl,'cubic');

