function [BestChart,BestFitness,farmlayout,farmlayout_NA] = GA(pop_size,wf,iterations,t,NA_type,tn,wt,algorithmDir,GA_parameters)
%% Parameter
BestChart=zeros(iterations,1);
BestFitness=zeros(iterations,1);
farmlayout = zeros(iterations,wf.rows*wf.cols);
farmlayout_NA = zeros(iterations,wf.rows*wf.cols);
pop_power = zeros(pop_size,wf.cols * wf.rows);

D = wf.turbine_num; down = 1;up = wf.rows *wf.cols;
%% Initialization
% Population initialization
[pop,lu] = windfarm_init(pop_size, wf.turbine_num,wf);
% Population evaluation

[popuFitness,power_order,lp_power_accum]= wf_fitness(wf,pop);
[popuFitness, sorted_index] = sort(popuFitness,'descend');
pop = pop(sorted_index,:);

FBest = popuFitness(1);

tournamentProbability = GA_parameters.tournamentProbability;
crossoverProbability  = GA_parameters.crossoverProbability;
mutationProbability   = GA_parameters.mutationProbability;
tournamentSize        = GA_parameters.tournamentSize;
%% loop
for iter = 1:iterations
%     pop = windfarm_constraint(pop, wf.NA_loc, D,down,up);
    [popuFitness,power_order,lp_power_accum]= wf_fitness(wf,pop);
    [popuFitness, sorted_index] = sort(popuFitness,'descend');
    pop = pop(sorted_index,:);
    power_order = power_order(sorted_index,:);
    lp_power_accum = lp_power_accum(sorted_index,:);

    if iter == 1
        BestFitness(iter,:) = popuFitness(1);
        BestChart(iter,:) = popuFitness(1) / wf.power_total;
        [best_farmlayout,best_farmlayout_NA]  = gene_layout_by_indices_one(wf,pop(1,:));
        farmlayout(iter,:) = best_farmlayout;
        farmlayout_NA(iter,:) = best_farmlayout_NA;
    else
        if popuFitness(1) > FBest
            FBest = popuFitness(1);
            BestFitness(iter,:) = FBest;
            BestChart(iter,:) = FBest / wf.power_total;
            [best_farmlayout,best_farmlayout_NA]  = gene_layout_by_indices_one(wf,pop(1,:));
            farmlayout(iter,:) = best_farmlayout;
            farmlayout_NA(iter,:) = best_farmlayout_NA;
        else
            BestFitness(iter,:) = BestFitness(iter-1);
            BestChart(iter,:) = BestChart(iter-1,:);
            farmlayout(iter,:) = farmlayout(iter-1,:);
            farmlayout_NA(iter,:) = farmlayout_NA(iter-1,:);
        end

    end

    pop = repalce_worst(pop, pop_size, power_order,wf.NA_loc,pop_power);
    
    %% Crossover
    temporaryPopulation = pop;
    for i = 1:2:pop_size
        i1 = TournamentSelect(popuFitness,tournamentProbability,tournamentSize);
        i2 = TournamentSelect(popuFitness,tournamentProbability,tournamentSize);
        r = rand;
        if (r < crossoverProbability)
            individual1 = pop(i1,:);
            individual2 = pop(i2,:);
            newIndividualPair = Cross(individual1, individual2);
            temporaryPopulation(i,:) = newIndividualPair(1,:);
            temporaryPopulation(i+1,:) = newIndividualPair(2,:);
        else
            temporaryPopulation(i,:) = pop(i1,:);
            temporaryPopulation(i+1,:) = pop(i2,:);
        end
    end
%     count = 0;
%     list = [];
%     for i =1:pop_size
%         if size(unique(pop(i,:)),2) <size(pop(i,:),2)
%             count = count +1;
%             list = [list;i];
%         end
%     end
%     fprintf('pop error: %d\n',count)
%     disp(list)
% 
% 
%     count = 0;
%     list = [];
%     for i =1:pop_size
%         if size(unique(temporaryPopulation(i,:)),2) <size(temporaryPopulation(i,:),2)
%             count = count +1;
%             list = [list;i];
%         end
%     end
%     fprintf('cross temporaryPopulation error: %d\n',count)
%     disp(list)

    %% Mutation
    temporaryPopulation(1,:) = pop(1,:);
    for i = 2:pop_size
        tempIndividual = Mutate(temporaryPopulation(i,:),mutationProbability);
        temporaryPopulation(i,:) = tempIndividual;
    end

%     count = 0;
%     list = [];
%     for i =1:pop_size
%         if size(unique(temporaryPopulation(i,:)),2) <size(temporaryPopulation(i,:),2)
%             count = count +1;
%             list = [list;i];
%         end
%     end
%     fprintf('Mutate temporaryPopulation error: %d\n',count)
%     disp(list)
    pop = temporaryPopulation;

    fprintf('%s NA %d Turbine Num:%d Wind %s run: %d iteration: %d  eta %f fitness %f\n',algorithmDir,NA_type,tn,wt,t,iter,(FBest / wf.power_total),FBest)

    %     [n_parent,parents_ind,parent_layouts,parent_layouts_NA,parent_pop_indices] = slection(pop,pop_NA,pop_indices,pop_size,ElitPercent,random_rate);
    %
    %
    %     [pop,pop_NA,pop_indices] = crossover(wf,pop,pop_NA,pop_indices,pop_size,n_parent,parents_ind,parent_layouts,parent_pop_indices);
    %
    %     [pop,pop_NA,pop_indices] = mutation(pop,pop_NA, pop_indices, pop_size, MutatPercent,wf);

    %     fprintf('%s NA %d Turbine Num:%d Wind %s run: %d iteration: %d  eta %f fitness %f\n',algorithmDir,NA_type,tn,wt,t,iter,BestChart(iter,:),BestFitness(iter,:))
end
disp('**')

end


