function result_status = proc_paid_images(B4,B5,B8,B8A,B11,TCI,JSONString)

pkg load image
JSONSentinel = JSONString

sizeB4 = [10980 10980];
format long;
red = 255;
green = 0;
stepsize = 25;
mat_step =1;
c_map = zeros(20,3);
while green<255
	green = green+stepsize
	if green > 255
		green = 255;
	end
	r_color = red/255;
	g_color = green/255;
	c_map(mat_step,1) = r_color;
	c_map(mat_step,2) = g_color;
	c_map(mat_step,3) = 0;
	mat_step = mat_step+1;
end
while red >0
	red = red-stepsize;
	if red<0
		red = 0;
	end
	r_color = red/255;
	g_color = green/255;
	c_map(mat_step,1) = r_color;
	c_map(mat_step,2) = g_color;
	c_map(mat_step,3) = 0;
	mat_step = mat_step+1;
end
%c_map

%finding field bounding coordinates
TileMaxLat =JSONSentinel.TileMaxLat
TileMinLat =JSONSentinel.TileMinLat
TileMaxLong = JSONSentinel.TileMaxLong
TileMinLong =JSONSentinel.TileMinLong
FieldMaxLat = JSONSentinel.FieldMaxLat
FieldMinLat = JSONSentinel.FieldMinLat
FieldMaxLong =JSONSentinel.FieldMaxLong
FieldMinLong =JSONSentinel.FieldMinLong
UID = JSONSentinel.UID
LatsPerPixel = (TileMaxLat - TileMinLat)/sizeB4(1)
LongsPerPixel = (TileMaxLong - TileMinLong)/sizeB4(1)
    PixelsPerLat = abs(sizeB4(1)/(TileMaxLat - TileMinLat));
PixelsPerLong = abs(sizeB4(2)/(TileMaxLong - TileMinLong));
StartPixelLat = PixelsPerLat*(1-abs(FieldMaxLat-TileMinLat));
StartPixelLat = round(abs(StartPixelLat))
EndPixelLat = PixelsPerLat*(1-abs(TileMinLat- FieldMinLat));
EndPixelLat =  round(abs(EndPixelLat))
if StartPixelLat > EndPixelLat
        temp = EndPixelLat;
        EndPixelLat = StartPixelLat;
        StartPixelLat = temp;
end

if StartPixelLat > 300
        StartPixelLat = StartPixelLat - 300
else
        StartPixelLat = 1
end

if EndPixelLat <sizeB4(1)-300
        EndPixelLat = EndPixelLat + 300
else
        EndPixelLat = sizeB4(1)
end
StartPixelLong = PixelsPerLong*(abs(TileMinLong - FieldMinLong));
StartPixelLong = round(abs(StartPixelLong))
EndPixelLong = PixelsPerLong*abs(TileMinLong - FieldMaxLong);
EndPixelLong = round(abs(EndPixelLong))
if StartPixelLong > EndPixelLong
        temp = EndPixelLong;
        EndPixelLong = StartPixelLong;
        StartPixelLong = temp;
end
if StartPixelLong > 300
        StartPixelLong = StartPixelLong - 300
else
        StartPixelLong =1
end
if EndPixelLong < sizeB4(1) - 300
        EndPixelLong = EndPixelLong+300
else
        EndPixelLong = sizeB4(1)
end


LatLength = EndPixelLat - StartPixelLat;
LongLength = EndPixelLong - StartPixelLong;
LatNDVI = LatLength
LongNDVI = LongLength
if StartPixelLat > sizeB4(1) | EndPixelLat > sizeB4(1) | StartPixelLong > sizeB4(1) | EndPixelLong > sizeB4(1)
	result_status2 = 'failed due to unbound start end'
	result_status = strcat(UID,"failed0")
else

%clip larger images
	clip = [StartPixelLat EndPixelLat StartPixelLong EndPixelLong];
	clip_tci = clip;
	B4 = B4(clip(1):clip(2),clip(3):clip(4));
	B8 = B8(clip(1):clip(2),clip(3):clip(4));
	TCI = TCI(clip(1):clip(2),clip(3):clip(4),:);

%find field coordinates for smaller resolution images
sizeB4 = [5490 5490];
rr = 11
LatsPerPixel = (TileMaxLat - TileMinLat)/sizeB4(1)
LongsPerPixel = (TileMaxLong - TileMinLong)/sizeB4(1)
    PixelsPerLat = abs(sizeB4(1)/(TileMaxLat - TileMinLat));
PixelsPerLong = abs(sizeB4(2)/(TileMaxLong - TileMinLong));
StartPixelLat = PixelsPerLat*(1-abs(FieldMaxLat-TileMinLat));
StartPixelLat = round(abs(StartPixelLat))
EndPixelLat = PixelsPerLat*(1-abs(TileMinLat- FieldMinLat));
EndPixelLat =  round(abs(EndPixelLat))
if StartPixelLat > EndPixelLat
        temp = EndPixelLat;
        EndPixelLat = StartPixelLat;
        StartPixelLat = temp;
end
if StartPixelLat > 150
        StartPixelLat = StartPixelLat - 150
else
        StartPixelLat = 1
end
if EndPixelLat <sizeB4(1)-150
        EndPixelLat = EndPixelLat + 150
else
        EndPixelLat = sizeB4(1)
end
StartPixelLong = PixelsPerLong*(abs(TileMinLong - FieldMinLong));
StartPixelLong = round(abs(StartPixelLong))
EndPixelLong = PixelsPerLong*abs(TileMinLong - FieldMaxLong);
EndPixelLong = round(abs(EndPixelLong))
if StartPixelLong > EndPixelLong
        temp = EndPixelLong;
        EndPixelLong = StartPixelLong;
        StartPixelLong = temp;
end
if StartPixelLong > 150
        StartPixelLong = StartPixelLong - 150
else
        StartPixelLong = 1
end
if EndPixelLong < sizeB4(1) - 150
        EndPixelLong = EndPixelLong+150
else
        EndPixelLong = sizeB4(1)
end
LatLength = EndPixelLat - StartPixelLat
LongLength = EndPixelLong - StartPixelLong

%clip smaller images
clip = [StartPixelLat EndPixelLat StartPixelLong EndPixelLong];
	B8A = B8A(clip(1):clip(2),clip(3):clip(4),:);
	B11 = B11(clip(1):clip(2),clip(3):clip(4),:);
	B5 = B5(clip(1):clip(2),clip(3):clip(4),:);
	rr = 12
	B8NDRE = B8;
	B4F_name = strcat(UID,"_B4.png");
	%imwrite(B4,B4F_name);
	B8NDRE = imresize(B8NDRE,[abs(clip(2)-clip(1)+1) abs(clip(4)-clip(3)+1)]);
	B8NDRE = single(B8NDRE);
	B4 = single(B4);
	B8 = single(B8);
	B8A = single(B8A);
	B11 = single(B11);
	B5 = single(B5);
	%calculate NDVI, NDRE, NDWI
	NDVI = (B8 - B4)./(B8 + B4);
	NDRE = (B8NDRE-B5)./(B8NDRE+B5);
	NDWI = (B11 - B8A)./(B8A + B11);
	ndvi_array = NDVI;
	NDSI = -NDWI;
	sNDSI = size(NDSI);
	
	%check if there is too much cloud
	sN = size(NDVI);
	pixelcount = sN(1)*sN(2);
	negcount = 0;
	negcount_cloud = 0;
	cloud_count = 0;
	for i = 1:sN(1)
		for j = 1:sN(2)
			if NDVI(i,j) < 0
				negcount=negcount+1;
			end
			%if B4(i,j)>0.65
			%	cloud_count=cloud_count+1;
			%end
		end
	end
	for i = 1:sNDSI(1)
		for j = 1:sNDSI(2)
			if NDSI(i,j) > 0.7
				negcount_cloud=negcount_cloud+1;
			end
			%if B4(i,j)>0.65
			%	cloud_count=cloud_count+1;
			%end
		end
	end
	negcount
	negcount_cloud
	pixelcount
	
	if(negcount>pixelcount*0.2 | negcount_cloud>pixelcount*0.5)
		result_status = strcat(UID,"failed1")
		result_status2 = 'Too much cloud in image'
		return
	end
	%check if the image has large black portion
	SizeNewB4 = size(B4);
	IsNaNMat = isnan(NDVI);
	IsNaNCount =sum(sum(IsNaNMat))
	if(IsNaNCount > SizeNewB4(1)*SizeNewB4(2)*0.5)
		result_status = strcat(UID,"failed1")
		result_status2 = 'Too much black in image'
		return
	else	

	%create ETCI image
		my_limit=0.5;
low_limit=0.008;
up_limit=0.992;
%----------------------------------------------------------------------
img=TCI;
[m1 n1 r1]=size(img);
%----------------------------------------------------------------------
if r1==3
    my_limit2=0.04;my_limit3=-0.04;
    a=rgb2ntsc(img);
    mean_adjustment=my_limit2-mean(mean(a(:,:,2)));
    a(:,:,2)=a(:,:,2)+mean_adjustment*(0.596-a(:,:,2));
    mean_adjustment=my_limit3-mean(mean(a(:,:,3)));
    a(:,:,3)=a(:,:,3)+mean_adjustment*(0.523-a(:,:,3));
else
    a=double(img)./255;
end
%----------------------------------------------------------------------
mean_adjustment=my_limit-mean(mean(a(:,:,1)));
a(:,:,1)=a(:,:,1)+mean_adjustment*(1-a(:,:,1));
if r1==3
    a=ntsc2rgb(a);
end
%----------------------------------------------------------------------
img=a.*255;
%--------------------calculation of vmin and vmax----------------------
for k=1:r1
    arr=sort(reshape(img(:,:,k),m1*n1,1));
    v_min(k)=arr(ceil(low_limit*m1*n1));
    v_max(k)=arr(ceil(up_limit*m1*n1));
end
%----------------------------------------------------------------------
if r1==3

trans = [ 0.299,  0.596,  0.211; 0.587, -0.274, -0.523; 0.114, -0.322,  0.312 ];

  v_min = v_min*trans;
	v_max = v_max*trans;
end
%----------------------------------------------------------------------
if abs(v_max(1)-v_min(1))<1e-3
	result_status = strcat(UID,"failed1")
                result_status2 = 'Too much black in image vmin vmax'
	return
end
%v_min
img=(img-v_min(1))/(v_max(1)-v_min(1));

startx=25;
endx = 900-25;
starty=181;
endy=1200-188;
TCIFilename = strcat(UID,"_ETCI.png");
%save ETCI image
imwrite(uint8(img.*255),TCIFilename);
TCIFilenameL = strcat(UID,"_ETCI_large.png");
%save ETCI image
imwrite(uint8(img.*255),TCIFilenameL);
%save TCI image
TCIFilename2 = strcat(UID,"_TCI.png");
imwrite(TCI,TCIFilename2);
TCIFilename2L = strcat(UID,"_TCI_large.png");
imwrite(TCI,TCIFilename2L);
%create color map NDVI image
NDVI_new = NDVI;
%ndvi_arr_name = strcat(UID,"_ndvi_arr.png");
%imwrite(NDVI_new,ndvi_arr_name);
%NDVI_new = imcomplement(NDVI);
colormap(c_map)
imagesc(NDVI_new, [0 1]);
axis square tight		
axis off
		ndviname = strcat(UID,"_ndvi")
		saveas(1,ndviname,'png')
%read NDVI image, crop, resize and mask on MNDVI
ndviname=strcat(ndviname,".png");
aa = imread(ndviname);
size_aa = size(aa)
sizeNDVI = size(NDVI)
%min(min(aa))
%max(max(aa))
aa_new = aa(startx:endx, starty:endy,:);

imwrite(aa_new,ndviname);

%create colormap NDRE image
NDRE_new = (1/(max(max(NDRE))))*NDRE;
%NDRE_new = imcomplement(NDRE);
%ndre_arr_name = strcat(UID,"_ndre_arr.png");
%imwrite(NDRE_new,ndre_arr_name);
colormap(c_map)
imagesc(NDRE_new, [0 1]);
axis square tight		
axis off
		NDRE_name = strcat(UID,"_ndre")
		saveas(1,NDRE_name,'png')
%read NDRE image, crop, resize and mask on MNDRE
i_name = strcat(UID,"_ndre.png");
aa = imread(i_name);
aa_new = aa(startx:endx, starty:endy,:);

imwrite(aa_new,i_name);

NDWI_new = (1/max(max(NDWI)))*NDWI;
NDWI_new = imcomplement(NDWI);
%ndwi_arr_name = strcat(UID,"_ndwi_arr.png");
%imwrite(NDWI_new,ndwi_arr_name);
colormap(c_map)
imagesc(NDWI_new);
axis square tight		
axis off
		imagename =strcat(UID,"_ndwi");
		saveas(1,imagename,'png')

%read NDRE image, crop, resize and mask on MNDRE	
i_name =strcat(UID,'_ndwi.png');
aa = imread(i_name);
aa_new = aa(startx:endx, starty:endy,:);

imwrite(aa_new,i_name);

result_status = strcat(UID,"successful")
end
end
