function [ average_x_values_jump, x_values_jump ] = quantile_boom_bust( num_var, num_jump, jump_matrix, x )


% PURPOSE: CALCULATE AVERAGE AND TWO STANDARD DEVIATION INTERVAL IN BOOM/BUST
% EPISODES IN THE PERIOD OF THE BOOM/BUST ONE AND TWO PERIODS BEFORE AND ONE 
% AND TWO PERIODS AFTER THE EPISODES
% 
% INPUT:
%     -num_var: NUMBER OF TIME SERIES FOR THE VARIABLE FOR WHICH WE CALCULATE THE AVERAGE
%     -num_jump: NUMBER OF BOOM OR BUST
%     -jump_matrix: MATRIX WITH ALL BOOM/BUST EPISODES, YEAR AND COUNTRY WHERE THEY HAPPEN
%     -x: VARIABLE FOR WHICH WE CALCULATE THE AVERAGE AND INTERVAL
%     
% OUTPUT:
%     -average_x_values_jump: AVERAGE FOR THE VARIABLE IN CASE OF A TYPOLOGY OF EPISODE
%     -std_x_values_jump: STANDARD DEVIATION FOR THE VARIABLE IN CASE OF A TYPOLOGY OF EPISODE
%     -upper_bound_jump: UPPER BOUND OF THE INTERVAL CONFIDENCE FOR THE VARIABLE IN CASE OF A TYPOLOGY OF EPISODE
%     -lower_bound_jump: LOWER BOUND OF THE INTERVAL CONFIDENCE FOR THE VARIABLE IN CASE OF A TYPOLOGY OF EPISODE
%     -x_values_jump: MATRIX WITH ALL VALUES OF THE VARIABLE IN ALL CASES OF A TYPOLOGY OF EPISODES DIVIDED BY COUNTRY


    %Create matrix for jumps (either booms or bustes)
    x_values_jump=nan(num_jump,num_var);
    for j=1:num_var
        for i=1:num_jump
            if jump_matrix(i,2)==1
                x_values_jump (i,:,j)=[jump_matrix(i,1); nan; nan; x(jump_matrix(i,1),jump_matrix(i,2),j); x(jump_matrix(i,1),jump_matrix(i,2)+1,j); x(jump_matrix(i,1),jump_matrix(i,2)+2,j)];
            elseif jump_matrix(i,2)==2
                x_values_jump (i,:,j)=[jump_matrix(i,1); nan; x(jump_matrix(i,1),jump_matrix(i,2)-1,j); x(jump_matrix(i,1),jump_matrix(i,2),j); x(jump_matrix(i,1),jump_matrix(i,2)+1,j); x(jump_matrix(i,1),jump_matrix(i,2)+2,j)];
            elseif jump_matrix(i,2)==16
                x_values_jump (i,:,j)=[jump_matrix(i,1); x(jump_matrix(i,1),jump_matrix(i,2)-2,j); x(jump_matrix(i,1),jump_matrix(i,2)-1,j); x(jump_matrix(i,1),jump_matrix(i,2),j); x(jump_matrix(i,1),jump_matrix(i,2)+1,j); nan];    
            elseif jump_matrix(i,2)==17
                x_values_jump (i,:,j)=[jump_matrix(i,1); x(jump_matrix(i,1),jump_matrix(i,2)-2,j); x(jump_matrix(i,1),jump_matrix(i,2)-1,j); x(jump_matrix(i,1),jump_matrix(i,2),j); nan; nan];    
            else
                x_values_jump (i,:,j)=[jump_matrix(i,1); x(jump_matrix(i,1),jump_matrix(i,2)-2,j); x(jump_matrix(i,1),jump_matrix(i,2)-1,j); x(jump_matrix(i,1),jump_matrix(i,2),j); x(jump_matrix(i,1),jump_matrix(i,2)+1,j); x(jump_matrix(i,1),jump_matrix(i,2)+2,j)];
            end
        end
    end
 
    %Classify by percentiles
    for i=1:size(x,1)
        if x(i,1)<=threshold_high
        quantile_indicator(i,1)=2;
        elseif aggregate_nfa(i,1)>=threshold_low
        quantile_indicator(i,1)=0;
        else
        quantile_indicator(i,1)=1;
        end
    end
    
    
    
    %Average for the different time horizons
    average_x_values_jump=nan(1,5,num_var);
    
    for j=1:num_var
        for i=1:5
            average_x_values_jump(1,i,j)=nanmean(x_values_jump(:,i+1,j),1);
        end
    end    

end