train_cluster_layer Subroutine

   subroutine train_cluster_layer(kohonen_map)
   !========================================================================================
!! Subroutine to train the cluster layer of a two_level self_organized_map
      class(two_level_self_organizing_map) :: kohonen_map
!! A `two_level_self_organizing_map` object
      integer :: ix,iy,iz,iepoch,ihit,ic,number_variables,iteration,ineigh,current_pos,ipattern,ipos,pos
      type(kohonen_prototype) :: current_prototype
      real(kind=wp) :: dist_hit,dist,distortion,maximum_radius,minimum_radius,current_radius
      real(kind=wp) :: alpha,h_neighborhood,sigma2,geometric_distance2
      real(kind=wp),dimension(kohonen_map%number_variables1,&
                              kohonen_map%number_variables2) :: current_values,prototype_values
      real(kind=wp),dimension(kohonen_map%number_variables1*kohonen_map%number_variables2) :: centers
      real(kind=wp),dimension(kohonen_map%number_clusters) :: distance_units
      real(kind=wp),dimension(kohonen_map%number_variables1*kohonen_map%number_variables2,&
                              kohonen_map%number_clusters) :: centers1
      logical :: testop
   !
      number_variables=kohonen_map%number_variables1*kohonen_map%number_variables2;
      !maximum_radius=real(max(kohonen_map%parameters%number_nodes_nx,kohonen_map%parameters%number_nodes_ny));
      maximum_radius=real(kohonen_map%number_clusters);
      minimum_radius=1.0_wp;
   
   !
      iteration = 0
      distortion=0.0_wp;
      if(kohonen_map%parameters(1)%view_flag) then
      write(*,*) 'TWO LEVEL SOM - Cluster Layer: Training starting...'
      endif
   !
      do iepoch=1,kohonen_map%parameters(2)%number_epochs
         if(kohonen_map%parameters(1)%view_flag) then
            write(6,*) ' Starting epoch -- distortion',iepoch,' -- ',distortion
         endif
         distortion=0.0_wp;
         do iz=1,size(kohonen_map%grid,3);
            do iy=1,size(kohonen_map%grid,2);
               do ix=1,size(kohonen_map%grid,1);
               iteration = iteration + 1
               ihit = 0;
               dist_hit = 100000.0;
               current_prototype=kohonen_map%grid(ix,iy,iz);
               call current_prototype%get_prototype(current_values);
               do ic=1,size(kohonen_map%cluster_layer)
                  dist=kohonen_map%cluster_layer(ic)%distance(current_prototype,kohonen_map%distance_function);
                  dist = dist/float(number_variables);                 
                  if (dist .lt. dist_hit) then
                        dist_hit = dist
                        ihit = ic;
                  endif 
               enddo!ic
               distortion = distortion + dist_hit;
               current_radius = max(maximum_radius*real(1001-iteration)/1000.0 + 0.9999999999,0.49_wp);
               sigma2=current_radius**2;
               !current_radius=0.0_wp;
               ! define learning rate
               alpha = max(kohonen_map%parameters(2)%learning_rate*&
                        (1.0_wp-real(iteration)/1000.0),0.01_wp);
               do ineigh = ihit-int(current_radius), ihit+int(current_radius)
                  if(ineigh .ge. 1 .and. ineigh .le. size(kohonen_map%cluster_layer)) then
                     select case(trim(kohonen_map%parameters(2)%neighborhood_type))
                        case('gaussian');
                        geometric_distance2=(ihit-ineigh)**2;
                        h_neighborhood=alpha*exp(-0.5*geometric_distance2/sigma2);
                        case('bubble');
                        h_neighborhood=alpha;
                     end select
                     call kohonen_map%cluster_layer(ineigh)%get_prototype(prototype_values);
                     prototype_values=prototype_values+h_neighborhood*(current_values-prototype_values);
                     call kohonen_map%cluster_layer(ineigh)%set_prototype(prototype_values);
                  endif
               enddo!ineigh
   
               enddo!iz
            enddo!iy
         enddo!ix
      enddo!iepoch
      if(kohonen_map%parameters(1)%view_flag) then
      write(*,*) 'TWO LEVEL SOM - Cluster Layer: Training finished';
   !
      write(6,*) 
      write(6,*) 'TWO LEVEL SOM - Cluster Layer: Find Best Match Unit...'
      write(6,*)
      endif
   !
      current_pos=0;
      do iz = 1,size(kohonen_map%grid,3);
         do iy=1,size(kohonen_map%grid,2);
            do ix=1,size(kohonen_map%grid,1);
               current_pos=current_pos+1;
               ihit = 0;
               dist_hit = 100000.0;
               current_prototype=kohonen_map%grid(ix,iy,iz);
               call current_prototype%get_prototype(current_values);
               do ic = 1, size(kohonen_map%cluster_layer);
                  dist=kohonen_map%cluster_layer(ic)%distance(current_prototype,kohonen_map%distance_function);
                  dist = dist/float(number_variables);
                  distance_units(ic)=dist;
                  if (dist .lt. dist_hit) then
                     dist_hit = dist
                     ihit = ic;
                     !write(*,*) 'hit= ',ix,iy,iz,ihit
                  endif
               enddo
               kohonen_map%cluster_cells_index(current_pos,1)=ix;
               kohonen_map%cluster_cells_index(current_pos,2)=iy;
               kohonen_map%cluster_cells_index(current_pos,3)=iz;
               kohonen_map%cluster_cells_index(current_pos,4)=ihit;
               !write(6,*) 'hit= ',ix,iy,iz,ihit
   !            if(kohonen_map%number_patterns(ix,iy,iz) .ge. 1) then
                  kohonen_map%grid_cluster(ix,iy,iz)=ihit;
   !            endif
               
   !             if(debug_option .gt. 0) then
   !               write(idbg,*) ix,iy,iz,ihit
   !             endif
            enddo !iz
         enddo!iy
      enddo!ix
   ! Print grid cluster
      inquire(unit=kohonen_map%parameters(1)%iclus,opened=testop);
      if(testop) then
      do iz=1,size(kohonen_map%grid_cluster,3);
            write(kohonen_map%parameters(1)%iclus,*) 'Layer ',iz
            do ix=1,size(kohonen_map%grid_cluster,1);
               write(kohonen_map%parameters(1)%iclus,'(100I5)') (kohonen_map%grid_cluster(ix,iy,iz),&
                     iy=1,size(kohonen_map%grid_cluster,2))
            enddo!ix
         enddo!iz
      endif
   !
   !     do ipattern=1,48
   !        write(*,*) kohonen_map%cells_index(ipattern,:)
   !     enddo
      do ipattern=1,kohonen_map%parameters(1)%number_patterns
         ix=kohonen_map%cells_index(ipattern,1);
         iy=kohonen_map%cells_index(ipattern,2);
         iz=kohonen_map%cells_index(ipattern,3);
         !write(*,*) 'pre= ',ipattern,ix,iy,iz,kohonen_map%parameters(1)%number_patterns
         do ipos=1,size(kohonen_map%cluster_cells_index,1)
            if(kohonen_map%cluster_cells_index(ipos,1) .eq. ix .and. &
               kohonen_map%cluster_cells_index(ipos,2) .eq. iy .and. &
               kohonen_map%cluster_cells_index(ipos,3) .eq. iz) then
               pos=ipos;
               exit
            endif
         enddo!ipos
         ic=kohonen_map%cluster_cells_index(pos,4);
         kohonen_map%cluster_samples(ipattern)=ic;
         !write(*,*) 'clust-train = ',pos,ipattern, ic
         inquire(unit=kohonen_map%parameters(1)%iclus1,opened=testop);
         if(testop) then
            write(kohonen_map%parameters(1)%iclus1,*) ipattern,ix,iy,iz,ic;
         endif
      enddo!ipattern
   !
      do ic=1,size(kohonen_map%cluster_layer);
         !write(unit1,*) 'Cluster= ',ic
         !call kohonen_map%cluster_layer(ic)%print(unit1);
         call kohonen_map%cluster_layer(ic)%get_prototype(centers)
         centers1(:,ic)=centers(:)
      enddo!ic
      inquire(unit=kohonen_map%parameters(1)%icen,opened=testop);
      if(testop) then
      do ix=1,size(centers1,1)
         write(kohonen_map%parameters(1)%icen,*) ix,(centers1(ix,ic),ic=1,size(centers1,2))
      enddo
      endif
   !
      if(kohonen_map%parameters(1)%view_flag) then
      write(6,*) 
      write(6,*) 'TWO LEVEL SOM - Cluster Layer: Find Best Match Unit finished'
      write(6,*)
      endif
   !
   end subroutine train_cluster_layer