From 174928d255e30000147c008c4c376c83cf63b239 Mon Sep 17 00:00:00 2001 From: Marc DeGraef Date: Mon, 20 Jan 2025 15:00:54 -0500 Subject: [PATCH] initial cleanup of DIC code --- Source/EMsoftOOLib/mod_DIC.f90 | 34 +++++++++---------- .../program_mods/mod_HREBSDDIC.f90 | 8 +++-- Source/TestPrograms/play.f90 | 12 +++++-- 3 files changed, 31 insertions(+), 23 deletions(-) diff --git a/Source/EMsoftOOLib/mod_DIC.f90 b/Source/EMsoftOOLib/mod_DIC.f90 index 83064eb..4c7650d 100644 --- a/Source/EMsoftOOLib/mod_DIC.f90 +++ b/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 diff --git a/Source/EMsoftOOLib/program_mods/mod_HREBSDDIC.f90 b/Source/EMsoftOOLib/program_mods/mod_HREBSDDIC.f90 index 344f825..eef70c9 100644 --- a/Source/EMsoftOOLib/program_mods/mod_HREBSDDIC.f90 +++ b/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 diff --git a/Source/TestPrograms/play.f90 b/Source/TestPrograms/play.f90 index 0ee5fa5..86c4c61 100644 --- a/Source/TestPrograms/play.f90 +++ b/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()