	program TriCS_CCL
c	-----------------
c
c	fits data from ccl-files and write *.col file for upals
c	and rafin.dat file for position refining
c

c	replaces "none"
c	8.5.99/21.6.99/15.7.99 Juerg Schefer
c	29.11.99: Lorentz-Korrektur SJ33
c	10.4.2000: flat background  SJ33
c       31.10.2000: transfered to UNIX ZO33
c       22.11.2000: problem with printing is fixed
c	02.02.01: h,k,l may be real ZM33
c     17.06.2010: Lorentz correction for tilt Geometry, SJ33
c
c	link ccl,fit4_shr/opt	! attention: check for last fit_n
c	

	implicit none

	integer         np_max
	parameter       (np_max=1000)
	real 		output(20),th2,omega_cal,omega_fit,chi,phi,yint,yy1,
     *			dint,delta,yoverall,h,k,l,s,temperature,range
	real            omega_err
	real 		output2(20), xx(np_max),lor_cor,fwhm,radian,degrees
        real            bgr, wmult
	logical 	plot,iexp,flat_backgr
	integer		lorentz,comens				! SJ33, June 2010
	integer         next_flag
	integer		iall,iincom,icom,icenter
	character*1	answer,hard
        character*128   input_files
        character*128   output_file1
        character*128   output_file2
        character*128	output_file3			! SJ33, May 2012
        character*128	output_file4			! SJ33, May 2012
	character*16 	out
	character*2	off,weak
	character*30    hkl
	character*12	hkl_int
	character*18	hkl_float
	character*24    datime
	integer         ldt

	integer ntyp,i,j,numor,ii,iout,np
	radian=3.14159265/180.00
	degrees=1./radian

        write(*,1000)
 1000   format (/,' FIT CCL Files: Version 2.08, update Jun. 21, 2010',
     *    //,' Data processing of files *.ccl from TriCS single'
     *      ,' detector mode')

	call sys_setenv('DAT_DEFSPEC', 'TRICS')
	call dat_ask_filelist(input_files, ' ')
	output_file1=input_files
	i=index(output_file1, '.')
	if (i .eq. 0) i=index(output_file1, '[')
	if (i .eq. 0) i=index(output_file1, ',')
	if (i .eq. 0) i=index(output_file1, ' ')
        do j=i-1,1,-1
          if (output_file1(j:j) .eq. '/') then
            goto 1002
	  endif
        enddo
        j=0
1002    if (i .gt. j+1) then
          output_file2=output_file1(j+1:i-1)//'_rafin.dat'
          output_file1=output_file1(j+1:i-1)//'.col'
          output_file3=output_file1(j+1:i-1)//'.comm'
          output_file4=output_file1(j+1:i-1)//'.incomm'
        else
          output_file2='rafin.dat'
          output_file1='upals.col'
	endif
c
	write(*,1004)
1004	format (/,' function type to be used (0=gauss[def],6=strange)',
     *     ' typ = ',$)
	read (5,'(i1)') ntyp
	if (ntyp .eq. 6) then
	  write(*,'(/,'' width factor [1.0] = '',$)')
	  read(5,'(f20.0)') wmult
	  if (wmult .eq. 0.0) wmult=1.0
	else
	  wmult=1.0
	  ntyp=0
	endif

	write(*,1013)
 1013	format (/,' Monitor to be standardized , monitor=',$)
	read (5,'(f20.0)') yoverall
	if (yoverall.le.1.1) yoverall = 100000.

	write(*,1015)
1015	format (/,
     *  ' plots (answer=n [no=def], y [yes,terminal], h [hard] ',$)
	read (5,'(a1)') answer
	plot = (answer.eq.'y' .or. answer .eq.'h')
	if (plot) then 
	  hard='y'
	  if (answer.eq.'y') hard='n'
	endif

	write(*,1019)
1019	format (' Flat Background (y/n, default=y) : ',$)
	read (5,'(a1)')answer
	flat_backgr = (answer.ne.'n' .and. answer .ne.'N')

	if (ntyp .eq. 0) then
	  write(*,1017)
1017	  format
     *	(//,' output intensities: 0=experimental (def), 1=gaussian ',$)
	  read (5,'(i1)') Iout
	  iexp= (iout.eq.0)
	else
	  iexp= .true.
	endif
	if (iexp) then
	  out='experimental    '
	else
	  out='gaussian        '
	end if

	write(*,1018)
1018	format (' Lorentz-Correction to be applied (0 = 4-circle-geometry(default), 1 = tilt-geometry, 2 = not applied) ?: ',$)
	read  (5,'(I1)') lorentz			! SJ33, June 2010, tilt option added (3 options now, def=0)
	if (lorentz.ne.1 .and. lorentz.ne.2) lorentz=0	! SJ33, June 2010
	write(*,1011) input_files,ntyp,wmult,
     *  output_file1,output_file2,out,yoverall,flat_backgr,lorentz
 1011	format (/,' ',78('-'),/,' input from file(s) :',A30,/,
     *  ' function type ',i3,'  width factor ',f6.3,//,
     *  ' output: for UPALS :',A50,/,
     *  '         for RAFIN :',A50,/,
     *  1X,A16,' intensities used'/,
     *   ' monitor: ',f12.1,/,
     *  ' Flat Background (yes/no): ',L33,/,
     *  ' Lorentz-Correction (options: 0 =  4-circle-mode,',
     *  ' 1 = tilt-mode, 2 = none)',/,
     *  ' applied option No. ',i1,/,
     *  ' No. 0 Icor=Iobs*sin(2theta)',/,
     *  ' No. 1 Icor=Iobs*sin(Gama)*cos(nu)',/,
     *  ' ',78('-'),/,' hit return to start ',$)
c
	read (5,'(i1)') ii
C	if (ii.eq.1) goto 504

c
	if (ntyp.ne.0 .and. ntyp.ne.6) ntyp = 0
c
c	this ends the input for the program, the rest goes automatic
c	------------------------------------------------------------
	if (plot) then 
		call sys_setenv('CHOOSER_PAN','9') ! 9 plots per page
		end if
c 
c	loop over all data following:
	iall    =0
	iincom  =0
	icom    =0
	icenter =0
	open (unit=1,file=output_file1,status='unknown')
	open (unit=2,file=output_file2,status='unknown')
	open (unit=3,file=output_file3,status='unknown')
	open (unit=4,file=output_file4,status='unknown')
c	write headers for fullprof
c	write (3,1027)9.99,0,1,1,0.00,0.00,0.00
c	write (4,1028)9.99,0,1,1,0.00,0.00,0.00
c1027	format (' title',/,22H(I6,3I10,2F10.2,4F8.2),  /,F7.4,i5,i6,/,i1,3f5.2)
c1028	format (' title',/,24H(I6,3F10.4,2F10.2,4F8.2),/,F7.4,I5,i6,/,i1,3f5.2)
	next_flag=0 ! start flag
	i=0
	call fit_init
c
c	------------------------------------------------------ start loop --
1500	continue

	call fit_dat_next(input_files, next_flag)
	if (next_flag .eq.0) goto 2222 ! last scan of last file read
	i=i+1

	call fit_mon(yoverall)

C	set limits if peak to narrow:
	call fit_get_array('P', output, 20, ii)
	call fit_get_array('X', xx, np_max, np)

	
c	get information out of title:
	call fit_get_real('two_theta',th2)
	call fit_get_real('omega',omega_cal)
	call fit_get_real('chi',chi)
	call fit_get_real('phi',phi)
	call fit_get_real('h',h)
	call fit_get_real('k',k)
	call fit_get_real('l',l)
	call fit_get_real('temp',temperature)	! SJ33, Aug.30,1999
	call fit_get_str('date', ldt, datime)

cjs	type *,' ***** temp=',temperature
	numor=i
	if (l.eq.88) then
	  write(*,1007)
 1007	  format (' data numor not found, exit')
	  goto 3000		
	endif
c
	if (ntyp .eq. 6) then
	  call fit_fun(6,1,wmult,0.0) ! select strange
	else
	  call fit_fun(ntyp,0,0.0,0.0) ! select other (guassian) 
	endif
	if (flat_backgr) call fit_set(2,0.,0.,0.,0.) ! flat background
c	range = scan with:
	range=xx(np)-xx(1)
cjs	type *,'******',output(5),output(6)
	if (ntyp .ne. 6) then
	  if (output(6) .lt. 5.0*range/np) then
c	    set position to center:
	    call fit_set(3, (xx(np)+xx(1))/2.0, range/np,
     &	      xx(1)+range/6,xx(np)-range/6)
c	    set fwhm 
            fwhm = max(0.2,range/5)
c           arbitrary limits: minimal 3 steps, maximal 1 /3 range
            call fit_set(6,fwhm,fwhm/10.,3*range/np, range/3)
	  else
	    fwhm=max(output(6),5*range/np)
	    call fit_set(6,fwhm,range/np
	1                ,3*range/np,range/3)
	  endif
	endif
	  call fit_fit(0)
c
	if (plot) call fit_plot(hard) 		! hardcopy / terminal
c
	call fit_get_array('P', output, 20, ii)
	call fit_get_array('E', output2, 20, ii)
cjs	type *,' output ',output
cjs	type *,' output2',output2
c
c
c	get results now:
	if (ntyp.eq.0) then
c	  function 1 (gauss: use only background, intensity=sum)
	  omega_fit = output (3)
	  omega_err = output2(3)
	  off  = '  '
	  if (iexp) then
	    yint  = output(8)
	    dint = output2(8)	
	  else
c	    modified 14.7.99 sj33
	    yint  = output(5)
	    dint = output2(5)	
	  endif
          bgr=output(1)
          fwhm=output(6)
	else
c	  function 6 (strange)
	  yint  = output(7)
	  dint = output2(7)	
	  omega_fit = output(3)
	  omega_err = output2(3)
          bgr=output(1)
          fwhm=output(4)
	end if
c
c	make notes for easy data inspections:
	delta = abs (omega_cal - omega_fit)
	if (delta.ge.0.3) off =' *'
	if (delta.ge.0.5) off ='**'
	weak = '  '
c
c	write on output files now:
	lor_cor = 1.00000
c	L = 1 / SIN (2-Theta), page 156, Schwarzenbach, EPFL:
c	L applies to the calculated value of the intensities
c	Ical*L = Iobs
c	(measured intensity is always bigger than corrected one)
	if (lorentz.eq.0) lor_cor = sin (abs(th2)*radian)
c
c	June 2010: SJ33
c	tilt angle nu is read-in as chi from file *.col
c	correction formula: not used (3.59) on page 154, Schwarzenbach/Chapuis Cristallographie, Presses EPFL	! SJ33, June 2010
c	used: Arndt and Wills, L=sin(gamma)*cos(nu), where gamma projection 2-theta, nu tilt
	yy1=lor_cor
	if (lorentz.eq.1) lor_cor = 
     *	sin(abs(th2)*radian)*cos(abs(chi)*radian)
c     *	  sqrt( (sin(abs(th2)*radian))**2 -(sin(abs(chi)*radian))**2 )     					! SJ33, June 2010
c	write (6,*)' chi=',chi,th2,radian,sin(chi*radian),sin(th2*radian)
c	write (6,*)' Lorentz=',yy1,lor_cor
c
c
        if (abs(nint(h)-h) .gt. abs(h)*1e-4 .or.
     &      abs(nint(k)-k) .gt. abs(k)*1e-4 .or.
     &      abs(nint(l)-l) .gt. abs(l)*1e-4) then
	  write(hkl, '(3(x,f9.4))')    h, k, l
	  write(hkl_float,'(3f6.2)')  h, k, l
	  iincom=iincom+1
	  iall=iall+1
	  comens=1			! SJ33-May 2012
	else
	  write(hkl, '(3(x,i9))') nint(h), nint(k), nint(l)
	  write(hkl_int,'(3i4)') nint(h), nint(k), nint(l)
	  iall=iall+1
	  icom=icom+1				! SJ33-May 2012
	  comens=0
	endif
c
        write (1, 1005) numor,hkl,lor_cor*yint,lor_cor*dint,
     &	  Th2/2.,omega_fit,chi,phi,temperature,yoverall/1000.,off,weak
     &                    ,bgr,fwhm,lor_cor,omega_err,datime(1:ldt)
c
     	if (comens.eq.0) then
     	write (3, 1025) numor,hkl_int,  lor_cor*yint,lor_cor*dint,
     *		Th2/2,omega_fit,chi,phi			! SJ33-May 2012
     	else
     	write (4, 1026) numor,hkl_float,lor_cor*yint,lor_cor*dint,
     *		Th2/2,omega_fit,chi,phi	
     	endif 		! SJ33-May 2012
1025	format (I6,A12,2F10.2,4f8.2)       	! SJ33, nuc
1026	format (I4,A18,2F10.2,4F8.2)		! SJ33 mag	              
1005	format (I6,a,1x,f9.2,1x,f9.3,f7.2,1x,f8.3,1x,f7.2,1x,f7.2,1x,
     &    F8.3,1x,F6.0,1x,2A2,1x,f6.1,1x,f5.3,1x,f6.4,1x,f6.3,1x,a)

	if (dint.ge.1) then
	  s = yint/dint
	  if (s.le.3) weak='w '	! weak
	  if (s.le.1) weak='vw'	! very weak
	  if (s.ge.15 .and. yint.ge.100) then ! write strong hkl to rafin.dat 
     	  write (2,1009) hkl,Th2,omega_fit,chi,phi,.2,temperature
     	  icenter=icenter+1
     	  end if
c
 1009	  format (a,f8.2,f8.3,2f8.2,1x,f5.3,f9.3)
	end if
C
 	goto 1500
c	-------------------------------------------------------- end loop --
C
 2222	write(*,1006)
 1006	format (' Bye bye, normal end of program.')
C
C	-------------------------------------------
C	error exits here, does not pass label 2222:
 3000	continue
	close (1)
	close (2)
	close (3)			! SJ33-May 2012
	close (4)			! SJ33-May 2012
c	
c	writing the end message for the user
	write (6,2040)
	write (6,2042) input_files(1:60)
	write (6,2041) iall,output_file1(1:30),icenter,output_file2(1:30),
     *	icom,output_file3(1:30),iincom,output_file4
2040	format(///,' trics_ccl prepared the following files for YOU:',//)
2042	format(' TriCS-input-file used:',A60,//)
2041	format(i6,' reflections ',a30,' for JANA2006',/
     *  ,i6,' positions   ',a30,' for rafin (to improve your UB using the dataset)',/
     *  ,I6,' reflections ',a30,' for Fullprof - commensurate data',/
     *  ,I6,' reflections ',a30,' for Fullprof - incommensurate data',//)
c
	call fit_exit
c
	end
