function madCyclePlot(inputData, startDay, endDay, dtn, errorBars, perBarSwitch,...
    outputDir, expID, yLabelStr, legendStr, Col, titleStr, forcedYTick, SW, inputDataScale, d, U, vWinds, nFreq, maxFreq, windYLabel)


      

        
% time units (in seconds)
minuteInSeconds = 60;                   % one minute in seconds
hourInSeconds   = 60*minuteInSeconds;   % one hour in seconds
dayInSeconds    = 24*hourInSeconds;     % one day in seconds
yearInSeconds   = 365*dayInSeconds;     % one year in seconds


if nargin > 11 && isempty(titleStr)>0
    titleSwitch = 0;
elseif nargin > 11 && isempty(titleStr)==0
    titleSwitch = 1;
else
    titleSwitch = 0;
end


if nargin > 12 && isempty(forcedYTick)==0
    niceYTickSwitch = 0;
elseif nargin > 12 && isempty(forcedYTick)>0
    niceYTickSwitch = 1;
else
    niceYTickSwitch = 1;
end

if nargin > 13 && isempty(SW)==0
    shortWaveSwitch = 1;
elseif nargin > 13 && isempty(SW)==1
    shortWaveSwitch = 0;
else
    shortWaveSwitch = 0;
end

if nargin > 14;
    fluxScaleSwitch = 1;
    fluxScaleData = {inputDataScale};
else
    fluxScaleSwitch = 0;
    fluxScaleData = {};
end
    
if nargin > 15;
    windRoseSwitch = 1;
    windData = {d, U};
    allData = [inputData;fluxScaleData;windData'];
else
    windRoseSwitch = 0;
    windData = {};
    allData = inputData;

end





%% 
nm  = 12;                 % number of months per year
nw = length(allData);   % number of datasets
startNum = datenum(startDay);
endNum   = datenum(endDay);
dayNum   = endNum - startNum + 1;
lenData  = cellfun(@length, allData);
lenData  = lenData(1);

dto = round(dayNum * dayInSeconds./lenData); % old time step in seconds
dtd = dayInSeconds/dtn;   % number of time steps per day
dtr = dtn./dto;            % fraction new/old timestep
  


if sum(mod(dtr,1)) > 0
    error('The new time step has to be a multiply of the old one.')
end


% calculate new mean
if sum(dtn-dto)>0
    inputData = cellfun( @(x,y) nanmean(reshape(x,y,[]))', allData,num2cell(dtr), 'UniformOutput',false );
elseif dtn<dto
    error('The new time step has to be greater than the old one (interpolation unsupported).')
end


dtf       = 1/dtd;             % timestep's fraction of a day
startDate = startNum;          % start date 
endDate   = endNum - dtf;  % end date

% create time vector
timeData   = cellfun( @(x) linspace(startDate,endDate,length(x))', allData, 'UniformOutput',false );

% allocate the respective month's number to each time step
monthCount = cellfun( @month, timeData,'UniformOutput',false);


% sort data by month
monthData = cell(nm,nw);
for i=1:nm
    monthData(i,:) = cellfun(@(x,y) x(y==i) , allData, monthCount, 'UniformOutput',false);
end

if windRoseSwitch == 1
    monthWindData = monthData(:,end-1:end);
    monthData = monthData(:,1:end-2);
    nw = nw - 2;
end




% number of days of the respective month
nday = cellfun( @(x) length(nanmean(reshape(x,dtd,[]))'), monthData, 'UniformOutput',false );
nday = cell2mat(nday);

% reshape into diurnal cycles
monthData = cellfun( @(x) reshape(x,dtd,[]), monthData, 'UniformOutput',false);

diuFraction = cell2mat(cellfun( @(x) sum(sum(isnan(x))==0)/size(x,2), monthData, 'UniformOutput',false ));


% average by diurnal cycles
diuCyc = cellfun( @(x) nanmean(x')', monthData, 'UniformOutput',false );
diuStd = cellfun( @(x) nanstd(x')', monthData, 'UniformOutput',false );
diuCyc = cell2mat(diuCyc);
diuStd = cell2mat(diuStd);

diuData = cell(nw,2);
for i=1:nw
    diuData{i,1} = reshape(diuCyc(:,i),dtd,[]);
    diuData{i,2} = reshape(diuStd(:,i),dtd,[]);
end

if fluxScaleSwitch == 1
    scaleData = diuData{end,1};
    %diuData   = diuData{1:end-1,:};
    nw = nw-1;
end


%% Plot Settings


%
monthStr   = {'Jan','Feb','Mar','Apr','May','Jun',...
               'Jul','Aug','Sep','Oct','Nov','Dec'};
diurnalStr = repmat({'00','06','12','18'}',nm,1);
percentStr = {'25%','50%','75%','100%'};

           
% time vector containing 12 randomly days for x axis
start12Days = datenum('01-01-2015');
end12Days   = datenum('01-13-2015');
time12Days  = reshape(linspace(start12Days,end12Days,nm*dtd)',dtd,[]);

%% Plot 

% Figure
if windRoseSwitch == 1
    % Wind Rose settings

    Options.colors       = rainbow(length(vWinds));
    Options.AngleNorth   = 0;
    Options.AngleEast    = 90;
    Options.Labels       = {'','','',''};
    Options.TitleString  = '';
    Options.MaxFrequency = maxFreq;
    Options.nFreq        = nFreq;
    Options.vWinds       = vWinds;
    Options.CenteredIn0  = false;
    Options.nDirections  = 12*2;
    Options.min_radius   = 0;

    F = [0.0 0.1 1.0 0.8];
else
    F = [0.0 0.1 1.0 0.6];
end
figure('Color','white','Units','normalized','Position',F)

if windRoseSwitch == 1
    lm = 0.05;      % left margin
    tm = 0.25;      % top margin
else
    if titleSwitch == 1
        lm = 0.03;      % left margin
        tm = 0.1;      % top margin
    else
        lm = 0.03;      % left margin
        tm = 0.06;      % top margin
    end
end
rm = 0.01;      % right margin
bm = 0.12;      % bottom margin
xw = 1-lm-rm;   % x-width
yh = 1-bm-tm;   % y-height
mw  = xw/nm;    % width  of month
P = [lm bm xw yh];
fontSize = 24;

%% Wind Rose Plot

if windRoseSwitch == 1;
    % wind rose scale
    wrs = 2;
    for i=1:nm
        box on
        xP = lm+(i-1)*mw;
        P1 = [xP bm+yh+0.01 mw 0.9*wrs*mw];
        subplot('Position',P1)
        Options.axes      = gca;
        [figure_handle,count,speeds,directions,Table] = WindRose(monthWindData{i,1},monthWindData{i,2},Options);
        legend(gca,'off');
        if i==1
            barPos = get(gca,'Position');
            numTicks = 6;
            c = colorbar('Position', [barPos(1)-0.015 barPos(2)+0.012 0.006  barPos(4)-0.02],...
                'FontSize',fontSize/2,'TickLabels',mat2cell(vWinds(1:2:end)',numTicks),...
                'Ticks',0:1/(numTicks-1):1);
                
            colormap(rainbow)
        end
    end
    annotation('line',[lm lm+xw], [bm+yh+wrs*mw bm+yh+wrs*mw],'LineWidth',1.1);
else
    wrs = 0;
end    
    
%% Monthly Averaged Diurnal Cycle Plot 

subplot('Position',P)
box on; grid on; hold on

if windRoseSwitch == 1;
    ylabel2 = text(-0.035,bm+yh+wrs*mw+0.1,windYLabel,'Rotation',90,...
        'units','normalized','Interpreter','Latex','FontSize',fontSize-10);
end

t = reshape(time12Days,nm*dtd,[]);

maxVal = nan(nw,1);
minVal = nan(nw,1);
for j=nw:-1:1
    if errorBars == 1
        for i=1:nm
        hline{j} = shadedErrorBar(time12Days(:,i),diuData{j,1}(:,i),diuData{j,2}(:,i),...
            {'LineStyle','-','LineWidth',2,'Color',Col{j}},1);
        end

        if fluxScaleSwitch==1 && j==1            
            for i=1:nm*dtd
                xCoord  = time12Days(i);
                yCoord1 = diuData{j,1}(i);
                yCoord2 = scaleData(i);
                
                line([xCoord,xCoord],[yCoord1,yCoord2],'LineWidth',1,'Color',Col{j})
                
            end
        end
        maxVal(j) = max(max(diuData{j,1} + diuData{j,2}));
        minVal(j) = min(min(diuData{j,1} - diuData{j,2}));
    else
        hline{j} = plot(time12Days,diuData{j},'LineWidth',2.5,'Color',Col{j});
        maxVal(j) = max(max(diuData{j,1}));
        minVal(j) = min(min(diuData{j,1}));
    end
end



maxVal = max(maxVal);
minVal = min(minVal);
if niceYTickSwitch == 1
    [yTicks,yMin,yMax,yInt,ynInt] = niceYTicks(minVal,maxVal);
else
    yTicks = forcedYTick;
    yMin   = forcedYTick(1);
    yMax   = forcedYTick(end);
    yInt   = forcedYTick(2)-forcedYTick(1);
    ynInt  = length(forcedYTick)-1;
end



% Axes
dtl = nm/dtd;   % time step length
% set x-axis
set(gca,'xTick',t(1):12/48:t(end), 'xTickLabel',diurnalStr, 'xLim',[t(1) t(end)])
% set y-axis
set(gca,'yTick',yTicks, 'yLim',[yMin yMax])
% set axes size
set(gca,'FontSize',fontSize-7,'FontName','Arial')


% Legend
legLines = nan(1,nw);
for j=1:nw
    if errorBars == 1
        legLines(j) = hline{j}.mainLine;
    else
        legLines(j) = hline{j}(1);
    end
end
% legend strign
leg1 = legend(gca,legLines,legendStr);
% legend design
if windRoseSwitch == 1
    legY   = 0.68;
else
    if perBarSwitch==0
        legY = 0.8;
    elseif numel(legLines)==4
        legY = 0.78;    
    elseif nw > 2 && mean(ismember('\ $G',legendStr{1}))<1
        legY = 0.78;
    elseif nw > 2 && mean(ismember('\ $G',legendStr{1}))==1
        legY = 0.2;
    else
        legY = 0.80;
    end        
end
legPos = [0.905,legY,0,0];
set(leg1,'Color','white','EdgeColor','white',...
    'FontSize',fontSize-3,'Interpreter','Latex','Position',legPos);

% generate zero line
if yMin < 0 && yMax > 0
    negInt = abs(yTicks(1)/yInt);
    PP4   = bm+negInt*(yh/ynInt);
    annotation('line',[lm lm+xw], [PP4 PP4],...
        'LineStyle','--','LineWidth',1);
end

%%% month strings and lines
if shortWaveSwitch==1
    FC     = [0.5 0.5 0.5];
    fa     = 0.2;
    SW = {SW};
    % sort data by month
    monthData = cell(nm,1);
    for i=1:nm
        monthData(i,:) = cellfun(@(x,y) x(y==i) , SW, {monthCount{1}}, 'UniformOutput',false);
    end
    % reshape into diurnal cycles
    monthData = cellfun( @(x) reshape(x,dtd,[]), monthData, 'UniformOutput',false);
    % average by diurnal cycles
    diuCyc = cellfun( @(x) nanmean(x')', monthData, 'UniformOutput',false );
    SW = cell2mat(diuCyc);
    SW    = reshape(SW,nm*dtd,[]);
    tneu = SW;
    tneu(SW<0.5) = 0;
    tneu(SW>0.5) = 1;
    ttt = [NaN;tneu(2:end)-tneu(1:end-1)];

    tNorm = t-t(1);
    tNorm = tNorm/nm*xw;
    tPos  = tNorm(ttt~=0);
    tPos(end+1) = 1; 
end

% percent string positions
perStrPos = [0 0.018 0.038 0.056]';
for i=1:nm
    % Lines
    if i==1
        xP = lm+(i-1)*mw;
        annotation('line',[xP xP],[bm bm+yh+wrs*mw],'LineWidth',1.1)
    end
    
    if i==nm-1
        xP = lm+i*mw;
        if perBarSwitch == 0
            annotation('line',[xP xP],[bm 0.88*legY],'LineWidth',1.1)
        elseif numel(legLines)==4
            annotation('line',[xP xP],[bm 0.88*legY],'LineWidth',1.1)    
        elseif mean(ismember('\ $G',legendStr{1}))<1
            annotation('line',[xP xP],[bm 0.9*legY],'LineWidth',1.1)
        elseif mean(ismember('\ $G',legendStr{1}))==1
            annotation('line',[xP xP],[1.45*legY bm+yh+wrs*mw],'LineWidth',1.1)
        end
        if windRoseSwitch == 1
            annotation('line',[xP xP],[bm+yh bm+yh+wrs*mw],'LineWidth',1.1)
        end
    else
        xP = lm+i*mw;
        annotation('line',[xP xP],[bm bm+yh+wrs*mw],'LineWidth',1.1)
    end
    
    
    % month String (Jan, Feb, Mar, ...)
    annotation('textbox',[xP-mw bm-0.11 mw mw],'String',monthStr{i},...
        'LineStyle','none','HorizontalAlignment','center','FontSize',fontSize-2,'FontName','Arial');
    
    if perBarSwitch == 1;
        % Number of days as percent bars
        if windRoseSwitch == 1
            hBox = 0.0333333333;    % height of the box
        else
            hBox = 0.0416666666;    % height of the box
        end
        
        pBox = bm+yh-hBox;  % y-position of the box
        % Background
        annotation('textbox',[xP-mw pBox mw hBox],'BackGroundColor','none')
        annotation('textbox',[xP-mw pBox mw hBox/2],'BackGroundColor',Col{1},...
            'FaceAlpha',0.2,'EdgeColor','none')
        annotation('textbox',[xP-mw pBox+hBox/2 mw hBox/2],'BackGroundColor',Col{2},...
            'FaceAlpha',0.2,'EdgeColor','none')
        % Percent Bar
        annotation('textbox',[xP-mw pBox diuFraction(i,1)*mw hBox/2],'BackGroundColor',Col{1},...
            'FaceAlpha',0.5,'EdgeColor','none')
        annotation('textbox',[xP-mw pBox+hBox/2 diuFraction(i,2)*mw hBox/2],'BackGroundColor',Col{2},...
            'FaceAlpha',0.5,'EdgeColor','none')
        % Percent Number
        if windRoseSwitch == 1
            for k=1:length(percentStr)
                annotation('textbox',[xP-mw+perStrPos(k) pBox mw hBox],'String',percentStr{k},...
                    'FontSize',fontSize-9,'FontName','Arial',...
                    'VerticalAlignment','bottom','LineStyle','none')
            end
        else
            for k=1:length(percentStr)
                annotation('textbox',[xP-mw+perStrPos(k) pBox mw hBox],'String',percentStr{k},...
                    'FontSize',fontSize-13,'FontName','Arial',...
                    'VerticalAlignment','bottom','LineStyle','none')
            end
        end
    end
end
if shortWaveSwitch==1
kLoop = 1:2:length(tPos)-1;
for k=kLoop
        if k==kLoop(end-1)
            if perBarSwitch == 0
                p2 = bm;
                p4 = 0.88*legY-bm;
            elseif numel(legLines)==4
                p2 = bm;
                p4 = 0.88*legY-bm;
            elseif mean(ismember('\ $G',legendStr{1}))<1
                p2 = bm;
                p4 = 0.9*legY-bm;
            else
                p2 = 1.45*legY;
                p4 = bm+yh+wrs*mw - 1.45*legY;
            end

        else
            p2 = bm;
            if perBarSwitch==1
                p4 = yh-hBox;
            else
                p4 = yh;
            end
        end
        if p4<0
            p4=0;
        end
        if k==length(tPos)-1
            recPos = [lm+tPos(k) p2 1-lm-tPos(k)-rm p4];
        else
            recPos = [lm+tPos(k) p2 tPos(k+1)-tPos(k) p4];
        end
            
            annotation('rectangle',recPos,'Color','none',...
                'FaceColor',FC,'FaceAlpha',fa,'units','normalized');
    end
end

% xLabel
%xLabelStr = {'Time of Day in Hours'};
xLabelStr = {''};

% xlabel(xLabelStr,'Interpreter','Latex','FontSize',fontSize);
% get(gca,''

xLab = xlabel(xLabelStr);
xLab.Interpreter = 'latex';
xLab.Units = 'normalized';
xLab.FontSize = fontSize;
xLabPos = xLab.Position;
xLabPos(2) = xLabPos(2)-0.05;
xLab.Position = xLabPos;

% yLabel
if windRoseSwitch == 1
    yLab = ylabel(yLabelStr);
    yLab.Interpreter = 'latex';
    yLab.FontSize = fontSize-2;
else
    text(0,1.035,yLabelStr,...
        'units','normalized','Interpreter','latex','Fontsize',fontSize-3)
end

% Titel
if titleSwitch == 1 
    set(gcf,'NextPlot','add');
    axes;
    h = title(titleStr,'Interpreter','latex','FontSize',fontSize+2);
    set(gca,'Visible','off');
    set(h,'Visible','on');
end


fid = fopen([expID '.txt'], 'w');
for j=1:nw
    % Ausgabe mit max. 6 Vorkomma- und 2 Nachkommastellen
    fprintf(fid, '%6.3f\t',mean(reshape(diuData{j,1},numel(diuData{j,1}),[])));
end
if fluxScaleSwitch==1
    fprintf(fid, '%6.3f\t',mean(reshape(scaleData,numel(scaleData),[])));
end
    
% Datei schliessen
fclose(fid);

if windRoseSwitch == 1
    % Create textbox
    annotation(gcf,'textarrow',[0.018 0.018],[0.9 0.9],'string','Wind Speed in m/s', ...
          'HeadStyle','none','LineStyle', 'none', 'TextRotation',90,'Interpreter','Latex','FontSize',fontSize-8);
end


% print
if errorBars == 1
    drucken([outputDir expID '.std.png']);
else
    drucken([outputDir expID '.png']);
end

end




