function generateDstAndWBCPatches(alreadySegmented)
      getTrainingSet32X32FcsStk_BFcsd_Feat_RGB('Type1', true, alreadySegmented); %This function get the fcsRGB patches, BFczdPathes and Features on RGB and stores it
      getTrainingSet32X32FcsStk_BFcsd_Feat_RGB('Type2', true, alreadySegmented);
end

function [FnlFeatSet, CellId] = getMyFeatures(cellsImg, imUint8RGB, slideName)
        %Lets do the processing on Green channel for the time being
        avMsk           = fspecial('average', 3);
        FnlFeatSet      = single(zeros(1, 14*3)); %3 Channels
        for ch = 1:3
            imGreen         = imUint8RGB(:, :, ch);
            imGreenDble     = im2double(imGreen);
            [grdMg, ~]      = imgradient(imGreenDble);
            meanLclMsk      = imfilter(imGreenDble, avMsk);
            varnLclMsk      = imfilter(imGreenDble.^2, avMsk) - meanLclMsk.^2;
            [SzM, SzN]      = size(imGreen);
            props           = regionprops(cellsImg, {'PixelIdxList', 'Centroid'});
            [numObjs, ~]    = size(props);
            myBgndImgI      = uint8(200* ones(SzM, SzN));
            FeatSet         = single(zeros(numObjs, 14));
            CellId          = [];
            for i = 1:numObjs
                pxlIdxLst      = props(i).PixelIdxList;
                [SR, SC]        = ind2sub(size(imGreen), pxlIdxLst);
                tmpI            = myBgndImgI;
                tmpI(pxlIdxLst) = imGreen(pxlIdxLst);
                %Get the patch
                minSR           = min(SR); maxSR = max(SR);  
                minSC           = min(SC); maxSC = max(SC);
                imgPatch        = tmpI(minSR:maxSR, minSC:maxSC);
                mskPatch        = cellsImg(minSR:maxSR, minSC:maxSC);
                FeatSet(i, 1:4) = getGLCMFeatPrPatch(imgPatch, mskPatch);

                imGrnDblPchPxls = imGreenDble(pxlIdxLst);
                meanGlblPxls    = mean(imGrnDblPchPxls);
                varGlblPxls     = var(imGrnDblPchPxls);
                minPxls         = min(imGrnDblPchPxls);
                maxPxls         = max(imGrnDblPchPxls);
                minGrdMag       = min(grdMg(pxlIdxLst));
                maxGrdMag       = max(grdMg(pxlIdxLst));   
                minLclMean      = min(meanLclMsk(pxlIdxLst));
                maxLclMean      = max(meanLclMsk(pxlIdxLst));
                minLclVarn      = min(varnLclMsk(pxlIdxLst));
                maxLclVarn      = max(varnLclMsk(pxlIdxLst));
                FeatSet(i,5:14) = [meanGlblPxls varGlblPxls minPxls maxPxls ...
                                   minGrdMag maxGrdMag minLclMean maxLclMean ...
                                   minLclVarn maxLclVarn];
                cellCentr       = round(props(i).Centroid);
                CellId(i).name  = [slideName '_' num2str(cellCentr(2)) '_' num2str(cellCentr(1))]; 
            end
            FnlFeatSet(1, (ch-1)*14+1:ch*14) = FeatSet;
        end
end

function GLCMFeat = getGLCMFeatPrPatch(imUint8Gray, Msk)
    im = imUint8Gray;
     %im has to be Single band [0 255]
    [SzR, SzC, SzO] = size(im);
    MskX = [Msk(:, 2:end) Msk(:, end)];
    MskY = [Msk(2:end, :); Msk(end, :)];
    MskD = Msk & MskX & MskY; %MskD = bwmorph(MskD, 'erode', 3); 
    %Compute the GLCM for each band for the reg and return the props. 
    Lvl = 32; LvlDiv = 256/Lvl; 
    regIntst = MskD; 
    GLCMMat = zeros(Lvl, Lvl);
    GLCMFeat = zeros(SzO, 4);
    for band = 1:SzO
        imBand = double(im);
        imgL = ceil((imBand+1)/LvlDiv);
        imgL(~regIntst) = -1;
        for i = 1:SzR
            for j = 1:SzC-1
                valLeft = imgL(i, j);
                valRight = imgL(i, j+1);
                if (valLeft ~= -1 && valRight ~= -1)
                    GLCMMat(valLeft, valRight) = GLCMMat(valLeft, valRight)+1;
                end
            end
        end
        GFeat = graycoprops(GLCMMat);
        GLCMFeat(band, :) = [GFeat.Contrast GFeat.Correlation GFeat.Energy GFeat.Homogeneity];
     end 
end

function getTrainingSet32X32FcsStk_BFcsd_Feat_RGB(Folder, dustExclFrTrain, alreadySegmented)
%This function saves the 32x32 focus stack patches of 5th, 21st, and 37th
%images of the stck FczdRGB, stote the 21st alone as BFczdRGB, and the 
%features of the BFcsdRGB as BFcsdFeatRGB. 
%
    getWBCDataSet = true; getDstDataSet = true;
    if (dustExclFrTrain)
        svPath = 'CNNPatches\DntUseDst4NgPatc\';
    else
        svPath = 'CNNPatches\DoUseDst4NgPatc\';
    end
    if (alreadySegmented)
        load ('CNNPatches\Segment\alreadySegImg');
        load ('CNNPatches\Segment\alreadyWBCImg');
        load ('CNNPatches\Segment\segIndxNames');
    end
    btchSz = 12000;
    BasePath = ['E:\Gopakumar\GopakumarIISTDrive\Dataset\' ...
                '4th IIST Visit\Mal01\Images\GndTruth\Separate\NewEasy\'];
    fnames = dir([BasePath '\ImgsMinAcrsStckTwoCls\' Folder '\*B.jpg']);
    FczStkPath = 'E:\Gopakumar\GopakumarIISTDrive\Dataset\Malaria\Patch32By32\DataSet32By32Color\Manual_But_Sel_Training\AllFcsStack\';
    
    numfids = length(fnames);  M = 32; N = 32; R = 1; Rad = 0.75*M;
    dskMskPstv = strel('disk', Rad); dskMskDust = strel('disk', Rad);
    posCnt = 0; negCnt = 0; negChnk = 1; prevVidIndx = ''; 
    for K = 1:numfids 
        K
        AbsFNme = [BasePath '\ImgsMinAcrsStckTwoCls\'  Folder '\' fnames(K).name]; 
        %Get details needed to acces the Ground Truth file
        [~, FileName, ~] = fileparts(AbsFNme);
        for i = 1:length(FileName); if (FileName(i) == '_'); vid = i-1; break; end; end;
        for j = vid+2:length(FileName); if (FileName(j) == '_'); stck = j-1; break; end; end;
        vidIndx = FileName(1:vid);  stckIndx = FileName(vid+2:stck);
        slideName = [vidIndx '_' stckIndx];
        %Read the min RGB image
        img     = imread(AbsFNme); [SizeR, SizeC, ~] = size(img);
        if(isTherWBC([vidIndx '_' stckIndx]))
            if (~alreadySegmented)
                [segmentedImg, segWithBndry, WBCs] = getSegmentation(img);
            else
                IndxPosSlctd = getSegmentedIndx([vidIndx '_' stckIndx], segIndxNames);
                segmentedImg = alreadySegImg(:, :, IndxPosSlctd);
                WBCs = alreadyWBCImg(:, :, IndxPosSlctd);
                if (sum(WBCs(:) > 0))
                    thereIsWBC = true;               
                end
            end
        else
            WBCs = false(SizeR, SizeC);
            thereIsWBC = false;
        end
        imF1 = im2single(imread([FczStkPath vidIndx '_' stckIndx '_5.jpg']));
        imF  = imread([FczStkPath vidIndx '_' stckIndx '_21.jpg']); imF2 = im2single(imF);
        imF3 = im2single(imread([FczStkPath vidIndx '_' stckIndx '_37.jpg']));
        imG(:, :, 1:3) = imF1; imG(:, :, 4:6) = imF2; imG(:, :, 7:9) = imF3;  
%         imG = 0;
%         imG = img; %(:, :, 2);
%         imG = im2single(imG);
        [SzM, SzN, SzO] = size(imG);
        %Get the positive patches
        load ([BasePath 'GndTrth\' vidIndx '_' stckIndx]);
        if (~strcmp(prevVidIndx, vidIndx)); 
            prevVidIndx = vidIndx;
            DstMsk = getDstLocnsFor(vidIndx);
        end 
        
        load ([BasePath 'GndTrth\Dst_' vidIndx '_' stckIndx]);
        [PlocR, PlocC] = find (Msk);
        falseMsk = false(size(Msk));
        mskdPstv = imdilate(Msk, dskMskPstv); 
        mskdDust = imdilate(DstMsk, dskMskDust); 
        
        regmin  = getMySuspectedRegion(img);
        if (thereIsWBC && getWBCDataSet)
            WBCMsks = WBCs & regmin; 
            Cntrids = regionprops(WBCMsks, 'centroid');
            [numWBCPts, ~] = size(Cntrids);
            %Add any missed location
            newIntstRegMn = ~WBCs & regmin;
            mYR = imF(:, :, 1); mYG = imF(:, :, 2);  mYB = imF(:, :, 3); 
            mYR(newIntstRegMn) = 255; mYG(newIntstRegMn) = 0;  mYB(newIntstRegMn) = 0; 
            showIm(:, :, 1) = mYR; showIm(:, :, 2) = mYG;  showIm(:, :, 3) = mYB; 
            imshow(showIm); [NC, NR, Btn] = ginput(1);
            while(Btn ~= 3) %until right button is pressed
                numWBCPts = numWBCPts + 1;
                Cntrids(numWBCPts).Centroid = [NC NR];
                [NC, NR, Btn] = ginput(1);
            end
            
            for i = 1:numWBCPts
                CCntrids = round(Cntrids(i).Centroid); 
                %For the time being, if there is a patch of the required size
                %around the point, then only we are considering it.
                minRw = CCntrids(2) - M/2;     minCl = CCntrids(1) - N/2;
                maxRw = CCntrids(2) + M/2-1;   maxCl = CCntrids(1) + N/2-1;
                if (minRw > 0 && minCl > 0 && maxRw <= SzM && maxCl <= SzN) 
                    currPatch = imG(minRw:maxRw, minCl:maxCl, :);
                    posCnt = posCnt + 1;
                    WBCStackFcsRGB(1:M, 1:N, 1:SzO, posCnt) = currPatch;
                    WBCBestFczdRGB(1:M, 1:N, 1:3, posCnt) = imG(minRw:maxRw, minCl:maxCl, 4:6);
                    falseMsk (minRw:maxRw, minCl:maxCl) = true;
                    [FeatSet, ~]    = getMyFeatures(falseMsk, imF, slideName);
                    falseMsk (minRw:maxRw, minCl:maxCl) = false;
                    WBCBstFczdFeat(posCnt, :) = FeatSet;
                end
            end
            thereIsWBC = false;
        end
        
        
        if (getDstDataSet)
            Cntrids = regionprops(DstMsk, 'centroid');
            [numDstPts, ~] = size(Cntrids);

            for i = 1:numDstPts
                CCntrids = round(Cntrids(i).Centroid); 
                %For the time being, if there is a patch of the required size
                %around the point, then only we are considering it.
                minRw = CCntrids(2) - M/2;     minCl = CCntrids(1) - N/2;
                maxRw = CCntrids(2) + M/2-1;   maxCl = CCntrids(1) + N/2-1;
                if (minRw > 0 && minCl > 0 && maxRw <= SzM && maxCl <= SzN) 
                    currPatch = imG(minRw:maxRw, minCl:maxCl, :);
                    posCnt = posCnt + 1;
                    DustStackFcsRGB(1:M, 1:N, 1:SzO, posCnt) = currPatch;
                    DustBestFczdRGB(1:M, 1:N, 1:3, posCnt) = imG(minRw:maxRw, minCl:maxCl, 4:6);
                    falseMsk (minRw:maxRw, minCl:maxCl) = true;
                    [FeatSet, ~]    = getMyFeatures(falseMsk, imF, slideName);
                    falseMsk (minRw:maxRw, minCl:maxCl) = false;
                    DustBstFczdFeat(posCnt, :) = FeatSet;
                end
            end
        end
    end  
    if (getDstDataSet)
        save ([svPath 'FczRGBDsStk_' Folder], 'DustStackFcsRGB');
        save ([svPath 'BFczdRGBDsStk_' Folder], 'DustBestFczdRGB');
        save ([svPath 'BFczdFeatRGBDsStk_' Folder], 'DustBstFczdFeat');
    end
    if (getWBCDataSet)
        save ([svPath 'FczRGBWcStk_' Folder], 'WBCStackFcsRGB');
        save ([svPath 'BFczdRGBWcStk_' Folder], 'WBCBestFczdRGB');
        save ([svPath 'BFczdFeatRGBWcStk_' Folder], 'WBCBstFczdFeat');
    end
    
%     negtvStackFcsRGB = negtvStackFcsRGB(1:M, 1:N, 1:SzO, 1:negCnt); %+1:end) = []; 
%     negtvBestFczdRGB = negtvBestFczdRGB(1:M, 1:N, :, 1:negCnt);
%     negtvBstFczdFeat = negtvBstFczdFeat(1:negCnt, :);
%     
%     save ([svPath 'FczRGBNgStk_' Folder '_' num2str(negChnk)], 'negtvStackFcsRGB'); 
%     save ([svPath 'BFczdRGBNgStk_' Folder '_' num2str(negChnk)], 'negtvBestFczdRGB');
%     save ([svPath 'BFczdFeatRGBNgStk_' Folder '_' num2str(negChnk)], 'negtvBstFczdFeat');
end
function  IndxPosSlctd = getSegmentedIndx(indxName, segIndxNames)
    [~, numNames] = size(segIndxNames);
    for i = 1:numNames
        if (strcmp(indxName, segIndxNames(i).name))
            IndxPosSlctd = i;
            break;
        end
    end
end
function stat = isTherWBC(fileName)
    stat = false;
    indx = {'1_9', '2_19', '3_3', '3_4', '4_8', '5_1', 'A_6', 'A_14', 'BB_3', 'BB_13', ...
            'BB_14', 'C_15', 'C_16', 'D_13', 'E_2', 'EE_4', 'F_3', 'F_5', 'F_9', 'FF_6', ...
            'GG_4', 'G_16', 'G_17', 'H_5', 'H_16', 'I_1', 'J_18', 'KK_9', 'L_6', ...
            'O_8', 'P_12', 'P_13', 'S_16', 'T_12', 'V_8', 'V_9', 'W_8', 'W_18', 'X_11', ...
            'X_13', 'Y_2', 'Z_11'};
    [~, num] = size(indx);
    for i = 1:num
        if (strcmp(fileName, indx(i)))
            stat = true;
            break;
        end
    end
end
function cDstLoc = getDstLocnsFor(vidIndx)
    Path = 'E:\Gopakumar\GopakumarIISTDrive\Dataset\4th IIST Visit\Mal01\Images\GndTruth\Separate\NewEasy\DstLocnsPgm\';
    load ([Path 'DstLocByPGM']); cnt = 0;
    vidIndxs = DstLocByPGM.vidIndx;
    [~, numVidIndxs] = size(vidIndxs);
    for i = 1:numVidIndxs
        if (strcmp(vidIndx, vidIndxs(i).name))
            cDstLoc = DstLocByPGM.DstLocByPgm(:, :, i);
            break;
        end
    end
end
function susReg = getMySuspectedRegion(img)
    MnSz = 900;
    hsv = rgb2hsv(img);
    vlue = hsv(:, :, 3);
    stDsk = strel('disk', 7);
    openimg = imopen(vlue, stDsk);
    mask = imregionalmin(openimg);
    %First Filtering Exclude all the Bgnd
    susReg = (bwareaopen((vlue < graythresh(vlue)), MnSz)) & mask;
end