if max(size(time_BG))>4

    plot_BG=0;


    BG_count=0;
    clear BG_cycle
    clear time_BG_cycle
    
    %This piece of code identifies the number of BG cycles    
    idx=find(diff(time_BG)>0.0035); %at least 5 min apart - refers to the last point of BG cycle, i.e., the point before the gap;
    
	%% Runs through BG data gaps and calculates averages between them
    for i=1:max(size(idx))
		BG_count=BG_count+1;
        if i==1 % Calculating the first BG cycle from 1 to the first gap 
			BG_cycle(BG_count,:)=nanmean(BG(1:idx(1),:),1); 
			time_BG_cycle(BG_count,:)=nanmean(time_BG(1:idx(1)));
		else % Ranging through the rest of the data gaps
			BG_cycle(BG_count,:)=nanmean(BG(idx(i-1)+1:idx(i),:),1);%averages skipping the first minute of BG check
			time_BG_cycle(BG_count,:)=nanmean(time_BG(idx(i-1)+1:idx(i)));%averages skipping the first minute of BG check 
			
			%This part I'll try to subtrack BG from vocs on-the-fly
			clear BG_dummy
			idx_ptrms=find(time_PTRMS>time_BG_cycle(BG_count-1)&time_PTRMS<time_BG_cycle(BG_count));
			for j=1:size(vocs,2)
				BG_dummy(:,j)=interp1(time_BG_cycle(BG_count-1:BG_count),BG_cycle(BG_count-1:BG_count,j),time_PTRMS(idx_ptrms));
            end
%             close all
%             circ(time_BG_cycle(BG_count-1:BG_count),BG_cycle(BG_count-1:BG_count,Idx79))
%             hold on
%             circ(time_BG,BG(:,Idx79),'g')
%             dynamicdateticks([],[],'dd/mm')
%             circ(time_PTRMS,vocs(:,Idx79),'r')
%             circ(time_PTRMS(idx_ptrms),BG_dummy(:,Idx79),'k')
			vocs(idx_ptrms,:)=vocs(idx_ptrms,:)-BG_dummy;
%             circ(time_PTRMS,vocs(:,Idx79),'c')
%             axis([time_BG_cycle(BG_count-1) time_BG_cycle(BG_count) 0 max(vocs(idx_ptrms,Idx79))])
%             legend('BG cycles','BG data points','measurement','interpolated BG','subtracted measurement')
        end     
    end
    
	%% Identifying the last cycle, after the last data gap
	BG_count=BG_count+1; 
	BG_cycle(BG_count,:)=nanmean(BG(idx(end)+1:end,:),1); 
	time_BG_cycle(BG_count,:)=nanmean(time_BG(idx(end)+1:end));
	
	clear BG_dummy
	idx_ptrms=find(time_PTRMS>time_BG_cycle(BG_count-1)&time_PTRMS<time_BG_cycle(BG_count));
	for j=1:size(vocs,2)
		BG_dummy(:,j)=interp1(time_BG_cycle(BG_count-1:BG_count),BG_cycle(BG_count-1:BG_count,j),time_PTRMS(idx_ptrms));
	end
% 	 close all
%     circ(time_BG_cycle(BG_count-1:BG_count),BG_cycle(BG_count-1:BG_count,Idx79))
%     hold on
%     circ(time_BG,BG(:,Idx79),'g')
%     dynamicdateticks([],[],'dd/mm')
%     circ(time_PTRMS,vocs(:,Idx79),'r')
%     circ(time_PTRMS(idx_ptrms),BG_dummy(:,Idx79),'k')
    vocs(idx_ptrms,:)=vocs(idx_ptrms,:)-BG_dummy;
%     circ(time_PTRMS,vocs(:,Idx79),'c')
%     axis([time_BG_cycle(BG_count-1) time_BG_cycle(BG_count) 0 max(vocs(idx_ptrms,Idx79))])
%     legend('BG cycles','BG data points','measurement','interpolated BG','subtracted measurement')
	
	%% Removing BG for measurement cycles before and after the first and last BG cycles (if existing)
	idx_ptrms=find(time_PTRMS<time_BG_cycle(1));
	if min(size(idx_ptrms))>0
        close all
%         circ(time_BG_cycle(1),BG_cycle(1,Idx79))
%         hold on
%         circ(time_BG,BG(:,Idx79),'g')
%         dynamicdateticks([],[],'dd/mm')
%         circ(time_PTRMS,vocs(:,Idx79),'r')
%         circ(time_PTRMS(idx_ptrms),vocs((idx_ptrms),Idx79),'k')
        for j=1:size(vocs,2)
            vocs(idx_ptrms,j)=vocs(idx_ptrms,j)-BG_cycle(1,j);
        end
%         circ(time_PTRMS,vocs(:,Idx79),'c')
%         axis([time_BG_cycle(BG_count-1) time_BG_cycle(BG_count) 0 max(vocs(idx_ptrms,Idx79))])
%         legend('BG cycles','BG data points','measurement','interpolated BG','subtracted measurement')
	end
	
	idx_ptrms=find(time_PTRMS>time_BG_cycle(end));
	if min(size(idx_ptrms))>0
        close all
        circ(time_BG_cycle(end),BG_cycle(end,Idx69))
        hold on
        circ(time_BG,BG(:,Idx79),'g')
        dynamicDateTicks([],[],'dd/mm')
        circ(time_PTRMS,vocs(:,Idx79),'r')
        circ(time_PTRMS(idx_ptrms),vocs((idx_ptrms),Idx69),'k')
        for j=1:size(vocs,2)
            vocs(idx_ptrms,j)=vocs(idx_ptrms,j)-BG_cycle(end,j);
        end
        circ(time_PTRMS,vocs(:,Idx69),'c')
        axis([time_BG_cycle(BG_count-1) time_BG_cycle(BG_count) 0 max(vocs(idx_ptrms,Idx69))])
        legend('BG cycles','BG data points','measurement','interpolated BG','subtracted measurement')
	end
	
	%% Plotting
     if plot_BG==1 && isunix==0
        for i=1:max(size(header))
            figure
            eval(['logy_circ(time_PTRMS,vocs(:,Idx',num2str(header(i)),'))']);
            hold on
            eval(['logy_circ(time_BG,BG(:,Idx',num2str(header(i)),'),''r'')']);
            eval(['logy_circ(time_BG_cycle,BG_cycle(:,Idx',num2str(header(i)),'),''g'')']);
            dynamicDateTicks([], [], 'dd/mm')
            eval(['title(label_',(num2str(header(i))),')'])
            legend('Measurement','BG','Avg BG')
%             eval(['print -djpeg 1_pics\BG_Meas_',num2str(header(i)),'.jpg'])
        end;
    end;

%{
    i=11;
    figure
    eval(['logy_circ(time_PTRMS,vocs(:,Idx',num2str(header(i)),'))']);
    hold on
    eval(['logy_circ(time_BG,BG(:,Idx',num2str(header(i)),'),''r'')']);
    eval(['logy_circ(time_BG_cycle,BG_cycle(:,Idx',num2str(header(i)),'),''g'')']);
    dynamicDateTicks([], [], 'dd/mm')
    eval(['title(label_',(num2str(header(i))),')'])

    count=0;
    %This piece of code removes the BG cycles from the measurement data
    if BG_count==1
        for i=1:max(size(header))
            vocs(:,i)=vocs(:,i)-BG_cycle(1,i);
        end;
    elseif BG_count>1
        for j=1:max(size(time_PTRMS))
            if time_PTRMS(j)<time_BG_cycle(1);
                vocs(j,:)=vocs(j,:)-BG_cycle(1,:);
            elseif time_PTRMS(j)>time_BG_cycle(end);
                vocs(j,:)=vocs(j,:)-BG_cycle(end,:);
            else
                [diff, idx]=min(abs(time_PTRMS(j) - time_BG_cycle));
                if time_PTRMS(j)>time_BG_cycle(idx)
                    idx_before=idx;
                    idx_after=idx+1;
                else
                    idx_before=idx-1;
                    idx_after=idx;
                end;

                tot_time=time_BG_cycle(idx_after)-time_BG_cycle(idx_before);
                count=count+1;
                BG_prop(count,:)=BG_cycle(idx_after,:).*(time_PTRMS(j)-time_BG_cycle(idx_before))./tot_time + BG_cycle(idx_before,:).*(time_BG_cycle(idx_after)-time_PTRMS(j))./tot_time;        
                vocs(j,:)=vocs(j,:)-BG_prop(count,:); %Simply subtracts the first BG check
                time_BG_prop(count)=time_PTRMS(j);
            end;                

        end;

  % logy_circ(time_BG_prop,BG_prop(:,11),'m');
  % legend('Measurement','BG','Avg BG','Prop BG')
    end;
%}
end;
