root/applications/dual/SIDP/transformace/sidp_transformace.m @ 1090

function H=sidp_transformace(sidp_parameters, rsss_parameters,system,apriori)
%transformace
eps=10^-3;
system.dim=2;
apriori.eta_tilda0=eps;
apriori.eta_tilda0_range=apriori.y0_range/(system.sigma+apriori.y0_range);
apriori.beta_tilda0=eps;
apriori.beta_tilda0_range=apriori.b0_range^2/(apriori.P0+apriori.b0_range^2);

%pomocne promenne a konstanty
step=[1 sidp_parameters.n_grid];

zero=zeros(sidp_parameters.horizont,2); %zero(1,:) nema smysl, ale takhle je to pekne
range=zeros(sidp_parameters.horizont,2);
eigvec=zeros(sidp_parameters.horizont,2,2);

[H,mi]=init_hyperstate(sidp_parameters,apriori);  %vytvori H a mi, nastavi H0 a mi*

for i=1:sidp_parameters.n_pass
    
    for j=1:sidp_parameters.n_iter
        [i j]
        search_region=sidp_parameters.gama^(j-1)*sidp_parameters.lambda^(i-1)*sidp_parameters.search_region_init;
        [H zero range eigvec mi]=update_hyperstate(H,zero, range, eigvec,mi,sidp_parameters.n_grid,1+(i+j>2));   %vygeneruje trajektorie, v zasazenem regionu rovnomerne rozmisti body a prekopiruje nalezena mi*
      %>1
        for k=sidp_parameters.horizont:-1:1    %pozor na meze
            level=sidp_parameters.horizont+1-k;          
            k
            
            for l=1:size(H,2)  
                
                %if mod(l,1000)==0
                 %   l
                %end
                
                candidates=H(k,l,2)/H(k,l,1)*generate_candidates(H(k,l,1)/H(k,l,2)*mi(k,l), search_region, sidp_parameters.num_of_candidates,sidp_parameters.generate_candidates_mode);
                %_H(k,l,2)/H(k,l,1)
                %compare candidates
                
                    %simple comparing
                    realization1=zeros(size(candidates,1),rsss_parameters.n0);
                    state(1,:)=H(k,l,:);
                    
                    for m=1:sidp_parameters.num_of_candidates
                        
                        best_control=candidates(m);
                        for n=1:rsss_parameters.n0
                            
                            %generate realization     
                            for o=1:level
                                
                                eta_0=state(1)/(1+eps-state(1));
                                beta_0=sqrt(state(2)/(1+eps-state(2)));
                                best_control=-best_control *eta_0/beta_0;
                                
                                s=randn;
                                pom=sqrt(1+best_control^2);
                                eta_1=abs(eta_0+beta_0*best_control+pom*s);                                
                                next_state(1)=eta_1/(1+eta_1);
                                
                                realization1(m,n)=realization1(m,n)+(eta_1)^2;

                                if o<level
                                    
                                    %kalman
                                    beta_1=pom*beta_0+(best_control*s);
                                    next_state(2)=beta_1^2/(1+beta_1^2);
                                    
                                    state=next_state;

                                    %find in H      
                                    stav=state-zero(k+o,:);
                                    %for p=1:system.dim
                                        %stav1(p)=(stav(1)*eigvec(k+o,1,p)+stav(2)*eigvec(k+o,2,p))*range(k+o,p)*(sidp_parameters.n_grid-1)+1;
                                    %end
                                    stav1=stav.*range(k+o,:)*(sidp_parameters.n_grid-1)+1;
                                    index=min(max(round(stav1),1),sidp_parameters.n_grid)-[0 1];
                                   % plot(H(k+o,:,1),H(k+o,:,2), 'r.' ); hold on; plot(state(1,1),state(1,2), 'b+' );
                                    index=sum(index.*step);
                                    
                                    %find_in_hyperstate5(state, H(k+o,:,:))
                                    
                                    %plot(H(k+o,index,1),H(k+o,index,2), 'g+' ); hold off
                                     %plot(H(k+o,index,1),H(k+o,index,2), 'g+' )
                                     
                      %  if imag(index)~=0
                       %     ;
                       % end    
                                    best_control=mi(k+o,index);
                                end

                            end

                        end
                    end
                    mean_values=mean(realization1,2);
                    [min_val min_index]=min(mean_values);
                    best_control=candidates(min_index);
                    mi(k,l)=best_control;
            end
        end     
       
        vypis(:,:)=H(:,:,1);
        save 'eta_tilda.txt' vypis -ASCII;
        vypis(:,:)=H(:,:,2);
        save 'beta_tilda.txt' vypis -ASCII;
        save 'mi.txt' mi -ASCII;        
    end
    
end

end