function [NodiFibase]=NodiSelez_k1k2(W2,Fi,NodiFi,NodiFibase,Geo,k1,k2)

%Scelta gdl aggiuntivi rispetto a 'NodiFibase' [ntbase 2] con stima geometrica-modale
%autovettori 'Fi' [nFi mFi] associati a corrispondenza 'NodiFi' [nFi 2],
%autovalori 'W2'=d1,d2 [mFi 1], geometria 'Geo' [nGeo 4]
% K1=2, k2=1
%function [NodiFibase]=NodiSelez(W2,Fi,NodiFi,NodiFibase,Geo)

if nargin==7
    %[nk2 mk2] = size(k2);
    [nFi mFi] = size(Fi);
    [nt mt] = size(NodiFi);
    [ntbase mtbase] = size(NodiFibase);
    [nW2 mW2] = size(W2);
    [nGeo mGeo] = size(Geo);
    
    if nFi==nt & mFi==nW2 & mW2==1 & mt==2 & mGeo==4% & mtbase==2
        if ntbase<nW2
            if rem(nW2,3)==0 & rem(ntbase,3)==0
                h = waitbar(0,'Wait ... ','Name','MoGeSeC selection');
                %Completa nodi di riduzione
                %con nodi 3D, cioč aggiungo 3 modi alla volta corrispondenti a nodo fisico (x,y,z)
                for cont2 = (ntbase+7) : 3 : nW2
                    waitbar(cont2/nW2,h,['Wait ... MoGeSeC selection: node ',num2str(cont2),' on ',num2str(nW2),' nodes'])

                    %ValNodo = zeros(nGeo,1);
                    ValNodoG = zeros(nGeo,1);
                    ValNodoM = zeros(nGeo,1);
                    
                    %Contributo spaziale
                    if ntbase==0
                        ValNodoG = ones(nGeo,1);
                    else
                        %approccio min(energia potenziale gravitazionale sul dominio)
                        for cont1 = 1 : 3 : ntbase
                            %k1 = 2;         %peso esponente funzione geometrica (migliore)
                            %k1 = 0.7;         %peso esponente funzione geometrica (simile al caso k1 = 0.5;)
                            %k1 = 0.5;       %peso esponente funzione geometrica (simile al caso k1 = 1;)
                            %k1 = 0.25;       %peso esponente funzione geometrica (peggiore)
                            temp1 = find(NodiFibase(cont1,1)==Geo(:,1));
                            if cont1==1
                                ValNodoG = 1./ ...
                                    (((Geo(:,2)-Geo(temp1,2)).^2 + ...
                                    (Geo(:,3)-Geo(temp1,3)).^2 + ...
                                    (Geo(:,4)-Geo(temp1,4)).^2) ./ ...
                                    (max((Geo(:,2)-Geo(temp1,2)).^2 + ...
                                    (Geo(:,3)-Geo(temp1,3)).^2 + ...
                                    (Geo(:,4)-Geo(temp1,4)).^2))).^k1;
                            else
                                ValNodoG = ValNodoG + 1./ ...
                                    (((Geo(:,2)-Geo(temp1,2)).^2 + ...
                                    (Geo(:,3)-Geo(temp1,3)).^2 + ...
                                    (Geo(:,4)-Geo(temp1,4)).^2) ./ ...
                                    (max((Geo(:,2)-Geo(temp1,2)).^2 + ...
                                    (Geo(:,3)-Geo(temp1,3)).^2 + ...
                                    (Geo(:,4)-Geo(temp1,4)).^2))).^k1;
                            end
                        end
                        ValNodoG = 1./ValNodoG;
                        [temp,indice] = max(ValNodoG);
                        ValNodoG = (1/temp(1,1))*ValNodoG;
                    end

                    %Contributo modale
                    %k2 = k2;       %peso esponente funzione modale
                    for cont1 = 1 : 1 : (ntbase+3)
                        ValNodoM = ValNodoM + ...
                                   ((Fi(1:3:nFi,cont1).^2 + ...
                                     Fi(2:3:nFi,cont1).^2 + ...
                                     Fi(3:3:nFi,cont1).^2) ./ ...
                                   max(Fi(1:3:nFi,cont1).^2 + ...
                                       Fi(2:3:nFi,cont1).^2 + ...
                                       Fi(3:3:nFi,cont1).^2)).^k2;
                    end
                    ValNodoM = ValNodoM/(ntbase+3);
                    
                    %Contributo totale
                    ValNodo = ValNodoG.*ValNodoM;
                    %[temp,indice] = max(ValNodo);
                    %ValNodo = (1/temp(1,1))*ValNodo;
                    [temp,indice] = max(ValNodo);
                    
                    %Grafici
                    % figure
                    % plot(ValNodo,'-k')
                    % hold on
                    % plot(ValNodoG,'-r')
                    % plot(ValNodoM,'-b')
                    % plot(indice,temp,'or')

                    %Grafici
%                     figure
%                     size_Fem = ((max(Geo(:,2))-min(Geo(:,2)))^2+(max(Geo(:,3))-min(Geo(:,3)))^2+(max(Geo(:,4))-min(Geo(:,4)))^2)^0.5;
%                     dist_num = 0.01;
%                     colorplot = [1,1,1;0.971,0.996,1;0.943,0.993,1;0.914,0.989,1;0.885,0.986,1;0.857,0.982,1;0.828,0.978,1;0.799,0.975,1;0.771,0.971,1;0.742,0.968,1;0.713,0.964,1;0.685,0.961,1;0.656,0.957,1;0.627,0.953,1;0.599,0.95,1;0.57,0.946,1;0.541,0.943,1;0.513,0.939,1;0.484,0.935,1;0.455,0.932,1;0.427,0.928,1;0.398,0.925,1;0.369,0.921,1;0.341,0.918,1;0.312,0.914,1;0.283,0.91,1;0.254,0.907,1;0.226,0.903,1;0.197,0.9,1;0.168,0.896,1;0.14,0.892,1;0.111,0.889,1;0.103,0.821,1;0.094,0.752,1;0.085,0.684,1;0.077,0.615,1;0.068,0.547,1;0.06,0.479,1;0.051,0.41,1;0.043,0.342,1;0.034,0.274,1;0.026,0.205,1;0.017,0.137,1;0.009,0.068,1;0,0,1;0,0,0.962;0,0,0.923;0,0,0.885;0,0,0.847;0,0,0.808;0,0,0.77;0,0,0.732;0,0,0.694;0,0,0.655;0,0,0.617;0,0,0.579;0,0,0.54;0,0,0.502;0,0,0.46;0,0,0.418;0,0,0.376;0,0,0.335;0,0,0.293;0,0,0.251];
%                     subplot(1,3,3)
% %                     discr = 1;
% %                     %discr = 10^(round(log10(size(Geo,1))-3))/4;
% %                     if discr<1
% %                         discr = 1;
% %                     end
% %                     for cont1 = 1 : discr : size(Geo,1)
% %                         %plot3(Geo(cont1,2),Geo(cont1,3),Geo(cont1,4),'o','Linewidth',1,'Color',colorplot(round(ValNodo(cont1,1)/temp*(size(colorplot,1)-1))+1,:));
% %                         %text(Geo(cont1,2)-0.025,Geo(cont1,3)+0.025,Geo(cont1,4),num2str(r),'Color','b');
% %                         plot3(Geo(cont1,2),Geo(cont1,3),Geo(cont1,4),'o','MarkerSize',4,...
% %                             'Color',colorplot(round(ValNodo(cont1,1)/temp*(size(colorplot,1)-1))+1,:),...
% %                             'MarkerFaceColor',colorplot(round(ValNodo(cont1,1)/temp*(size(colorplot,1)-1))+1,:));
% %                         hold on
% %                     end
%                     for cont1 = 1 : size(colorplot,1)
%                         cont2 = find(ValNodo/temp>=(cont1-1)/size(colorplot,1));
%                         if cont1==size(colorplot,1)
%                             cont3 = find(ValNodo/temp<=cont1/size(colorplot,1));
%                         else
%                             cont3 = find(ValNodo/temp<cont1/size(colorplot,1));
%                         end
%                         cont4 = intersect(cont2,cont3,'rows');
% %                         plot3(Geo(cont4,2),Geo(cont4,3),Geo(cont4,4),'o','MarkerSize',4,...
% %                             'Color',colorplot(cont1,:),'MarkerFaceColor',colorplot(cont1,:));
%                         hold on
%                     end
%                     cont1 = [];
%                     for cont3 = 1:fix(size(NodiFibase,1)/3)
%                         cont2 = find(NodiFibase(3*(cont3-1)+1,1)==Geo(:,1));
%                         cont1(cont3,1) = cont2(1,1);
%                     end
%                     col = 'or';
%                     colsize = 6;
% %                     plot3(Geo(cont1,2),Geo(cont1,3),Geo(cont1,4),col,'Linewidth',4,'MarkerSize',colsize);
%                     axis equal
%                     xlabel('\itx \rmaxis [mm]')%,'FontSize',12,'FontName','TimesNewRoman')
%                     ylabel('\ity \rmaxis [mm]')%,'FontSize',12,'FontName','TimesNewRoman')
%                     zlabel('\itz \rmaxis [mm]')%,'FontSize',12,'FontName','TimesNewRoman')
%                     grid on
%                     title('\itW')
%                     view([-75 5])
%                     [az,el] = view;
%                     axislim = axis;
% %                     caxis([0 100])
% %                     colorbar('Ylim',[0 100])
% %                     colormap(colorplot)
%                     set(gcf,'color',[1 1 1])
% 
%                     subplot(1,3,1)
%                     for cont1 = 1 : size(colorplot,1)
%                         cont2 = find(ValNodoG/max(ValNodoG)>=(cont1-1)/size(colorplot,1));
%                         if cont1==size(colorplot,1)
%                             cont3 = find(ValNodoG/max(ValNodoG)<=cont1/size(colorplot,1));
%                         else
%                             cont3 = find(ValNodoG/max(ValNodoG)<cont1/size(colorplot,1));
%                         end
%                         cont4 = intersect(cont2,cont3,'rows');
% %                         plot3(Geo(cont4,2),Geo(cont4,3),Geo(cont4,4),'o','MarkerSize',4,...
% %                             'Color',colorplot(cont1,:),'MarkerFaceColor',colorplot(cont1,:));
%                         hold on
%                     end
%                     axis equal
%                     xlabel('\itx \rmaxis [mm]')%,'FontSize',12,'FontName','TimesNewRoman')
%                     ylabel('\ity \rmaxis [mm]')%,'FontSize',12,'FontName','TimesNewRoman')
%                     zlabel('\itz \rmaxis [mm]')%,'FontSize',12,'FontName','TimesNewRoman')
%                     grid on
%                     title('\itW_g')
%                     view([az,el])
%                     axis(axislim)
%                     set(gcf,'color',[1 1 1])
%                     
%                     subplot(1,3,2)
%                     for cont1 = 1 : size(colorplot,1)
%                         cont2 = find(ValNodoM/max(ValNodoM)>=(cont1-1)/size(colorplot,1));
%                         if cont1==size(colorplot,1)
%                             cont3 = find(ValNodoM/max(ValNodoM)<=cont1/size(colorplot,1));
%                         else
%                             cont3 = find(ValNodoM/max(ValNodoM)<cont1/size(colorplot,1));
%                         end
%                         cont4 = intersect(cont2,cont3,'rows');
% %                         plot3(Geo(cont4,2),Geo(cont4,3),Geo(cont4,4),'o','MarkerSize',4,...
% %                             'Color',colorplot(cont1,:),'MarkerFaceColor',colorplot(cont1,:));
%                         hold on
%                     end
%                     axis equal
%                     xlabel('\itx \rmaxis [mm]')%,'FontSize',12,'FontName','TimesNewRoman')
%                     ylabel('\ity \rmaxis [mm]')%,'FontSize',12,'FontName','TimesNewRoman')
%                     zlabel('\itz \rmaxis [mm]')%,'FontSize',12,'FontName','TimesNewRoman')
%                     grid on
%                     title('\itW_m')
%                     view([az,el])
%                     axis(axislim)
%                     set(gcf,'color',[1 1 1])

                    NodiFibase = [NodiFibase; ...
                                  Geo(indice(1,1),1) 1; ...
                                  Geo(indice(1,1),1) 2; ...
                                  Geo(indice(1,1),1) 3];
                    [temp,indice] = sort(NodiFibase(:,1)+0.1*NodiFibase(:,2));
                    NodiFibase = NodiFibase(indice,:);
                    ntbase = ntbase+3;
                end
                close(h);

            else
                error('Numero modi o numero nodi base non multipli di 3 ... non posso continuare ...')
            end
        elseif ntbase==nW2
            NodiFibase = NodiFibase;
        else
            error('Troppi nodi (n > r) ... non posso continuare ...')
        end
    end
    
else
    error('Dati input / output errati ... non posso continuare ...')
end
