program loops
integer::i,i1,i2,i3,n1,n2,n3,n(4),l,m
integer,allocatable,dimension(:,:)::a,b,c
real::d1,d2
real,allocatable::smoother(:,:)
n1=5
n2=4
n3=3
n4=2
n=(/n1,n2,n3,n4/)
i=product(n)
allocate(a(i,4),b(i,4),c(i,4))



i=0
do i4=1,n4
do i3=1,n3
 do i2=1,n2
  do i1=1,n1
   i=i+1
   a(i,1)=i1;a(i,2)=i2;a(i,3)=i3;a(i,4)=i4
!   write(0,*) i,i1,i2,i3,i4
!   write(0,*) i-n1*(i2-1)-(i3-1)*n2*n1,i1,1+ (i-1)/n1 - (i3-1)*n2,i2,1+ (i-1)/n1/n2,i3
end do;end do; end do;end do

!write(0,*)' '
!do i=1,n1*n2*n3
! i3=1+ (i-1)/n1/n2
! i2=1+ (i-1)/n1 - (i3-1)*n2
! i1=i            -(i3-1)*n2*n1- (i2-1)*n1
!write(0,*) i,i1,i2,i3
!end do

write(0,*)' '
do i=1,n1*n2*n3*n4
 i4=1+ (i-1)/n1/n2/n3
 i3=1+ (i-1)/n1/n2 -(i4-1)*n3
 i2=1+ (i-1)/n1    -(i4-1)*n3*n2   - (i3-1)*n2
 i1=i              -(i4-1)*n3*n2*n1- (i3-1)*n2*n1- (i2-1)*n1
   b(i,1)=i1;b(i,2)=i2;b(i,3)=i3;b(i,4)=i4
!   write(0,*) i,i1,i2,i3,i4
end do
!ii=0
!do ii=1,product(n)
! ihy=1+ (ii-1)/n(1)/n(2)/n(3)
! ihx=1+ (ii-1)/n(1)/n(2) -(ihy-1)*n(3)
! imy=1+ (ii-1)/n(1)    -(ihy-1)*n(3)*n(2)   - (ihx-1)*n(2)
! imx=ii              -(ihy-1)*n(3)*n(2)*n(1)- (ihx-1)*n(2)*n(1)- (imy-1)*n(1)
!   c(ii,1)=imx;c(ii,2)=imy;c(ii,3)=ihx;c(ii,4)=ihy
!!   write(0,*) ii,imx,imy,ihx,ihy
!end do







write(0,*)n,product(n)
write (0,*)'  '
ii=0
do ii=1,product(n)
    ihy=1+(ii-1)/n(1)/n(2)/n(3)
    ihx=1+(ii-1)/n(1)/n(2)     -(ihy-1)*n(3)
    imy=1+(ii-1)/n(1)          -(ihy-1)*n(3)*n(2)     -(ihx-1)*n(2)
    imx=ii                     -(ihy-1)*n(3)*n(2)*n(1)-(ihx-1)*n(2)*n(1)-(imy-1)*n(1)
   c(ii,1)=imx;c(ii,2)=imy;c(ii,3)=ihx;c(ii,4)=ihy
!!    write(0,*) i,ihy,ihx,imy,imx
end do
!write(0,*) sum(a-b)
!write(0,*) (a-b)
!write (0,*)'  '
!write(0,*) sum(a-c)
!write(0,*) (a-c)
!write (0,*)'  '
!write(0,*) sum(c-b)
!write(0,*) (c-b)

l=5;m=1
cx=1+l/2
cy=1+m/2
write (0,*)'smoother  ',cx,cy
allocate(smoother(l,m))

          do i2=1,m
             do i1=1,l
		d1=real(1+l/2-i1)
		d2=real(1+m/2-i2)
                smoother(i1,i2)=sqrt(real(cx-i1)**2+real(cy-i2)**2)+1
             end do
          end do
!          smoother(l/2+1,m/2+1)=1.
write(0,*) 1./smoother
!write(0,*) smoother
end program












