initial cleanup of DIC code

Source/EMsoftOOLib/mod_DIC.f90

@@ -51,6 +51,8 @@ module mod_DIC
     integer(ip)                     :: ky = 5     ! spline order
     integer(kind=irg)               :: nx         ! pattern x-size
     integer(kind=irg)               :: ny         ! pattern y-size
+    integer(kind=irg)               :: nxc        ! actual size of x-coordinate array
+    integer(kind=irg)               :: nyc        ! actual size of y-coordinate array
     real(wp)                        :: rnxi       ! pattern x scale factor
     real(wp)                        :: rnyi       ! pattern y scale factor
     real(wp)                        :: aspectratio! smallest n over largest n
@@ -65,6 +67,7 @@ module mod_DIC
     logical                         :: verbose = .FALSE.   ! useful for debugging
     logical                         :: normalizedcoordinates = .TRUE.
 
+
 ! the arrays that are tied to a given reference/target pattern are defined here
     real(wp),allocatable            :: x(:)
     real(wp),allocatable            :: y(:)
@@ -138,24 +141,26 @@ module mod_DIC
 contains
 
 !--------------------------------------------------------------------------
-recursive type(DIC_T) function DIC_constructor( nx, ny, normalize) result(DIC)
+recursive type(DIC_T) function DIC_constructor( nx, ny, normalize ) result(DIC)
 !DEC$ ATTRIBUTES DLLEXPORT :: DIC_constructor
  !! author: MDG
  !! version: 1.0
  !! date: 12/01/24
  !!
  !! constructor for the DIC_T Class;
+ !!
+ !! only normalized coordinates have been tested so far...
 
 use bspline_kinds_module, only: wp, ip 
 
 IMPLICIT NONE
 
-integer(kind=irg), INTENT(IN)     :: nx
-integer(kind=irg), INTENT(IN)     :: ny
-logical,INTENT(IN),OPTIONAL       :: normalize
+integer(kind=irg), INTENT(IN)               :: nx
+integer(kind=irg), INTENT(IN)               :: ny
+logical,INTENT(IN),OPTIONAL                 :: normalize
 
-integer(kind=irg)                 :: i, j
-real(wp)                          :: ratio=1.0_wp 
+integer(kind=irg)                           :: i, j
+real(wp)                                    :: ratio=1.0_wp 
 
 DIC%normalizedcoordinates = .FALSE.
 
@@ -164,24 +169,18 @@ recursive type(DIC_T) function DIC_constructor( nx, ny, normalize) result(DIC)
 DIC%ny = ny 
 
 ! allocate and initialize the normalized coordinate arrays
-allocate( DIC%x(0:nx-1), DIC%y(0:ny-1) )
+allocate( DIC%x(0:DIC%nx-1), DIC%y(0:DIC%ny-1) )
+
 DIC%x = (/ (real(i,wp),i=0,nx-1) /)
 DIC%y = (/ (real(j,wp),j=0,ny-1) /)
 
 if (present(normalize)) then  
   if (normalize.eqv..TRUE.) then ! we use normalized coordinates
-    DIC%rnxi = 1.0_wp/real(nx-1,wp)
-    DIC%rnyi = 1.0_wp/real(ny-1,wp)
+    DIC%rnxi = 1.0_wp/real(DIC%nx-1,wp)
+    DIC%rnyi = 1.0_wp/real(DIC%ny-1,wp)
     DIC%x = DIC%x * DIC%rnxi
     DIC%y = DIC%y * DIC%rnyi
     DIC%normalizedcoordinates = .TRUE.
-    ! if (nx.gt.ny) then 
-    !   ! ratio = real(ny,wp) / real(nx,wp)
-    !   DIC%y = DIC%y * ratio
-    ! else
-    !   ! ratio = real(nx,wp) / real(ny,wp)
-    !   DIC%x = DIC%x * ratio
-    ! end if
   end if
 end if
 
@@ -488,8 +487,6 @@ recursive subroutine defineSR_(self, nbx, nby, PCx, PCy)
     self%xiY = self%y(nby:self%ny-nby-1) - PCy
   end if
 end if
-write (*,*) 'defineSR : ',minval(self%xiX), maxval(self%xiX)
-write (*,*) 'defineSR : ',minval(self%xiY), maxval(self%xiY)
 if (self%verbose) call Message%printMessage(' defineSR_::xiX, xiY arrays allocated')
 
 ! allocate array for the product of the gradient and the Jacobian
@@ -587,6 +584,7 @@ recursive subroutine applyHomography_(self, h, PCx, PCy, dotarget)
 
 lnx = self%nx
 lny = self%ny
+
 if (self%aspectratio.eq.1.0_wp) then 
   lPCx = PCx 
   lPCy = PCy 

Source/EMsoftOOLib/program_mods/mod_HREBSDDIC.f90

@@ -636,13 +636,17 @@ subroutine HREBSD_DIC_(self, EMsoft, progname, HDFnames)
     end if
   end do
 
-  if (mod(ii,25).eq.0) then 
+  if (mod(ii,250).eq.0) then 
     io_int(1) = ii
     io_int(2) = numpats
     call Message%WriteValue(' completed # patterns/total ',io_int,2)
   end if  
   hg = DIC%getHomography(W)
-  homographies(1:8,ii) = dble(hg)
+  if (jj.eq.enl%maxnumit+1) then  ! zero solution if no convergence is reached
+    homographies(1:8,ii) = (/ (0.0_wp, i=1,8) /)
+  else
+    homographies(1:8,ii) = dble(hg)
+  end if
   normdp(ii) = dble(ndp)
   residuals(ii) = CIC
   nit(ii) = jj

Source/TestPrograms/play.f90

@@ -70,7 +70,7 @@ program EMplay
                                 minx, miny, xi1max, xi2max, normdp, oldnorm, oldW(3,3), horiginal(8), CIC, sol(8), &
                                 homographies(8,1000), hpartial(8), scalingfactor
 real(kind=dbl)               :: Wnew(3,3), Winv(3,3), dx, dy, p2(3), Woriginal(3,3), alp, srt(3,3), srtrot(3,3)
-integer(kind=irg)            :: nx, ny, nxy, nbx, nby, i, ii, j, NSR, cnt, nxSR, nySR, jj, recordsize, ierr  
+integer(kind=irg)            :: nx, ny, nxy, nbx, nby, i, ii, j, NSR, cnt, nxSR, nySR, jj, recordsize, ierr, maxnumit  
 real(wp)                     :: tol
 integer(kind=4)              :: hnstat
 character(fnlen)             :: fname, gname, hostname 
@@ -156,6 +156,8 @@ program EMplay
 horiginal = (/ (0.0_wp, i=1,8) /)
 call DIC%applyHomography(horiginal, PCx, PCy)
 
+maxnumit = 50 
+
 do jj=1, 1000
     ! call Message%printMessage(' ---------------------- ')
     if (mod(jj,100).eq.0) write (*,*) 'starting pattern ', jj
@@ -189,7 +191,7 @@ program EMplay
 scalingfactor = 1.5D0  
 
 ! and here we start the loop 
-do ii=1,50 
+do ii=1, maxnumit
     ! write (*,*) ' iteration # ',ii
     if (ii.eq.1) then  ! initialize to identity homography in first cycle
         hpartial = (/ (0.0_wp, i=1,8) /)
@@ -244,7 +246,11 @@ program EMplay
 ! write (*,*) horiginal-hg
 
 ! write results to data file (single precision because IDL has a bug for double precision)
-write (unit=28,FMT='(9(F12.8,","),I4)') real(hg), real(normdp), ii
+if (ii.eq.maxnumit+1) then 
+    write (unit=28,FMT='(9(F12.8,","),I4)') (/ (0.0, i=1,8) /), real(normdp), ii
+else
+    write (unit=28,FMT='(9(F12.8,","),I4)') real(hg), real(normdp), ii
+end if
 
 ! if (jj.eq.3) stop
   call DIC%cleanup()