function [C,Res,rL,cL,cGrn,rGrn,rNum,cNum,wghN,unwN,SIGC,waveGr,CD] = bigtab(row,col,statistics,layout,ltperc,srtvl,options,boxes,marks,valid_values, index,select, footer, label, decimal)
%ovo je preimenovana fukcija bigtabSG3
% sada radi sa mean(varijabla)
%v5 - dodata i populacija
%v4 - popravljen i ncase koji neradi kako treba sa filterom
%TODO - labele za blok multiple -npr.  nivo ispitanika, nivo domacinstva

%statistics koje su definisane
%count(1) rowpct(2) colpct(3) totpct(4) mean(5) trimmean(5.1) median(6)
%min(7) max(8) range(9) sd(10) var(11) ciXX(12) trimciXX(12.1)
%index(13) colshare(14) rowshare(15) trimcolshare(14.1) trimrowshare(15.1)
%poppct(17) popXY(16)-gde je XY broj u hiljadama ili globalna varijabla popsize

%summarize
%statistics
%layout
%ltperc
%options
    SIGC = [];
    
    if ~def(boxes)
        boxes = "";
    end
    
	global gWave
	if def(wave)
		gWave = ucase(wave(1));
	else
		gWave = "WAVE";
	end
	if def(valid_values)
		if size(valid_values,1)>size(valid_values,2)
			validValues = valid_values';
		else
			validValues = valid_values;
		end
	else
		validValues = [];
	end

	if def(row)
		row=resolveSet(ucase(row));
		indx = find(~strlen(row));
		if any(indx)
			row(indx)=[];
		end
		vars_1 = max(size(row));
		if vars_1
			rowDonje = nestpars(row);
			rowBM = DataGetBMVar(rowDonje);
		end
	else
		vars_1 = 0;
	end

	nivoaBM = 0;
	col_1 = 1;
	global g_col;
	if def(col)
	  i = 0;
	elseif g_col~=[]
	  col = g_col;
	else
		col_1 = 0;
    end

	if col_1
	  col=resolveSet(ucase(col));
		indx = find(~strlen(col));
		if any(indx)
			col(indx)=[];
		end
		col_1 = max(size(col));
		colDonje = nestpars(col);
		colBM = DataGetBMVar(colDonje);
		%provera da li su sa istog nivoa
		if vars_1
			nivoaBM = proveraBM(rowDonje, rowBM', "redovima", colDonje, colBM, "kolonama");
		end
    end


	global g_col1,g_col2,g_col3;
	colt = 0;
	for cj=1:3
		if def(col#cj)
			col#cj = col#cj;
		elseif size(g_col#cj,1)
		  col#cj = g_col#cj;
		else
			col#cj = [];
		end

		if size(col#cj,1)==0
			col#cj="";
		else
			if colt<cj-1;error("Argument 'col"+num2str("%d",cj-1)+"' is not specified");end
			col#cj=resolveSet(ucase(col#cj));
			colt = colt+1;
		end
	end

	printTab = 1;
	sampleSize = 0;
	ownbase = 0;
	outzero = 0;
	plus = 0;
	indexcolor = 0;
	rlabvals = 0;
	clabvals = 0;
	sysmis = 0;
	whole = 0;
    sigc = 0;
    nobonf = 0;
    sigw = 0;
    nototal = 0;
    color3 = 0;
	abc = 0;
	opt_nps = 0;
	
	if def(options)
		options=ucase(options);
		if any(options=="NOPRINT")
			printTab = 0;
		end
		if any(options=="SAMPLESIZE")
			sampleSize = 1;
		end
		if any(options=="OWNBASE")
			ownbase = 1;
		end
		if any(options=="NOZEROCOLUMNS")
			outzero = 1;
		end
		if any(options=="PLUS")
			plus = 1;
		end
		if any(options=="INDEXCOLOR")
			indexcolor = 1;
		end
		if any(options=="LABVALS")
			rlabvals = 1;
			clabvals = 1;
		end
		if any(options=="RLABVALS")
			rlabvals = 1;
		end
		if any(options=="CLABVALS")
			clabvals = 1;
		end
		if any(options=="SYSMIS")
			sysmis = 1;
		end
		if any(options=="WHOLE")
			whole = 1;
		end
        if any(options=="SIGC")
			sigc = 1;
		end
        if any(options=="NOBONF")
			nobonf = 1;
		end
        if any(options=="SIGW")
			sigw = 1;
		end
        if any(options=="NOTOTAL")
			nototal = 1;
		end
        if any(options=="COLOR3")
			color3 = 1;
		end
		if any(options=="ABC")
			abc = 1;
		end
		if any(options=="NPS")
			opt_nps = 1;
		end
	end

	if def(ltperc)
		if size(ltperc,1)
			if isstr(ltperc)
				ltperc = str2num(ltperc);
			end
			if prod(size(ltperc))~=2
				error("Argument 'ltperc' must have two parameters")
			end
		else
			ltperc = [];
		end
	else
		ltperc = [];
	end

	if def(marks)
		if max(size(marks))>1
			error("Marks have to be scalar value")
		else
			marks = str2num(marks);
			if marks<=0
				error("Marks have to be greater then zero")
			end
		end
	else
		marks = 0;
	end

	if def(index) %%%% ne znam sta ovo size         & size(index,1)
		indx = find(~strlen(index));
		if any(indx)
			index(indx)=[];
		end
		index = ucase(index);
		if prod(size(index))==1
			if ~strcmp(index,"$TOTAL")
				index_num = Inf;
			else
				error("Argument 'index' have to have to parameters, variable name and value")
			end
		elseif prod(size(index))==0
			index_num = NaN;
		elseif prod(size(index))>2
			error("Argument 'index' have to have to parameters, variable name and value")
		else
			index_num = str2num(index(2));
			index = index(1);
			indexBM = DataGetBMVar(index);
			%provera da li su sa istog nivoa
			if vars_1
				i = proveraBM(rowDonje, rowBM', "rows", index, indexBM, "index");
				if nivoaBM
					if i~=nivoaBM
						error("Different variable level in rows and in index from those in columns")
					end
				end
				nivoaBM = i;
			end
			if col_1
				i = proveraBM(colDonje, colBM', "kolonama", index, indexBM, "indexu");
				if nivoaBM
					if i~=nivoaBM
						error("Different variable level in columns and in index from those in rows")
					end
				end
				nivoaBM = i;
			end
		end
	else
		index_num = NaN;
	end




	if def(statistics)
		statistics = ucase(statistics);
		statistics = statistics(:);
		if size(statistics,2)>size(statistics,1)
			statistics = statistics';
		end
		if index_num ~= NaN
			if ~any(statistics=="INDEX")
				statistics = [statistics; "INDEX"];
			end
		end
		[statistics,summarize] = parseS(statistics);

		mInd = make(size(summarize),0);
		[sum_vars,indx]=freq(summarize);
		sumv_1 = size(sum_vars,1);
		if sum_vars(sumv_1)==""
			sum_vars(sumv_1) = [];
			sumv_1 = sumv_1-1;
		end
		for i=1:sumv_1
			indx = find(summarize==sum_vars(i));
			mInd(indx)=i;
		end
		stat_1 = size(statistics,1);
	else
		if index_num ~= NaN
			statistics = "INDEX";
		else
			statistics = "COLPCT";
		end
		stat_1 = 1;
		sumv_1 = 0;
		mInd = 0;
		sum_v = 0;
	end
	stat_type =  make(size(statistics),0);
	if ~def(srtvl)
		srtvl = Inf;
	end


	if sumv_1
		sumBM = DataGetBMVar(sum_vars);
		%provera da li su sa istog nivoa

		if vars_1
			i = proveraBM(rowDonje, rowBM', "redovima", sum_vars, sumBM', "sumiranim");
			if nivoaBM
				if i~=nivoaBM
					error("Different variable level in rows and summerized variables from those in columns")
				end
			end
			nivoaBM = i;
		end
		if col_1
			i = proveraBM(colDonje, colBM', "kolonma", sum_vars, sumBM', "sumiranim");
			if nivoaBM
				if i~=nivoaBM
					error("Different variable level in columns and summerized variables from those in rows")
				end
			end
			nivoaBM = i;
		end
	end


	%ovo je broj statistika za summarize varijable
	%sumstats =  make(size(statistics),0);
	indx = find(~strcmp(statistics,"COUNT"));
	if any(indx)
		stat_type(indx) = 1;
		statistics(indx) = "Count";
		if any(mInd(indx)~=0)
			error("Count and % do not have variables is parantheses'")
		end
	end
	indx = find(~strcmp(statistics,"ROWPCT"));
	if any(indx)
		stat_type(indx) = 2;
		statistics(indx) = "Row %";
		if any(mInd(indx)~=0)
			error("Count and % do not have variables is parantheses'")
		end
	end
	indx = find(~strcmp(statistics,"COLPCT"));
	if any(indx)
		stat_type(indx) = 3;
		statistics(indx) = "Col %";
		if any(mInd(indx)~=0)
			error("Count and % do not have variables is parantheses'")
		end
	end
	indx = find(~strcmp(statistics,"TOTPCT"));
	if any(indx)
		stat_type(indx) = 4;
		statistics(indx) = "Tot %";
		if any(mInd(indx)~=0)
			error("Count and % do not have variables is parantheses'")
		end
	end
	indx = find(~strcmp(statistics,"MEAN"));
	if any(indx)
		stat_type(indx) = 5;
		statistics(indx) = "Mean";
		if any(mInd(indx)==0)
			error("Variable for summerzied statistic is not defined")
		end
	end
	indx = find(~strcmp(statistics,"TRIMMEAN"));
	if any(indx)
		stat_type(indx) = 5.1;
		statistics(indx) = "TMean";
		if any(mInd(indx)==0)
			error("Variable for summerzied statistic is not defined")
		end
	end
	indx = find(~strcmp(statistics,"MEDIAN"));
	if any(indx)
		stat_type(indx) = 6;
		statistics(indx) = "Median";
		if any(mInd(indx)==0)
			error("Variable for summerzied statistic is not defined")
		end
	end
	indx = find(~strcmp(statistics,"MIN"));
	if any(indx)
		stat_type(indx) = 7;
		statistics(indx) = "Min";
		if any(mInd(indx)==0)
			error("Variable for summerzied statistic is not defined")
		end
	end
	indx = find(~strcmp(statistics,"MAX"));
	if any(indx)
		stat_type(indx) = 8;
		statistics(indx) = "Max";
		if any(mInd(indx)==0)
			error("Variable for summerzied statistic is not defined")
		end
	end
	indx = find(~strcmp(statistics,"RANGE"));
	if any(indx)
		stat_type(indx) = 9;
		statistics(indx) = "Range";
		if any(mInd(indx)==0)
			error("Variable for summerzied statistic is not defined")
		end
	end
	indx = find(~strcmp(statistics,"SD"));
	if any(indx)
		stat_type(indx) = 10;
		statistics(indx) = "Std";
		if any(mInd(indx)==0)
			error("Variable for summerzied statistic is not defined")
		end
	end
	indx = find(~strcmp(statistics,"VAR"));
	if any(indx)
		stat_type(indx) = 11;
		statistics(indx) = "Var";
		if any(mInd(indx)==0)
			error("Variable for summerzied statistic is not defined")
		end
	end
	indx = sortd(find(~strcmp(left(statistics,2),"CI")));
	if any(indx)
		if any(strlen(statistics(indx))>5 | strlen(statistics(indx))<3)
			error("Incorrect definition of percentage values in CI")
		end
		ci = make(size(statistics),NaN);
		if any(strlen(statistics(indx))==4)
			ci(indx)=-str2num(right(statistics(indx),2));
		end
		if any(strlen(statistics(indx))==3)
			ci(indx)=-str2num(right(statistics(indx),1));
		end
		indx1 = find(ci==NaN);
		if any(indx1)
			ci(indx1)=[];
		end
		stat_type(indx) = 12;
		if any(mInd(indx)==0)
			error("Variable for summary statistic is not defined")
		end
		for i=1:size(indx,1)
			if indx(i)<size(statistics,1)
				if outputType
					statistics = [statistics(1:indx(i));"|";statistics(indx(i)+1:NaN)];
				else
					statistics = [statistics(1:indx(i));"";statistics(indx(i)+1:NaN)];
				end
				stat_type = [stat_type(1:indx(i));ci(i);stat_type(indx(i)+1:NaN)];
				mInd = [mInd(1:indx(i));0;mInd(indx(i)+1:NaN)];
			else
				if outputType
					statistics = [statistics;"|"];
				else
					statistics = [statistics;""];
				end
				stat_type = [stat_type;ci(i)];
				mInd = [mInd;0];
			end
		end
		stat_1 = stat_1 + size(indx,1);
	end

	indx = sortd(find(~strcmp(left(statistics,6),"TRIMCI")));
	if any(indx)
		if any(strlen(statistics(indx))>9 | strlen(statistics(indx))<7)
			error("Incorrect definition of percntage values in CI")
		end
		ci = make(size(statistics),NaN);
		if any(strlen(statistics(indx))==8)
			ci(indx)=-str2num(right(statistics(indx),2));
		end
		if any(strlen(statistics(indx))==7)
			ci(indx)=-str2num(right(statistics(indx),1));
		end
		indx1 = find(ci==NaN);
		if any(indx1)
			ci(indx1)=[];
		end
		stat_type(indx) = 12.1;
		if any(mInd(indx)==0)
			error("Variable for summerzied statistic is not defined")
		end
		for i=1:size(indx,1)
			if indx(i)<size(statistics,1)
				if outputType
					statistics = [statistics(1:indx(i));"|";statistics(indx(i)+1:NaN)];
				else
					statistics = [statistics(1:indx(i));"";statistics(indx(i)+1:NaN)];
				end
				stat_type = [stat_type(1:indx(i));ci(i);stat_type(indx(i)+1:NaN)];
				mInd = [mInd(1:indx(i));0;mInd(indx(i)+1:NaN)];
			else
				if outputType
					statistics = [statistics;"|"];
				else
					statistics = [statistics;""];
				end
				stat_type = [stat_type;ci(i)];
				mInd = [mInd;0];
			end
		end
		stat_1 = stat_1 + size(indx,1);
	end

	indx = find(~strcmp(statistics,"INDEX"));
	if any(indx)
		stat_type(indx) = 13;
		statistics(indx) = "Index";
		if any(mInd(indx)~=0)
			error("Count and % do not have variables is parantheses'")
		end
	end

	indx = find(~strcmp(statistics,"COLSHARE"));
	if any(indx)
		stat_type(indx) = 14;
		statistics(indx) = "Col Share";
		%if sumv_1
		%	sumstats(indx) = 1;
		%end
	end

	indx = find(~strcmp(statistics,"TRIMCOLSHARE"));
	if any(indx)
		stat_type(indx) = 14.1;
		statistics(indx) = "TColShare";
		if any(mInd(indx)==0)
			error("Variable for summerzied statistic is not defined")
		end
	end

	indx = find(~strcmp(statistics,"ROWSHARE"));
	if any(indx)
		stat_type(indx) = 15;
		statistics(indx) = "Row Share";
		%if sumv_1
		%	sumstats(indx) = 1;
		%end
	end

	indx = find(~strcmp(statistics,"TRIMROWSHARE"));
	if any(indx)
		stat_type(indx) = 15.1;
		statistics(indx) = "TRowShare";
		if any(mInd(indx)==0)
			error("Variable for summerzied statistic is not defined")
		end
	end

	indx = find(~strcmp(statistics,"POPPCT"));
	if any(indx)
		stat_type(indx) = 17;
		statistics(indx) = "Pop %";
		if any(mInd(indx)~=0)
			error("Count and % do not have variables is parantheses'")
		end
	end



	indx = find(~strcmp(left(statistics,3),"POP"));
	if any(indx)
		stat_type(indx) = 16;
		a = strlen(statistics(indx));
		indx1 = find(a>3);
		if any(indx1)
			popSize = str2num(right(statistics(indx(indx1(1))),	a(indx1(1))-3))/100;
		else
			global popsize
			if max(size(popsize))
				popSize = popsize(1)/100;
			else
				error("Size of population is not defined either in statistics 'pop' nor in the global variable 'popsize'")
			end
		end
		statistics(indx) = "in 000";
		if any(mInd(indx)~=0)
			error("Count and % do not have variables is parantheses'")
		end
	end

	%sumstats_ = msum(sumstats);
	indx = find(~stat_type);
	if any(indx)
		error("Unknown statistics: "+ statistics(indx(1)))
	end

	if def(layout)
		layout=ucase(layout);
		if any(~strcmp(layout,"STATISTICSINCOL"))
			statsinrow=0;
		else
			statsinrow=1;
		end
		%if any(~strcmp(layout,"SUMMARIZEINCOL"))
		%	sumsinrow=0;
		%else
		%	sumsinrow=1;
		%end
	else
		%sumsinrow=1;
		statsinrow=1;
	end
	sumsinrow = statsinrow;


	if size(filter,1)
		_ncase = sum(filter);
	else
		_ncase = ncase;
	end

	lRF = NaN;
	lDK = NaN;
	tot_lok_mis = [];
	global g_missing;
	if size(g_missing,1)
		lok_mis = g_missing;
		bt_miss;
	end
	if def(missing)
		lok_mis = missing;
		bt_miss;
	end
	%%%%%%%%%%%%%%%%%%%%%%%% Selector
	%selval = NaN;
	%if def(select)
	%	select=ucase(select);
	%	if max(size(select))>1
	%		selval = str2num(select(2:NaN));
	%	end
	%	select = select(1);
	%
	%	global Selector
	%	if nivoaBM>1
	%		name1 = GetVarName(select);
	%		Selector = make(ncase*nivoaBM,1,1);
	%		for i=1:nivoaBM
	%			if ~vext(name1(i))
	%				error("Varijabla "+name1(i)+" ne postoji, koja je deo "+select)
	%			end
	%			Selector(((i-1)*ncase+1):i*ncase)=Data(:,name1(i));
	%		end
	%	else
	%		Selector = Data(:,select);
	%	end
	%	if selval(1)==NaN
	%		Selector = Selector>0;
	%	else
	%		if size(selval,1)>size(selval,2)
	%			selval = selval';
	%		end
	%		Selector = rsum(Selector==make(size(Selector,1),1,selval))>0;
	%	end
	%end
	if def(select)
		global Selector
		if nivoaBM>1
			SEL = calcFilter(select);
			Selector = make(ncase*nivoaBM,1,1);
			for i=1:nivoaBM
				Selector(((i-1)*ncase+1):i*ncase)=SEL;
			end
			clear SEL;
		else
			Selector = calcFilter(select);
		end
	end

	if col_1
		varss = col;
		ukljuciMiss;
		[cD, cL, cGrn, cGrnM, cNum, cT, cTT, waveGr] = nestpars(col, [sysmis,0,ownbase,nivoaBM,plus,clabvals,whole], [], [],[]);
		if def(missing)
			SetMissVal(varss, missT, missS)
		end
	else
		[cD,cL,num,cT,cTT] = dummy("$TOTAL", [sysmis,0,ownbase,nivoaBM,plus,clabvals,whole],"",index_num);
		cL = "Total";
		cGrn = [1;1];
		cGrnM= [1;1];
	end
	
	%%%%%%%%%%%%%%%%%%%%%%%%%%
	global Alph, Alph2;
	if abc == 1
		Alph2 = [];
		zstart = 0;
		for z = 1: size(cGrn, 2)
			zind = ( (cGrn(1,z) - zstart) : (cGrn(2,z) - zstart) )';
			zstart = cGrn(2,z);
			xx = Alph(1, zind);
			Alph2 = [Alph2, xx];
		end
	else
		Alph2 = Alph;
	end
	%%%%%%%%%%%%%%%%%%%%%%%%%%%

	if index_num ~= NaN %treba i indx
		[bin,lab,num,binB,binBB] = dummy(index, [sysmis,0,ownbase,nivoaBM,plus,0],"",index_num);
		cD = [bin, cD];
		cT = [binB, cT];
		cTT = [binBB, cTT];
		cL = [make(size(cL,1),1,"Index"),cL];
		cGrn = [[1;1],cGrn+1];
		cGrnM = [[1;1],cGrnM+1];
	end

	if vars_1 | opt_nps
		varss = row;
		ukljuciMiss;
		[rD, rL, rGrn, rGrnM, rNum, rT, rTT] = nestpars(row, [sysmis,marks,ownbase,nivoaBM,plus,rlabvals,whole],ltperc,srtvl,boxes,validValues);
		if def(missing)
			SetMissVal(varss, missT, missS)
		end
		rowTotal = 0;
		if opt_nps
			C = rD;
			Res = rNum;
			return 
		end
	else
		[rD,rL,num,rT,rTT] = dummy("$TOTAL", [sysmis,0,ownbase,nivoaBM,plus,rlabvals,whole]);
		%rD = make(_ncase*nivoaBM,1,1);
		%rT = make(_ncase*nivoaBM,1,0);
		%rT = make(_ncase*nivoaBM,1,0);
		rL = "Total";
		rNum = 1;
		rGrn = [1;1];
		rGrnM = [1;1];
		vars_1 = 1;
		rowTotal = 1;
	end
	
	
	%%%%%%%%%%%%%%%%%%%%%%%
	global g_reScaleWeight;
	if wgh~=[]
		wghOriginal = wgh;
		if size(g_reScaleWeight,1)
			rfndR = find(rNum == g_reScaleWeight);
			if (rfndR ~= [])
				rfnd = rsum(rD(:, rfndR));
			else
				rfnd = 0;
			end
			
			cfndC = find(cNum == g_reScaleWeight);
			if (cfndC ~= [])
				cfnd = rsum(cD(:, cfndC));
			else
				cfnd = 0;
			end
			
			rcfnd = find(rsum(rfnd + cfnd) > 0);
			if (rcfnd ~= [])
				rcff = find((rsum(rD(rcfnd,:)) == 1) | (rsum(cD(rcfnd,:)) == 1));
				if (size(rcff,1))
					rcfnd = rcfnd(rcff);
				else
					rcfnd = [];
				end
			end
			if size(rD,2)>1
				if rfndR~=[]
					rD(:,rfndR) = [];
					rT(:,rfndR) = [];
					rTT(:,rfndR) = [];
					rNum(rfndR) = [];
					rL(:,rfndR) = [];
				end
			end
			if size(cD,2)>1
				if cfndC~=[]
					cD(:,cfndC) = [];
					cT(:,cfndC) = [];
					cTT(:,cfndC) = [];
					cNum(cfndC) = [];
					cL(:,cfndC) = [];
					cci = find (cGrn>=cfndC);
					if (cci~=[])
						cGrn(cci)=cGrn(cci)-1;
					end
				end
			end
			
			if (size(rcfnd, 1))
				tg = (rfnd + cfnd) > 0;
				
				[val, fre] = freq(Data(:,"NW"));
				for zzz = 1:size(val,1)
					nw3 = find(Data(:,"NW") == val(zzz));
					rcfnd = find(Data(:,"NW") == val(zzz) & tg);
					if rcfnd~=[] 
						tmp_wgh = wgh;
						totalwgh = sum(tmp_wgh(nw3));
						izbaci = sum(tmp_wgh(rcfnd));
						tmp_wgh(nw3) = tmp_wgh(nw3) * totalwgh / (totalwgh - izbaci);
						tmp_wgh(rcfnd) = 0;
						wgh = tmp_wgh;
					end
				end
			end
		end
	end

	rT = rT - rD;
	rTT = rTT - rD;
	vars_1 = size(rD,2);
	
	cT = cT - cD;
	cTT = cTT - cD;
	col_1 = size(cD,2);
 
 
	if size(rL,1)
		if outputType
			indx1 = find(~(strcmp(rL,"~")));
			indx2 = find(~(strcmp(rL,"|")));
			if any(indx1)
				rL(indx1)="|";
			end
			if any(indx2)
				rL(indx2)="~";
			end
		end
		rL = rL';
	end

	%irm=1;
	%icm=1;
	if sumv_1
		rLS = summarize;
		i = find(mInd~=0);
		rLnew(i) = GetVarLabel(summarize(i));
		indx = find(strlen(rLnew(i)));
		if any(indx)
			rLS(i(indx)) = rLnew(i(indx));
		end
	end

	sumstats_ = sum(mInd~=0);
	tInd = stat_type;
	%mInd je vec napravljen

	if statsinrow
		if sumstats_
			if outputType
				rLP = make(max(size(statistics)), 2, "|");
			else
				rLP = make(max(size(statistics)), 2, "");
			end
			rLP(1:NaN,1)=statistics(1:NaN);
			rLP(1:NaN,2)=rLS(1:NaN);
			cLP = [];
		else
			rLP = statistics;
			mInd = 0;
			cLP = [];
		end
	else
		if sumstats_
			rLP = [];
			if outputType
				cLP = make(2,max(size(statistics)), "~");
			else
				cLP = make(2,max(size(statistics)), "");
			end
			cLP(1,1:NaN)=statistics(1:NaN);
			cLP(2,1:NaN)=rLS(1:NaN);
		else
			rLP = [];
			cLP =  statistics';
			mInd = 0;
		end
	end

	if size(rLP,1)
		if outputType
			rLnew = make(size(rLP,1)*size(rL,1),size(rLP,2)+size(rL,2),"|");
		else
			rLnew = make(size(rLP,1)*size(rL,1),size(rLP,2)+size(rL,2),"");
		end
		rLnew(1:NaN:size(rLP,1),1:size(rL,2)) = rL;
		if statsinrow
			if sumstats_
				mIndnew = make(size(rLP,1)*size(rL,1),1,NaN);
				for i = 1:size(rLnew,1):size(rLP,1)
					mIndnew(i:(i+size(rLP,1)-1)) = mInd;
				end
				mInd = mIndnew;
			end

			mIndnew = make(size(rLP,1)*size(rL,1),1,NaN);
			for i = 1:size(rLnew,1):size(rLP,1)
				mIndnew(i:(i+size(rLP,1)-1)) = tInd;
			end
			tInd = mIndnew;
		end
		for i = 1:size(rLnew,1):size(rLP,1)
			rLnew(i:(i+size(rLP,1)-1), (size(rLnew,2)-size(rLP,2)+1):size(rLnew,2)) = rLP;
		end
		rL = rLnew;
	end
	if size(cLP,1)
		if outputType
			rLnew = make(size(cLP,1)+size(cL,1),size(cLP,2)*size(cL,2),"~");
		else
			rLnew = make(size(cLP,1)+size(cL,1),size(cLP,2)*size(cL,2),"");
		end
		rLnew(1:size(cL,1),1:NaN:size(cLP,2)) = cL;
		if ~statsinrow
			if sumstats_
				mIndnew = make(size(cLP,2)*size(cL,2),1,NaN);
				for i = 1:size(rLnew,2):size(cLP,2)
					mIndnew(i:(i+size(cLP,2)-1)) = mInd;
				end
				mInd = mIndnew;
			end
			mIndnew = make(size(cLP,2)*size(cL,2),1,NaN);
			for i = 1:size(rLnew,2):size(cLP,2)
				mIndnew(i:(i+size(cLP,2)-1)) = tInd;
			end
			tInd = mIndnew;
		end
		for i = 1:size(rLnew,2):size(cLP,2)
			rLnew((size(rLnew,1)-size(cLP,1)+1):size(rLnew,1), i:(i+size(cLP,2)-1)) = cLP;
		end
		cL = rLnew;
	end
	%if sumstats_
	%	indx = [find(~sumstats);find(sumstats)];
	%	statistics = statistics(indx);
	%	stat_type = stat_type(indx);
	%	sumstats = sumstats(indx);
	%end
	C = make(size(rL,1),size(cL,2),NaN);

	ir = size(cL,2)/max(max(cGrn));
	if ir~=1
		cGrn(1,:) = (cGrn(1,:)-1)*ir + 1;
		cGrn(2,:) = cGrn(2,:)*ir;
		cGrnM(1,:) = (cGrnM(1,:)-1)*ir + 1;
		cGrnM(2,:) = cGrnM(2,:)*ir;
	end
	ir = size(rL,1)/max(max(rGrn));
	if ir~=1
		rGrn(1,:) = (rGrn(1,:)-1)*ir + 1;
		rGrn(2,:) = rGrn(2,:)*ir;
		rGrnM(1,:) = (rGrnM(1,:)-1)*ir + 1;
		rGrnM(2,:) = rGrnM(2,:)*ir;
	end
	%ovde jos za nivoaBM
	if sumv_1
		varss = sum_vars;
		ukljuciMiss;

		a = DataIsMultiple(sum_vars);
		if any(a)
			s_D = make(size(rD,1),sumv_1,NaN);
			if nivoaBM>1 %ovo je sada blok multipla
				for ir = 1:sumv_1
					a = GetVarName(sum_vars);
					for i=1:nivoaBM
						s_D(((i-1)*_ncase+1):i*_ncase,ir) = Data(:,a(i));
					end
				end
			else
				for ir = 1:sumv_1
					if a(ir)
						b = Data(:,sum_vars(ir));
						b1 = b==NaN;
						indx = find(b1);
						if any(indx)
							b(indx)=0;
						end
						s_D(:,ir)=rsum(b);

						indx = find(rprod(b1));
						if any(indx)
							s_D(indx,ir) = NaN;
						end
					else
						s_D(:,ir)=Data(:,sum_vars(ir));
					end
				end
			end
		else
			s_D = Data(:,sum_vars);
		end

		if def(missing)
			SetMissVal(varss, missT, missS)
		end
	else
		s_D = rD;
	end
	%ovde se razuna wgh za average, i 'whole' opcija bi trebala ovde da uzme
	%s_DSumm = ones(size(s_D,1),1);
	s_DSumm = rsum(s_D~=NaN)>0;
	if rowTotal
		rD = s_DSumm;
	end
	C = make(size(rL,1),size(cL,2),NaN);
	Res = make(size(C),0);
	Cn1 = make(size(C),0);
	Cn2 = make(size(C),0);
    
    SD = make(size(C),0);

	CBck = C;

	unwN = (rD(:,1)+rT(:,1))'*cD;
	unwNShr = s_DSumm'*cD;
	if wgh~=[]
		if nivoaBM>1
			indx = make(size(wgh,1)*nivoaBM,1,0);
			for i=1:nivoaBM
				indx((i-1)*size(wgh,1)+1:(i*size(wgh,1))) = wgh;
			end
			wgh = indx;
		end
		rD = rD.*wgh;
		rT = rT.*wgh;
		rTT = rTT.*wgh;
		s_DSumm = s_DSumm .*wgh;
		if sampleSize
			sampleSize = _ncase/sum(wgh);
		end
	end
	CD = rD'*cD;

	CA = rT'*cT;
	CA_= rTT'*cTT;

	CB = rT'*cD;
	CB_ = rTT'*cD;

	CC = rD'*cT;
	CC_= rD'*cTT;

	wghN = (rD(:,1)+rT(:,1))'*cD;
	wghNShr = s_DSumm'*cD;

	EC = ((CC+CD).*(CB+CD))./(CD+CA+CB+CC);

	if sampleSize
		RES = (CD - EC)./sqrt(EC) * sqrt(sampleSize);
		indx = find(EC*sampleSize<5);
	else
		RES = (CD - EC)./sqrt(EC);
		indx = find(EC<5);
	end

	%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% sve ocekivane frekvence manje od 5 nisu znacajne

	if any(indx)
		RES(indx) = 0;
	end

	if any(stat_type == 2)
		ROWPCT = CD./(CD+CC);
		
		PN = CD+CC;
		
		if sigc
			if sampleSize
				sSize = sampleSize;
			else
				sSize = 1;
			end
			
            SIGC = make(size(ROWPCT),"");
            for i3 = 1:size(PN,2)
				N = PN(:,i3) * sSize;
				n = CD(:,i3) * sSize;
				P = ROWPCT(:,i3);
				nc = size(P,1);
				
				ord = make(1,nc,(1:nc));
				tSIGC = make(size(ord),"");
				
				
				Pt = (make(nc,1, n')+make(1,nc, n))./(make(nc,1, N')+make(1,nc, N));
				S = sqrt(Pt.*(1-Pt).*(1./make(nc,1, N') + 1./make(1,nc, N)));
				Z = (make(nc,1, P') - make(1,nc, P)) ./ S;
				p = 2*normcdf(-abs(Z));
				
				indx = find(Z<0);
				if indx~=[]
					p(indx) = 1;
				end
				
				indx = find(P==0 | P==1);
				if indx~=[]
					p(indx,:) = 1;
					p(:,indx) = 1;
				end
				
				indx = find(PN(1,:) < 2);
				if indx~=[]
					p(indx,:) = 1;
					p(:,indx) = 1;
				end
				
				if nobonf
					indx = find(p< 0.1);
				else
					indx = find(p< 0.1/(nc*(nc-1)/2));
				end
				if indx~=[] & color3==1
					tSIGC(indx) = getAlpha(ord(indx),1)+"*";
				end

				if nobonf
					indx = find(p< 0.05);
				else
					indx = find(p< 0.05/(nc*(nc-1)/2));
				end
				if indx~=[]
					tSIGC(indx) = getAlpha(ord(indx),1);
				end

				if nobonf
					indx = find(p< 0.01);
				else
					indx = find(p< 0.01/(nc*(nc-1)/2));
				end
				if indx~=[]
					tSIGC(indx) = getAlpha(ord(indx),0);
				end
				
				
				
				for i2=1:nc
					SIGC(i2,i3) = implode(tSIGC(:,i2),"");
				end
			end
			%indx = find(SIGC==" ");
			%if indx~=[]
			%    SIGC(indx) = ".";
                    %end
        end
		
		
		
		ROWPCT = ROWPCT * 100;
	end
	if any(stat_type == 3) | any(stat_type == 13)
		COLPCT = CD./(CD+CB);
		PN = CD+CB;
		if sigc
			SIGC = make(size(COLPCT),"");
			if sampleSize
				sSize = sampleSize;
			else
				sSize = 1;
			end
			
			for i3 = 1:size(cGrn ,2)
				if cGrn(1,i3)~=cGrn(2,i3)
					rng = cGrn(1,i3):cGrn(2,i3);
					from = cGrn(1,i3)-1;
					for i1 = 1:size(PN,1)
						N = PN(i1,rng) * sSize;
						n = CD(i1,rng) * sSize;
						P = COLPCT(i1,rng);
						nc = size(P,2);
						
						ord = make(1,nc,(1:nc));
						tSIGC = make(size(ord),"");
						
						
						Pt = (make(nc,1, n)+make(1,nc, n'))./(make(nc,1, N)+make(1,nc, N'));
						S = sqrt(Pt.*(1-Pt).*(1./make(nc,1, N) + 1./make(1,nc, N')));
						Z = (make(nc,1, P) - make(1,nc, P')) ./ S;
						p = 2*normcdf(-abs(Z));
						
						indx = find(Z<0);
						if indx~=[]
							p(indx) = 1;
						end
						
						indx = find(P==0 | P==1);
						if indx~=[]
							p(indx,:) = 1;
							p(:,indx) = 1;
						end
						
						indx = find(PN(1,rng) < 2);
						if indx~=[]
								p(indx,:) = 1;
								p(:,indx) = 1;
						end
						if abc
							if nototal
								dodatak = 0;
							else
								dodatak = 0;
							end
						else
							if nototal
								dodatak = 1;
							else
								dodatak = 0;
							end
						end
						
						%***********************************************
						%		 Ovde je treci stepen za sigc, PI = 90%
						%
						if nobonf
							indx = find(p< 0.1);
						else
							indx = find(p< 0.1/(nc*(nc-1)/2));
						end
						if indx~=[] & color3==1
							tSIGC(indx) = getAlpha(ord(indx)+from-dodatak,1)+"*";
						end
						
						%************************************************

						if nobonf
							indx = find(p< 0.05);
						else
							indx = find(p< 0.05/(nc*(nc-1)/2));
						end
						if indx~=[]
							tSIGC(indx) = getAlpha(ord(indx)+from-dodatak,1);
						end

						if nobonf
							indx = find(p< 0.01);
						else
							indx = find(p< 0.01/(nc*(nc-1)/2));
						end
						
						if indx~=[]
							tSIGC(indx) = getAlpha(ord(indx)+from-dodatak,0);
						end
						for i2=1:nc
							SIGC(i1,from+i2) = implode(tSIGC(:,i2),"");
						end
					end
					%indx = find(SIGC==" ");
					%if indx~=[]
					%    SIGC(indx) = ".";
					%end
				end
			end
		end
		COLPCT = COLPCT*100;
	end
	if any(stat_type == 4)
		TOTPCT = CD./(CD+CB+CA+CC)*100;
	end
	if any(stat_type == 16) | any(stat_type == 17)
		TOTPCT_ = CD./(CD+CB_+CA_+CC_)*100;
	end
    
    
    
    %ovo je za JTI kalkulacija reziduala 
	%CB+CD
    if waveGr~=[] & sigw
        JtiCols = [];
        JtiRefs = [];
        for i=1:size(waveGr,2)
            x = waveGr(1,i):waveGr(1,i)+waveGr(2,i)*waveGr(3,i)-1:waveGr(2,i);
            x = x';
            JtiCols = [JtiCols, x];
            %JtiRefs = [JtiRefs, make(size(x), x(size(x,2)))];
            JtiRefs = [JtiRefs, x(1), x(1:size(x,2)-1)];
        end
        
        if any(stat_type == 3)
            JtiN = CD+CB;
            JtiP = ((COLPCT(:, JtiCols) .* JtiN(:, JtiCols)) + (COLPCT(:, JtiRefs) .* JtiN(:, JtiRefs))) ./ (JtiN(:, JtiCols) + JtiN(:, JtiRefs)) / 100;
            JtiSE = sqrt(JtiP.*(1-JtiP).* (1./JtiN(:, JtiCols) + 1./JtiN(:, JtiRefs)));
            RES(:,JtiCols) = (COLPCT(:, JtiCols) - COLPCT(:, JtiRefs))./JtiSE./100;
        end
	end
    %%%%%%%%%%%%%%%%%%%%%%%%%%%% JTI End
    
    
    BioMean = 0;
	irow = 1;
	icol = 1;
	if ~sumv_1
		mInd=0;
	end
	ir = 1;
	ic = 1;
	im = 0;
	%aa = [rD,cD,s_D];
	%aa(1:30,:)
	%[rT,cT,s_D]
	irow_ = size(C,1)/vars_1;
	icol_ = size(C,2)/col_1;
	for irow = 1:size(C,1)
		ir = ceil(irow/irow_);
		for icol = 1:size(C,2)
			ic = ceil(icol/icol_);

			if statsinrow
				it = tInd(irow);
				if sumstats_
					im = mInd(irow);
				end
			else
				it = tInd(icol);
				if sumstats_
					im = mInd(icol);
				end
			end

			if im==0
				if it == 1
					C(irow,icol) = CD(ir,ic);
					Res(irow,icol) = RES(ir,ic);
				elseif it == 2
					C(irow,icol) = ROWPCT(ir,ic);
					Res(irow,icol) = RES(ir,ic);
				elseif it == 3
					C(irow,icol) = COLPCT(ir,ic);
					Res(irow,icol) = RES(ir,ic);
                    

				elseif it == 4
					C(irow,icol) = TOTPCT(ir,ic);
					Res(irow,icol) = RES(ir,ic);
				elseif it == 16
					C(irow,icol) = TOTPCT_(ir,ic) * popSize;
					Res(irow,icol) = RES(ir,ic);
				elseif it == 17
					C(irow,icol) = TOTPCT_(ir,ic);
					Res(irow,icol) = RES(ir,ic);
				elseif it == 13
					p1 = COLPCT(ir,ic);
					p2 = COLPCT(ir,1);
					C(irow,icol) = p1/p2*100;
					if indexcolor
						if (C(irow,icol)>=150)
							Res(irow,icol) = 3;
						elseif (C(irow,icol)>=130)
							Res(irow,icol) = 2;
						elseif (C(irow,icol)>=115)
							Res(irow,icol) = 1.95;
						elseif (C(irow,icol)<=50)
							Res(irow,icol) = -3;
						elseif (C(irow,icol)<=70)
							Res(irow,icol) = -2;
						elseif (C(irow,icol)<=85)
							Res(irow,icol) = -1.95;
						end
					else
						Res(irow,icol) = (p1 - p2)/sqrt(p1*(100-p1)/CD(ir,ic) + p2*(100-p2)/CD(ir,1));
					end
				end

				if it == 14 | it == 15 | it == 14.1 | it == 15.1
					C(irow,icol) = CD(ir,ic);
					Cn1(irow,icol) = CD(ir,ic);
					Cn1Unw(irow,icol) = CD(ir,ic);
					if it == 15 | it == 15.1
						CBck(irow,icol) = CC(ir,ic);
						Cn2(irow,icol) = CC(ir,ic);
					else
						CBck(irow,icol) = CB(ir,ic);
						Cn2(irow,icol) = CB(ir,ic);
					end
				end
			else
				if it>4 & it~=13
					%ovde se razuna average, i 'whole' opcija bi trebala ovde da uzme sve
					%i treba odluciti sta sa sysmis u nulu ili u mean i kako onda SD

					%ovo neradi jer ima ~ ispred ???? PROVERITI OBAVEZNO - greska u KAL parseru
					%icol = (ic-1)*icm + (im-1)*imm_c + (~statsinrow)*(it-1) + 1;
					% zagrade mu smetaju oko ~statsinrow
					if def(valid_values) & max(size(validValues)) > 0
						Dnan = s_D(:,im)~=NaN & s_D(:,im) rin validValues; %ako su ocene odnda prosek od 1,2,3,4,5
					elseif marks
						Dnan = s_D(:,im)~=NaN & s_D(:,im)<=marks; %ako su ocene odnda prosek od 1,2,3,4,5
					else
						if sumv_1
							Dnan = s_D(:,im)~=NaN;
						else
							Dnan = s_D(:,ir)~=NaN;
						end
					end
					%Dnan = ones(size(s_D,1),1);
					col = cD(:,ic) & rD(:,ir) & Dnan;
					%printf("%d %d %d %d %d \n",ic,  sum(cD(:,ic)),  ir, sum(rD(:,ir)), sum(Dnan))
					%sum(col)
					svecol = cT(:,ic) & rD(:,ir) & Dnan;
					%printf("%d %d %d %d %d \n",ic,  sum(cT(:,ic)),  ir, sum(rD(:,ir)), sum(Dnan))
					%sum(svecol)
					%sverow = cD(:,ic) & rT(:,ir) & Dnan;
					if (any(col) & any(svecol))
						fcol = find(col);
						fsvecol = find(svecol);
						%fsverow = find(sverow);
						if sumv_1
							rb1 = s_D(fcol,im);
							rb2 = s_D(fsvecol,im);
						else
							rb1 = s_D(fcol,ir);
							rb2 = s_D(fsvecol,ir);
						end
						%rb3 = s_D(fsverow,im);
						if wgh~=[]
							if sampleSize
								pcol = wgh(fcol)*sampleSize;
								psvecol = wgh(fsvecol)*sampleSize;
								%psverow = wgh(fsverow)*sampleSize;
							else
								pcol = wgh(fcol);
								psvecol = wgh(fsvecol);
								%psverow = wgh(fsverow);
							end
							n1 = sum(pcol);
							n2 = sum(psvecol);
							%n3 = sum(psverow);
						else
							pcol = [];
							psvecol = [];
							psverow = [];
							n1 = max(size(fcol));
							n2 = max(size(fsvecol));
							%n3 = max(size(fsverow));
						end
						%print([irow,icol])
						Cn1(irow,icol) = n1;
						Cn2(irow,icol) = n2;
						Cn1Unw(irow,icol) = max(size(fcol));
						%%%%%%%%%%%%%%%%%%%%%%%%% sada statistik po statistik
						if it == 5 | it == 12
						%MEAN IC
                            BioMean = 1;
							as1=mean(rb1,pcol);
							as2=mean(rb2,psvecol);
							std1=sd(rb1,pcol);
							std2=sd(rb2,psvecol);
							dfh = round(n1+n2-2);
							s = (n1*std1^2 + n2*std2^2)/dfh;
							sedh = sqrt(s*(n1+n2)/n1/n2);
							dif=as1-as2;
							th=dif/sedh;
							if th==NaN | th==-Inf | th==Inf
								th=0;
							end
							th = norminv(tcdf(th,dfh));
							if it == 5
								C(irow,icol) = as1;
								Res(irow,icol) = th;
                                SD(irow,icol) = std1;
							else % CI
								if n1>1 & std1
									if statsinrow
										it = tInd(irow+1);
									else
										it = tInd(icol+1);
									end
									t95= -std1/sqrt(n1-1) * tinv((1+it/100)/2, n1-1);
									C(irow,icol) = as1 - t95;
									Res(irow,icol) = th;
									if statsinrow
										C(irow+1,icol) = as1 + t95;
										Res(irow+1,icol) = th;
									else
										C(irow,icol+1) = as1 + t95;
										Res(irow,icol+1) = th;
									end
								end
							end
						elseif it == 5.1 | it == 12.1
						%MEAN IC
							[as1,std1]=trim5(rb1,pcol);
							[as2,std2]=trim5(rb2,psvecol);
							dfh = round(n1+n2-2);
							s = (n1*std1^2 + n2*std2^2)/dfh;
							sedh = sqrt(s*(n1+n2)/n1/n2);
							dif=as1-as2;
							th=dif/sedh;
							if th==NaN | th==-Inf | th==Inf
								th=0;
							end
							th = norminv(tcdf(th,dfh));
							if it == 5.1
								C(irow,icol) = as1;
								Res(irow,icol) = th;
							else % CI
								if n1>1 & std1
									if statsinrow
										it = tInd(irow+1);
									else
										it = tInd(icol+1);
									end
									t95= -std1/sqrt(n1-1) * tinv((1+it/100)/2, n1-1);
									C(irow,icol) = as1 - t95;
									Res(irow,icol) = th;
									if statsinrow
										C(irow+1,icol) = as1 + t95;
										Res(irow+1,icol) = th;
									else
										C(irow,icol+1) = as1 + t95;
										Res(irow,icol+1) = th;
									end
								end
							end
						elseif it == 6
							mdn = percentile([rb1;rb2],50,[pcol;psvecol]);
							fcol = find(rb1<=mdn);
							fsvecol = find(rb2<=mdn);
							if wgh~=[]
								if any(fcol)
									n1ispod = sum(pcol(fcol));
								else
									n1ispod = 0;
								end
								if any(fsvecol)
									n2ispod = sum(psvecol(fsvecol));
								else
									n2ispod = 0;
								end
							else
								n1ispod = max(size(fcol));
								n2ispod = max(size(fsvecol));
							end

							est = n1*(n1ispod+n2ispod)/(n1+n2);
							Res(irow,icol) = -(n1ispod - est)./sqrt(est);
							mdn1=percentile(rb1,50,pcol);
							C(irow,icol) = mdn1;
						elseif it == 7
							C(irow,icol) = min(rb1);
						elseif it == 8
							C(irow,icol) = max(rb1);
						elseif it == 9
							C(irow,icol) = max(rb1) - min(rb1);
						elseif it == 10
							C(irow,icol) = sd(rb1,pcol)*sqrt(n1/(n1-1));
						elseif it == 11
							C(irow,icol) = var(rb1,pcol)*n1/(n1-1);
						elseif it == 14 | it == 15 %colshare %rowshare
							C(irow,icol) = sum(rb1,pcol);
							CBck(irow,icol) = sum(rb2,psvecol);
						elseif it == 14.1 | it == 15.1 %trimcolshare %trimrowshare
							C(irow,icol) = trim5(rb1,pcol) * n1; %0.9*n1;
							CBck(irow,icol) = trim5(rb2,psvecol)*n2; %0.9*n2;
						end
					else
						fcol = find(col);
						if fcol~=[]
							if sumv_1
								rb1 = s_D(fcol,im);
							else
								rb1 = s_D(fcol,ir);
							end
							if wgh~=[]
								if sampleSize
									pcol = wgh(fcol)*sampleSize;
								else
									pcol = wgh(fcol);
								end
								n1 = sum(pcol);
							else
								pcol = [];
								n1 = max(size(fcol));
							end
							%print([irow,icol])
							Cn1(irow,icol) = n1;
							Cn2(irow,icol) = 0;
							Cn1Unw(irow,icol) =  max(size(fcol));
							if it == 5 %MEAN
								as1=mean(rb1,pcol);
								C(irow,icol) = as1;
                                SD(irow,icol) = sd(rb1,pcol);
							elseif it == 5.1 %TRIMMEAN
								as1 = trim5(rb1,pcol);
								C(irow,icol) = as1;
							elseif it == 12
								as1=mean(rb1,pcol);
								std1 = sd(rb1,pcol);
								if n1>1 & std1
									if statsinrow
										it = tInd(irow+1);
									else
										it = tInd(icol+1);
									end
									t95= -std1 /sqrt(n1-1) * tinv((1+it/100)/2, n1-1);
									C(irow,icol) = as1 - t95;
									if statsinrow
										C(irow+1,icol) = as1 + t95;
									else
										C(irow,icol+1) = as1 + t95;
									end
								end
							elseif it == 12.1
								[as1,std1]=trim5(rb1,pcol);
								n1=0.9*n1;
								if n1>1 & std1
									if statsinrow
										it = tInd(irow+1);
									else
										it = tInd(icol+1);
									end
									t95= -std1 /sqrt(n1-1) * tinv((1+it/100)/2, n1-1);
									C(irow,icol) = as1 - t95;
									if statsinrow
										C(irow+1,icol) = as1 + t95;
									else
										C(irow,icol+1) = as1 + t95;
									end
								end
							elseif it == 6
								mdn1 = percentile(rb1,50,pcol);
								C(irow,icol) = mdn1;
								Res(irow,icol) = 0;
							elseif it == 7
								C(irow,icol) = min(rb1);
							elseif it == 8
								C(irow,icol) = max(rb1);
							elseif it == 9
								C(irow,icol) = max(rb1) - min(rb1);
							elseif it == 10
								C(irow,icol) = sd(rb1,pcol)*sqrt(n1/(n1-1));
							elseif it == 11
								C(irow,icol) = var(rb1,pcol)*n1/(n1-1);
							elseif it == 14 | it == 15 %rowshare i colshare
								C(irow,icol) = sum(rb1,pcol);
								CBck(irow,icol) = 0;
							elseif it == 14.1 | it == 15.1 %trimrowshare i trimcolshare
								C(irow,icol) = trim5(rb1,pcol) * n1; %0.9*n1;
								CBck(irow,icol) = 0;
							end
						else
							Cn1(irow,icol) = NaN;
							Cn2(irow,icol) = NaN;
							Cn1Unw(irow,icol) =  NaN;
							C(irow,icol) = NaN;
							CBck(irow,icol) = NaN;
						end
					end
				end
			end
		end
	end
    
    if waveGr~=[] & BioMean & sigw
        dfh = round(Cn1(:, JtiCols)+Cn1(:, JtiRefs)-2);
        s = (Cn1(:, JtiCols).*SD(:, JtiCols).^2 + Cn1(:, JtiRefs).*SD(:, JtiRefs).^2)./dfh;
        sedh = sqrt(s.*(Cn1(:, JtiCols)+Cn1(:, JtiRefs))./Cn1(:, JtiCols)./Cn1(:, JtiRefs));
        dif=C(:, JtiCols)-C(:, JtiRefs);
        th=dif./sedh;
        indx = find(th==NaN | th==-Inf | th==Inf);
        if indx~=[]
            th(indx)=0;
        end
        RES = norminv(tcdf(th',dfh'));
		indx = find(RES==NaN | RES==-Inf | RES==Inf);
		RES(indx) = 0;
		Res(:,2:Inf) = RES';
    end
    
    if sigc & size(SIGC,1)==0
        global Alph;
        smAlph = lcase(Alph);
        SIGC = make(size(C),"");
        for i3 = 1:size(cGrn ,2)
            if cGrn(1,i3)~=cGrn(2,i3)
                rng = cGrn(1,i3):cGrn(2,i3);
                from = cGrn(1,i3)-1;
                for i1 = 1:size(C,1)
                    nc = max(size(rng));
                    
                    N = make(1,nc, Cn1(i1,rng)');
                    AS = make(1,nc,C(i1,rng)');
                    STD = make(1,nc, SD(i1,rng)');
                    
                    ord = make(1,nc,(1:nc));
                    tSIGC = make(size(ord),"");
                    
                    
                    dfh = round(N'+N-2);
                    s = (N.*STD.^2 + N'.*(STD').^2)./dfh;
                    sedh = sqrt(s.*(N+N')./N./(N'));
                    dif=AS'-AS;
                    th=dif./sedh;
                    indx = find(th==NaN | th==-Inf | th==Inf);
                    if indx~=[]
                        th(indx)=0;
                    end
                    p = tcdf(-abs(th),dfh)*2;
                    
                    indx = find(th<0);
                    if indx~=[]
                        p(indx) = 1;
                    end
                    
                    %indx = find(P==0 | P==1);
                    %if indx~=[]
                    %    p(indx,:) = 1;
                    %    p(:,indx) = 1;
                    %end
                    
                    %indx = find(PN(1,rng) < 2);
                    %if indx~=[]
                    %    p(indx,:) = 1;
                    %    p(:,indx) = 1;
                    %end
                    if color3
                        if nobonf
                            indx = find(p< 0.10);
                        else
                            indx = find(p< 0.10/(nc*(nc-1)/2));
                        end
                        if indx~=[]
                            tSIGC(indx) = "*"+getAlpha(ord(indx)+from-nototal, 1);
                        end
                    end
                    
                    %***********************************************
                    %		 Ovde je treci stepen za sigc, PI = 90%
                    %
                    if nobonf
                        indx = find(p< 0.1);
                    else
                        indx = find(p< 0.1/(nc*(nc-1)/2));
                    end
                    if indx~=[] & color3==1
                        tSIGC(indx) = getAlpha(ord(indx)+from-nototal,1)+"*";
                    end
                        
                    
                    if nobonf
                        indx = find(p< 0.05);
                    else
                        indx = find(p< 0.05/(nc*(nc-1)/2));
                    end
                    if indx~=[]
                        tSIGC(indx) = getAlpha(ord(indx)+from-nototal, 1);
                    end

                    if nobonf
                        indx = find(p< 0.01);
                    else
                        indx = find(p< 0.01/(nc*(nc-1)/2));
                    end
                    if indx~=[]
                        tSIGC(indx) = getAlpha(ord(indx)+from-nototal, 0);
                    end
                    for i2=1:nc
                        SIGC(i1,from+i2) = implode(tSIGC(:,i2),"");
                    end
                end
            end
        end
    end
	%a sada jos jednom da bi napravili sharove

	if any(stat_type == 14 | stat_type == 15 | stat_type == 14.1 | stat_type == 15.1)
		if Cn1~=[]
			if mmin(Cn1==Cn1(1,:))
				wghN = Cn1(1,:);
				unwN = Cn1Unw(1,:);
			else
				wghN = make(size(wghN),NaN);
				unwN = wghN;
			end
		end

		indx = find(C==NaN);
		CnoNan = C;
		if any(indx)
			CnoNan(indx)=0;
		end

		indx = find(CBck==NaN);
		if any(indx)
			CBck(indx)=0;
		end

		indx = find(tInd==14);
		if any(indx)
			for i = 1:size(rGrnM,2)
				if statsinrow
					rgrnmt = indx(find(indx>=rGrnM(1,i) & indx<=rGrnM(2,i)));
					cgrnmt = 1:size(C,2);
				else
					rgrnmt = rGrnM(1,i):rGrnM(2,i);
					cgrnmt = indx;
				end

				Crowshr = C(rgrnmt,cgrnmt);
				CSumCshr = csum(CnoNan(rgrnmt,cgrnmt));
				C(rgrnmt,cgrnmt) = Crowshr ./ CSumCshr *100;

				Crowshr = CBck(rgrnmt,cgrnmt);
				CSumCshr = csum(CBck(rgrnmt,cgrnmt));
				CBck(rgrnmt,cgrnmt) = Crowshr ./ CSumCshr *100;

				a = C(rgrnmt,cgrnmt)/100;
				b = CBck(rgrnmt,cgrnmt)/100;
				Res(rgrnmt,cgrnmt) = (a - b) ./ sqrt(a.*(1-a)./Cn1(rgrnmt,cgrnmt) + b.*(1-b)./Cn2(rgrnmt,cgrnmt));
			end
		end

		indx = find(tInd==14.1);
		if any(indx)
			for i = 1:size(rGrnM,2)
				if statsinrow
					rgrnmt = indx(find(indx>=rGrnM(1,i) & indx<=rGrnM(2,i)));
					cgrnmt = 1:size(C,2);
				else
					rgrnmt = rGrnM(1,i):rGrnM(2,i);
					cgrnmt = indx;
				end

				Crowshr = C(rgrnmt,cgrnmt);
				CSumCshr = csum(CnoNan(rgrnmt,cgrnmt));
				C(rgrnmt,cgrnmt) = Crowshr ./ CSumCshr *100;

				Crowshr = CBck(rgrnmt,cgrnmt);
				CSumCshr = csum(CBck(rgrnmt,cgrnmt));
				CBck(rgrnmt,cgrnmt) = Crowshr ./ CSumCshr *100;

				a = C(rgrnmt,cgrnmt)/100;
				b = CBck(rgrnmt,cgrnmt)/100;
				Res(rgrnmt,cgrnmt) = (a - b) ./ sqrt(a.*(1-a)./Cn1(rgrnmt,cgrnmt) + b.*(1-b)./Cn2(rgrnmt,cgrnmt));
			end
		end

		indx = find(tInd==15);
		if any(indx)
			for i = 1:size(cGrnM,2)
				if statsinrow
					rgrnmt = indx;
					cgrnmt = cGrnM(1,i):cGrnM(2,i);
				else
					rgrnmt = 1:size(C,1);
					cgrnmt = indx(find(indx>=cGrnM(1,i) & indx<=cGrnM(2,i)));
				end

				Crowshr = C(rgrnmt,cgrnmt);
				CSumCshr = rsum(CnoNan(rgrnmt,cgrnmt));
				C(rgrnmt,cgrnmt) = Crowshr ./ CSumCshr *100;

				Crowshr = CBck(rgrnmt,cgrnmt);
				CSumCshr = rsum(CBck(rgrnmt,cgrnmt));
				CBck(rgrnmt,cgrnmt) = Crowshr ./ CSumCshr *100;

				a = C(rgrnmt,cgrnmt)/100;
				b = CBck(rgrnmt,cgrnmt)/100;
				Res(rgrnmt,cgrnmt) = (a - b) ./ sqrt(a.*(1-a)./Cn1(rgrnmt,cgrnmt) + b.*(1-b)./Cn2(rgrnmt,cgrnmt));
			end
		end

		indx = find(tInd==15.1);
		if any(indx)
			for i = 1:size(cGrnM,2)
				if statsinrow
					rgrnmt = indx;
					cgrnmt = cGrnM(1,i):cGrnM(2,i);
				else
					rgrnmt = 1:size(C,1);
					cgrnmt = indx(find(indx>=cGrnM(1,i) & indx<=cGrnM(2,i)));
				end

				Crowshr = C(rgrnmt,cgrnmt);
				CSumCshr = rsum(CnoNan(rgrnmt,cgrnmt));
				C(rgrnmt,cgrnmt) = Crowshr ./ CSumCshr *100;

				Crowshr = CBck(rgrnmt,cgrnmt);
				CSumCshr = rsum(CBck(rgrnmt,cgrnmt));
				CBck(rgrnmt,cgrnmt) = Crowshr ./ CSumCshr *100;

				a = C(rgrnmt,cgrnmt)/100;
				b = CBck(rgrnmt,cgrnmt)/100;
				Res(rgrnmt,cgrnmt) = (a - b) ./ sqrt(a.*(1-a)./Cn1(rgrnmt,cgrnmt) + b.*(1-b)./Cn2(rgrnmt,cgrnmt));
			end
		end
	end

	if index_num ~= NaN
		a = cGrn(2,1);
		Res(:,1:a)=[];
		C(:,1:a)=[];
		cL(:,1:a)=[];
		cGrn(:,1)=[];
		cGrn = cGrn-a;
		cGrnM(:,1)=[];
		cGrnM = cGrnM-a;
	end

	a=find(~unwN);
    if outzero & a~=[]
    	C(:,a)=[];
    	Res(:,a)=[];
    	allout = 0;
	    for i=1:size(cGrn,2)
	    	alloutl = allout;
	    	for j=cGrn(1,i):cGrn(2,i)
	    		if ~unwN(j)
	    			alloutl=alloutl+1;
	    			if outputType
	    				a=find(cL(:,j)~="~");
	    			else
	    				a=find(cL(:,j)~="");
	    			end
	    			if a~=[] %ima neki koji nije ~ pa zato to treba prebaciti u sledecu kolonu
	    				if j<size(unwN,2)
	    					if outputType
	    						a1 = find(cL(a,j+1)=="~");
	    					else
	    						a1 = find(cL(a,j+1)=="");
	    					end
	    					if a1~=[]
	    						cL(a(a1),j+1) = cL(a(a1),j);
	    					end
	    				end
	    			end
	    		end
	    	end
	    	cGrn(1,i)=cGrn(1,i)-allout;
    		cGrn(2,i)=cGrn(2,i)-alloutl;
    		allout=alloutl;
	    end
	    a=find(cGrn(1,:)>cGrn(2,:));
	    if any(a)
	    	cGrn(:,a)=[];
	    end
	    a=find(~unwN);
	    allout = 0;
	    for i=1:size(cGrnM,2)
	    	alloutl = allout;
	    	for j=cGrnM(1,i):cGrnM(2,i)
	    		if ~unwN(j)
	    			alloutl=alloutl+1;
	    		end
	    	end
	    	cGrnM(1,i)=cGrnM(1,i)-allout;
    		cGrnM(2,i)=cGrnM(2,i)-alloutl;
    		allout=alloutl;
	    end
	    unwN(:,a)=[];
	    wghN(:,a)=[];
	    cL(:,a)=[];
	    cNum(a)=[];
	    a=find(cGrnM(1,:)>cGrnM(2,:));
	    if any(a)
	    	cGrnM(:,a)=[];
	    end
	end

	if printTab

		global g_decimal;
		defdecimal = 0;
		if def(decimal)
			if max(size(decimal))
				decimal = str2num(decimal);
				if max(size(decimal))>3
					error("Argument 'decimal' can have maximum of 3 parameters")
				elseif max(size(decimal))==1
					decimal = make(1,2,decimal);
				end
				defdecimal = 1;
			else
				decimal = [1 0];
			end
		elseif g_decimal~=[]
			decimal = g_decimal;
		else
			decimal = [1 1];
		end
		
		boje;
		color =3;
		BojaFC = make(size(C),NaN);      %Boja Fonta za tabelu
		BojaBC = make(size(C),NaN);      %Boja Fonta za celiju
		if color>0
			indx = find(Res>2.58);
			if indx~=[];BojaBC(indx) = plava(1,1);BojaFC(indx) = plava(1,2);end;
			if color>1
				indx = find(Res<=2.58 & Res>1.96);
				if indx~=[];BojaBC(indx) = plava(2,1);BojaFC(indx) = plava(2,2);end;
				if color>2
					indx = find(Res<=1.96 & Res>1.64);
					if indx~=[];BojaBC(indx) = plava(3,1);BojaFC(indx) = plava(3,2);end;
				end
			end
		end

		if color>0
			indx = find(Res<-2.58);
			if indx~=[];BojaBC(indx) = crvena(1,1);BojaFC(indx) = crvena(1,2);end;
			if color>1
				indx = find(Res>=-2.58 & Res<-1.96);
				if indx~=[];BojaBC(indx) = crvena(2,1);BojaFC(indx) = crvena(2,2);end;
				if color>2
					indx = find(Res>=-1.96 & Res<-1.64);
					if indx~=[];BojaBC(indx) = crvena(3,1);BojaFC(indx) = crvena(3,2);end;
				end
			end
		end

		CStr = num2str("%0.1f",C);
		if 1
			strfor = "%0." + num2str("%0.0f", decimal(2)) + "f";
			CStr(:,(2- nototal):size(C,2)) = num2str(strfor,C(:,(2- nototal):NaN));
			m1 = 10^decimal(2);
			C(:,(2- nototal):NaN) = round(C(:,(2- nototal):NaN)*m1)/m1;
		end
		if ~nototal
			strfor = "%0." + num2str("%0.0f", decimal(1)) + "f";
			CStr(:,1) = num2str(strfor,C(:,1));
			m1 = 10^decimal(1);
			C(:,1) = round(C(:,1)*m1)/m1;
		end

		%%%%%%%% OVO SAMO ZA KWS
		%if size(cL,1)==2
		%	cL = [cL(1,:);make(1,size(cL,2),"");cL(2,:)];
		%end

		if statsinrow
			if sumv_1
				indx = find(rL(:,size(rL,2)-1)=="Index" | rL(:,size(rL,2)-1)=="Count");
			else
				indx = find(rL(:,size(rL,2))=="Index" | rL(:,size(rL,2))=="Count");
			end
			if any(indx)
				CStr(indx,:)=num2str("%0.0f",C(indx,:));
			end
		else
			if sumv_1
				indx = find(cL(size(cL,1)-1,:)=="Index" | cL(size(cL,1)-1,:)=="Count");
			else
				indx = find(cL(size(cL,1),:)=="Index" | cL(size(cL,1),:)=="Count");
			end
			if any(indx)
				CStr(:,indx)=num2str("%0.0f",C(:,indx));
			end
		end

		indx = find(C==0 | C==NaN | C==Inf | C==-Inf);
		if any(indx)
			CStr(indx)="";
		end


		LinCol = make(1, max(max(cGrn)), NaN); %pravljenje linija
		LinCol(cGrn(2,:)) = 2;                 %za kolone
		LinRow = make(max(max(rGrn)), 1, NaN);
		LinRow(rGrn(2,:)) = 2;

		LinCol = [make(1,size(rL,2),NaN), 2, LinCol];
		LinRow = [make(size(cL,1),1,NaN); 2; LinRow];
		LinCol(1)=2;
		LinRow(1)=2;
		af = make(size(cL),NaN);
		bf = make(size(rL),NaN);
		cf = make(size(cL,1),size(rL,2),NaN);
		ab = make(size(cL),NaN);
		bb = make(size(rL),NaN);
		cb = make(size(cL,1),size(rL,2),NaN);

		af(1,:) = LFColor(3);
		ab(1,:) = LBColor(3);

		%bf(:,1) = LFColor(3);
		%bb(:,1) = LBColor(3);

		%cf(:,1) = LFColor(3);
		%cb(:,1) = LBColor(3);
		cf(1,:) = LFColor(3);
		cb(1,:) = LBColor(3);

		labVar = "";
		if def(label) & label(1)~=""
			for i=1:max(size(label))
				labVar = labVar+" "+label(i);
			end
		end
		heading(3,labVar)
		print(CStr,rL,cL,[cf,af;bf,BojaFC],[cb,ab;bb,BojaBC], LinRow, LinCol);

		global g_footer;
		if def(footer)
			foot = footer(1);
			for i=2:max(size(footer))
				foot = foot+" "+footer(i);
			end
		elseif size(g_footer,1)
			foot = g_footer;
		else
			foot = "";
		end
		if foot~=""
			echo(foot);
		end
		print("")
	end
	clear global Selector
	if wgh~=[]
		Data(:,wgh_name) = wghOriginal;
	end
%C
%Res
%rL
%cL
%cGrn
%rGrn
%rNum
%cNum
%wghN
%unwN
%SIGC
%waveGr	
	
end

function xM = proveraBM(row, rowBM, rowL, col, colBM, colL)
%funkija vraca broj varijabli u bloku
	x = (make(size(rowBM),colBM)==rowBM).*colBM;
	%[mmax(x), mmin(x)==mmax(x)]
	xM = mmax(x);
	if xM
		indx = find(x);

		xm = mmin(x(indx));
		if (xm~=xM)
			if any(colBM==xM) & any(colBM~=xM)
				indx1 = find(colBM==xM);
				indx2 = find(colBM~=xM);
				rez1="";
				rez2="";
				for i=1:size(indx1,1)
					rez1 = rez1 + " " + col(indx1(i));
				end
				for i=1:size(indx2,1)
					rez2 = rez2 + " " + col(indx2(i));
				end
				error("Variables in  "+ colL + " have different levels: " + rez1 + " and " + rez2)
			elseif any(rowBM==xM) & any(rowBM~=xM)
				indx1 = find(rowBM==xM);
				indx2 = find(rowBM~=xM);
				rez1="";
				rez2="";
				for i=1:size(indx1,1)
					rez1 = rez1 + " "  + row(indx1(i));
				end
				for i=1:size(indx2,1)
					rez2 = rez2 + " " + row(indx2(i));
				end
				error("Variables in "+ rowL + " have different levels: " + rez1 + " and " + rez2)
			end
		end

		indx = find(colBM==xM);
		if any(indx)
			xM = max(size(GetVarName(col(indx(1)))));
		else
			indx = find(rowBM==xM);
			xM = max(size(GetVarName(row(indx(1)))));
		end
	else
		xM = 1;
	end
end

function bigtabMiss(M', setInd, RFDK)
% M je Mx3 matrica sa vrednostima za missing
% setInd je 0/1 da li se vrednost setuje ili iskljucuje
% RFDK je vrednost
	indx = find(M(:,3)==RFDK);
	if any(indx)
		M(indx,3)=NaN;
	end
	indx = find(M(:,2)==RFDK);
	if any(indx)
		M(indx,2)=M(indx,3);
		M(indx,3)=NaN;
	end
	indx = find(M(:,1)==RFDK);
	if any(indx)
		M(indx,1:2)=M(indx,2:3);
		M(indx,3)=NaN;
	end
	if setInd
		indx = find(M(:,3)==NaN & M(:,1)~=NaN & M(:,2)~=NaN);
		if any(indx)
			M(indx,3)=RFDK;
		end
		indx = find(M(:,2)==NaN & M(:,1)~=NaN);
		if any(indx)
			M(indx,2)=RFDK;
		end
		indx = find(M(:,1)==NaN);
		if any(indx)
			M(indx,1)=RFDK;
		end
	end
end

macro bt_miss
	global RF, DK

	lok_mis = ucase(lok_mis);

	Indx = [];
	indx = find(lok_mis=="NONE");
	if any(indx)
		if max(size(lok_mis))>1
			error("If specified NONE has to be only value")
		end
		lNONE = 1;
		Indx = [Indx;indx];
	else
		lNONE = 0;
	end

	indx = find(lok_mis=="RF");
	if any(indx)
		lRF = 1;
		Indx = [Indx;indx];
	end
	indx = find(lok_mis=="DK");
	if any(indx)
		lDK = 1;
		Indx = [Indx;indx];
	end

	indx = find(lok_mis=="NORF");
	if any(indx)
		lRF = 0;
		Indx = [Indx;indx];
	end
	indx = find(lok_mis=="NODK");
	if any(indx)
		lDK = 0;
		Indx = [Indx;indx];
	end
	if any(Indx)
		lok_mis(Indx)=[];
	end
	if size(lok_mis,1)
		lok_mis = str2num(lok_mis);
		if (lok_mis anyin tot_lok_mis)
			for lm=1:max(size(lok_mis))
				lm_i = find(tot_lok_mis == lok_mis(lm));
				if (lm_i~=[])
					tot_lok_mis(lm_i) = [];
				end
			end
		end
		if (lRF~=NaN & lRF~=0)+(lDK~=NaN & lDK~=0)+max(size(tot_lok_mis))+max(size(lok_mis))>3
			error("More then 3 missing values are specified")
		end
		tot_lok_mis = [tot_lok_mis,lok_mis];
		clear lok_mis;
	end
end
macro ukljuciMiss
	if def(missing) | size(g_missing,1)
		varss = GetVarName(varss);
		[missT, missS] = GetMissVal(varss);
		misT = missT;
		misS = missS;
		if sum(misT<0)
			error("Variable has defined 'range' for missing values, which is not supported yet")
		end

		if lNONE==1
			misT = make(size(missT),NaN);
			misS = make(size(missS),NaN);
		end
		if lRF~=NaN
			bigtabMiss(misS, lRF, RF);
		end
		if lDK~=NaN
			bigtabMiss(misS, lDK, DK);
		end
		if size(tot_lok_mis,1)
			for i=1:max(size(tot_lok_mis))
				bigtabMiss(misS, 1, tot_lok_mis(i))
			end
		end
		misT = rsum(misS~=NaN);
		SetMissVal(varss, misT, misS)
	end
	
end

function y = getAlpha(indx, small) 
    global Alph2;
    maxI = max(indx);
    if maxI>size(Alph2,2)
        alph = Alph2;
        i = 1;
        while maxI > size(alph,2)
            %error (sprintf("Imas %d kolonu u tabeli, a samo %d obelezije kolona\notvori stable.s u text editoru i dodaj obelezija u Alph2 definiciju u 19.redu",size(cL,2),size(Alph2,2)))
            alph = [alph, Alph2+num2str("%d",i)];
            i = i+1;
        end
        y = alph(indx);
    else
        y = Alph2(indx);
    end
    if small
        y = lcase(y);
    end
end