part 1 of forced matching (34a6fd8c) · Commits · Zolnierczuk, Piotr / DrSPINE

sources/matching.f90

+128 −63

Original line number	Diff line number	Diff line
		@@ -4,6 +4,7 @@
		!! With these the proper resolution and backgrounds shall be found
		!!
		module matching
		use, intrinsic :: iso_fortran_env, only : output_unit
		use drspine_parameters
		use constants_module
		use data_types
		@@ -34,7 +35,9 @@ module matching
		0.1*ANGSTROEM &
		)

		public :: taupoint_match, taupoint_mismatch, match_scans, cformat_match_accuracy
		public :: taupoint_match, taupoint_mismatch, cformat_match_accuracy
		public :: clear_matches, list_matching_ids
		public :: match_scans, match_show
		public :: matching_tolerance

		contains
		@@ -80,14 +83,16 @@ contains

		taupoint_match = .false.

		! field integral
		if( diff_values (a%physics%field_integral, b%physics%field_integral) > matching_tolerance%J0_abs) return
		if( diff_relative(a%physics%field_integral, b%physics%field_integral) > matching_tolerance%J0_rel) return

		! scattering and sample angles
		if( diff_values(a%physics%scattering_angle, b%physics%scattering_angle) > matching_tolerance%scattering_angle) return
		if( diff_values(a%physics%sample_angle, b%physics%sample_angle) > matching_tolerance%sample_angle) return

		! wavelength spectrum
		if( a%spectrum(0)%no_lambda_bins .ne. b%spectrum(0)%no_lambda_bins ) return

		if(abs( (spectrum_min_lambda(a%spectrum(0))-spectrum_min_lambda(b%spectrum(0))) ) > matching_tolerance%lambda_min) return
		if(abs( (spectrum_max_lambda(a%spectrum(0))-spectrum_max_lambda(b%spectrum(0))) ) > matching_tolerance%lambda_max) return

		@@ -136,6 +141,75 @@ contains



		!moved from data_types.f90, they really belong here
		!>
		!
		subroutine add_matching_ptr( matches, ptr_to_add )
		implicit none
		type(scan_data_struct_pointer) , intent(inout) :: matches(:)
		type(scan_data_struct), pointer, intent(inout) :: ptr_to_add
		!
		integer :: i

		if(.not.associated(ptr_to_add)) then
		call msg_warning('add_matching_ptr', 'invalid pointer (not associated)')
		return
		endif

		d1: do i=1, size(matches)
		if(associated(matches(i)%ptr,ptr_to_add)) cycle
		if(.not.associated(matches(i)%ptr) ) then
		matches(i)%ptr => ptr_to_add
		return
		endif
		enddo d1
		call msg_error('matching', 'to many matches to one point ', ERROR_DATA_PROCESSING)

		end subroutine add_matching_ptr


		subroutine clear_matches( matches )
		implicit none
		type(scan_data_struct_pointer) , intent(inout) :: matches(:)

		integer :: i

		do i=1, size(matches)
		matches(i)%ptr => null()
		enddo

		end subroutine clear_matches


		function list_matching_ids(matches, full) result(clist)
		implicit none
		integer, parameter :: idfmt_width = 10
		type(scan_data_struct_pointer) , intent(in) :: matches(:)
		logical, optional :: full
		character(len=(idfmt_width+1)*size(matches)) :: clist
		!
		character(len=idfmt_width) :: idbuf
		logical :: full_match
		integer :: i

		clist = " "
		full_match = .false.
		if (present(full)) full_match = full
		d1: do i=1, size(matches)
		if(.not. associated(matches(i)%ptr) ) then
		cycle d1
		endif
		if (full_match) then
		write(idbuf,'(i0,":",i0.2)') matches(i)%ptr%id, matches(i)%ptr%point
		clist = trim(clist)//trim(idbuf)//" "
		else
		write(idbuf,'(i0)') matches(i)%ptr%id
		clist = trim(clist)//" "//idbuf
		end if
		enddo d1

		end function list_matching_ids




		@@ -150,16 +224,16 @@ contains

		cresult = ' '
		del = diff_values(apoint%physics%field_integral, bpoint%physics%field_integral)
		write(cmsg,'(a,e10.3)') 'dJ=', del
		write(cmsg,'(a,e10.3)') 'dJ/Tm=', del
		cresult = trim(cresult)//' '//trim(cmsg)

		del = diff_values(apoint%physics%scattering_angle, bpoint%physics%scattering_angle)
		write(cmsg,'(a,e10.3)') 'dphi/deg=', RAD2DEG(del)
		cresult = trim(cresult)//' '//trim(cmsg)
		cresult = trim(cresult)//', '//trim(cmsg)

		del = diff_values(apoint%physics%sample_angle, bpoint%physics%sample_angle)
		write(cmsg,'(a,e10.3)') 'dpsi/deg=', RAD2DEG(del)
		cresult = trim(cresult)//' '//trim(cmsg)
		cresult = trim(cresult)//', '//trim(cmsg)

		end function cformat_match_accuracy

		@@ -172,12 +246,8 @@ contains
		type(scan_struct), intent(in), target :: target_data ! resolution, background
		!
		integer :: i, j
		character(len=MAX_NAME_LENGTH) :: cmsg

		!! integer :: l !>>test !! (paz) unused variables

		source_loop: do i=1, source_data%number_of_points
		! >>new>> mm
		if(has_role(target_data%role, ROLE_REFERENCE)) then
		if(.not. associated( source_data%scan_point(i)%matching_ref)) &
		call clear_matches(source_data%scan_point(i)%matching_ref_arr)
		@@ -185,55 +255,33 @@ contains
		if(.not. associated( source_data%scan_point(i)%matching_bgr)) &
		call clear_matches(source_data%scan_point(i)%matching_bgr_arr)
		endif
		!!TEST!! write(,)"TESTTAUMATCH 1 : clear ", i !!TEST!!
		! <<new <<
		target_loop: do j=1, target_data%number_of_points
		! FIXME: implicit assumption that the target role is either ROLE_REFERENCE or ROLE_BACKGROUND
		if (taupoint_match(source_data%scan_point(i), target_data%scan_point(j))) then
		! FIXME: implicit assumption that role is either ROLE_REFERENCE or ROLE_BACKGROUND HERE
		if(has_role(target_data%role, ROLE_REFERENCE)) then
		source_data%scan_point(i)%matching_ref => target_data%scan_point(j)
		! >>new>> mm
		call add_matching_ptr(source_data%scan_point(i)%matching_ref_arr, &
		source_data%scan_point(i)%matching_ref )
		!!TEST!! write(,)"TESTTAUMATCH 2 : ref ", i, j, target_data%scan_point(j)%id, & !!TEST!!
		!!TEST!! source_data%scan_point(i)%matching_ref_arr(1)%ptr%id !!TEST!!
		!!!TEST!!>TEST
		!!TEST!! do l=1,size(source_data%scan_point(i)%matching_ref_arr)
		!!TEST!! if(.not.associated(source_data%scan_point(i)%matching_ref_arr(l)%ptr)) then
		!!TEST!! write(,) l, " not assoc"
		!!TEST!! else
		!!TEST!! write(,) l,source_data%scan_point(i)%matching_ref_arr(l)%ptr%id
		!!TEST!! endif
		!!TEST!! enddo
		!!!TEST!!>TEST
		! <<new<<
		cmsg = 'res'
		else if(has_role(target_data%role, ROLE_BACKGROUND)) then
		source_data%scan_point(i)%matching_bgr => target_data%scan_point(j)
		! >>new>> mm
		call add_matching_ptr(source_data%scan_point(i)%matching_bgr_arr, &
		source_data%scan_point(i)%matching_bgr )
		! <<new<<
		cmsg = 'bgr'
		else
		cycle target_loop
		end if
		call msg_debug('tau_match_data', trim(cmsg)//' match found '//&
		'source/'//trim(cformat_role(source_data%role))//&
		trim(msg_fmtint('(i5)',source_data%id))//'/'//&
		adjustl(trim(msg_fmtint('(i5)',i)))//' => '//&
		'target/'//trim(cformat_role(target_data%role))//&
		trim(msg_fmtint('(i5)',target_data%id))//'/'//&
		adjustl(trim(msg_fmtint('(i5)',j)))//' '//&
		call msg_debug('tau_match_data','match found: '//&
		'source/'//trim(cformat_role(source_data%role))//' '//&
		trim(msg_fmt("(i0,':',i0)",[source_data%id, i]))//' => '//&
		'target/'//trim(cformat_role(target_data%role))//' '//&
		trim(msg_fmt("(i0,':',i0)",[target_data%id, j]))//' '//&
		trim(cformat_match_accuracy(source_data%scan_point(i), target_data%scan_point(j))))
		cycle source_loop ! found match for this signal point
		cycle source_loop ! found match for this source point
		end if
		end do target_loop
		end do source_loop
		end subroutine tau_match_data



		!> match scans based on role
		!! This subroutine matches scans in scan_data array based on match_role.
		!! @param scan_data - array of scan_struct data
		@@ -244,35 +292,52 @@ contains
		integer, intent(in) :: match_role
		!
		integer :: i, j
		!!integer :: l,lp ! TEST !! (paz), unused variables
		!
		source_loop: do i=1, scan_size
		if ( .not. is_valid_scan(scan_data(i)) ) cycle source_loop
		if ( has_role(scan_data(i)%role, ROLE_REFERENCE) ) cycle source_loop ! resolution cannot be a match source (PAZ: FIXME?)
		if ( has_role(scan_data(i)%role, ROLE_REFERENCE) ) cycle source_loop ! resolution cannot be a matched to source (PAZ: FIXME?)
		target_loop: do j=1, scan_size
		if (i/=j .and. is_valid_scan(scan_data(j), mode=scan_data(i)%mode, role=match_role)) then
		call tau_match_data(scan_data(i), scan_data(j))
		end if
		end do target_loop

		!>TEST
		!!TEST!! do l=1,size(scan_data(i)%scan_point(1)%matching_ref_arr)
		!!TEST!! do lp=1,scan_data(i)%number_of_points
		!!TEST!! if(.not.associated(scan_data(i)%scan_point(lp)%matching_ref_arr(l)%ptr)) then
		!!TEST!! write(,) "TEST TPMS",i, lp, l," not assoc"
		!!TEST!! else
		!!TEST!! write(,) "TEST TPMS", i, lp, l, scan_data(i)%scan_point(lp)%matching_ref_arr(l)%ptr%id
		!!TEST!! endif
		!!TEST!! enddo
		!!TEST!! enddo
		!>TEST

		end do source_loop
		end subroutine match_scans


		!> show current matches
		!! @param scan_data - array of scan_struct data
		!! @param scan_size - size of scan_data
		subroutine match_show(scan_data, scan_size)
		type(scan_struct), dimension(:), target :: scan_data
		integer, intent(in) :: scan_size
		!
		integer :: i, j
		type(scan_data_struct), pointer :: src, tgt



		!
		source_loop: do i=1, scan_size
		if ( .not. is_valid_scan(scan_data(i)) ) cycle source_loop
		if ( has_role(scan_data(i)%role, ROLE_REFERENCE) ) cycle source_loop ! resolution cannot be a matched to source (PAZ: FIXME?)
		tau_loop: do j=1, scan_data(i)%number_of_points
		src => scan_data(i)%scan_point(j)
		write(output_unit,'(1x,a)', advance='no') trim(cformat_role(scan_data(i)%role, short=.true.))//':'
		write(output_unit,'(1x,a)', advance='no') trim(msg_fmt("(i0,':',i0.2)",[src%id, src%point]))
		if(associated( src%matching_ref)) then
		tgt => src%matching_ref
		write(output_unit,'(a)', advance='no') ' => ref: '
		write(output_unit,'(a)', advance='no') trim(msg_fmt("(i0,':',i0.2)",[tgt%id,tgt%point]))//' ['
		write(output_unit,'(a)', advance='no') trim(list_matching_ids(src%matching_ref_arr,.true.))//']'
		end if
		if(associated( src%matching_bgr)) then
		tgt => src%matching_bgr
		write(output_unit,'(a)', advance='no') ' => bgr: '
		write(output_unit,'(a)', advance='no') trim(msg_fmt("(i0,':',i0.2)",[tgt%id,tgt%point]))//' ['
		write(output_unit,'(a)', advance='no') trim(list_matching_ids(src%matching_bgr_arr,.true.))//']'
		end if
		write(output_unit,*)
		end do tau_loop
		end do source_loop
		end subroutine match_show

		end module matching