update_weights Subroutine

   subroutine update_weights(kohonen_map,current_values,ihit,jhit,khit,&
               maximum_radius,iteration) 
   !========================================================================================
!!    Subroutine to update the weights
   class(two_level_self_organizing_map) :: kohonen_map
!! A `two_level_self_organizing_map` object
   real(kind=wp),dimension(:,:),intent(inout) :: current_values
!! A real array
   integer,intent(inout) :: ihit,jhit,khit,iteration
!! Integer variables
   real(kind=wp),intent(inout) :: maximum_radius
!! A real variable with the maximum radius
   real(kind=wp),dimension(size(current_values,1),size(current_values,2)) :: prototype_values
   real(kind=wp),dimension(size(current_values,1),size(current_values,2)) :: winner_values,term1,term2
   integer :: nx,ny,nz,debug_option,ic,current_pos,ineigh,jneigh,kneigh,idbg
   real(kind=wp) :: time_factor,current_radius,alpha,sigma2,h_neighborhood,real_distance,term3
   real(kind=wp) :: distance_ratio,geometric_distance2,eps,current_distance,lambda
   !
   nx=kohonen_map%parameters(1)%number_nodes_nx;
   ny=kohonen_map%parameters(1)%number_nodes_ny;
   nz=kohonen_map%parameters(1)%number_nodes_nz;
   debug_option=kohonen_map%parameters(1)%debug_level;
   idbg=kohonen_map%parameters(1)%idbg;
   lambda=2.0_wp*(1.0_wp/maximum_radius);
   time_factor=1.0_wp-dble(iteration)/&
               dble(kohonen_map%parameters(1)%number_epochs*kohonen_map%parameters(1)%number_patterns);
   !current_radius = max(maximum_radius*real(1001-iteration)/1000.0 + 0.9999999999,4.0_wp);
   current_radius = max(maximum_radius*time_factor,4.0_wp);
   !alpha = max(kohonen_map%parameters%learning_rate*(1.0_wp-real(iteration)/1000.0),0.01_wp);
   alpha = max(kohonen_map%parameters(1)%learning_rate*time_factor,0.01_wp);
   sigma2=current_radius**2
   !
   do ic=1,size(kohonen_map%coordinates,1)
      current_pos=position2index(ihit,jhit,khit,nx,ny);
      current_distance=kohonen_map%cells_distances(current_pos,ic)
      if(current_distance .lt. current_radius) then
         geometric_distance2=current_distance**2;
         call index2position(ic,nx,ny,nz,ineigh,jneigh,kneigh);
                  !write(*,*) ic,ineigh,jneigh,kneigh,ihit,jhit,khit
         select case(trim(kohonen_map%parameters(1)%neighborhood_type))
            case('gaussian')
               h_neighborhood=alpha*exp(-0.5*geometric_distance2/sigma2);
            case('bubble')
               h_neighborhood=alpha;
            end select
            if(debug_option .gt. 0) then
               write(idbg,*) ihit,jhit,khit,ineigh,jneigh,kneigh
            endif
            select case(trim(kohonen_map%parameters(1)%som_type))
            case('normal_som')                      
               call kohonen_map%grid(ineigh,jneigh,kneigh)%get_prototype(prototype_values);
               prototype_values=prototype_values+h_neighborhood*(current_values-prototype_values);
               call kohonen_map%grid(ineigh,jneigh,kneigh)%set_prototype(prototype_values)
            case('visom')
               !write(*,*) trim(kohonen_map%parameters%som_type)
               call kohonen_map%grid(ineigh,jneigh,kneigh)%get_prototype(prototype_values);
               call kohonen_map%grid(ihit,jhit,khit)%get_prototype(winner_values);
               real_distance=sum((winner_values-prototype_values)**2);
               if( (ineigh .eq. ihit) .and. (jneigh .eq. jhit) .and. (kneigh .eq. khit) ) then                           
                  prototype_values=prototype_values+h_neighborhood*(current_values-prototype_values);
               else
                  distance_ratio=dsqrt(real_distance)/(dsqrt(geometric_distance2)*lambda);
                  term1=(current_values-winner_values);
                  term2=(winner_values-prototype_values);
                  eps=max(1.0_wp*time_factor,0.0_wp);
                  term3=1.0_wp;!((1.0_wp-eps)+eps)
                  prototype_values=prototype_values+h_neighborhood*(term1+term2*(distance_ratio-1.0_wp)*term3);
                  endif
                        !write(*,*) iteration,dsqrt(real_distance),dsqrt(geometric_distance2)*lambda,distance_ratio
                  call kohonen_map%grid(ineigh,jneigh,kneigh)%set_prototype(prototype_values); 
            end select
      endif
      enddo!ic
   
   end subroutine update_weights