!
! Optimize choice of two parameters $KSOL and $BSOL for use in 
! bulk solvent correction using Babinet principal. 
! The correction to be applied is:
!
!		Fc = (1 - Ksol exp(-Bsol*s**2))Fc
!
! We will loop over several possible values of $BSOL, and dump the resulting
! Luzzati plots to a series of output files defined by the name
! $KBSOLNAME.  These can be plotted using gnuplot, and the "flattest"
! result taken as indicating optimal choice of $BSOL for this value
! of $KSOL.
!
eval ($TRIAL_KSOL = 0.7)
!
! We just set the value of KSOL to be used in this run.
! You could write another script to loop over different values of
! $KSOL, or just edit this one and run it again with a new $TRIAL_KSOL.
!
! Ethan A Merritt - 1995

eval	($CHAINID_IS_SEGID = 1)

structure	@ generate.psf end
coordinates	@ current.protein

{* Add in discrete solvent model *}
 		@ solvent_model.inp
coordinates	@ current.solvent

delete		selection=(not known) end

		@xprepare.inp

xrefine
	resolution 1.8 100.
end

!
! Each run on the bulk solvent correction will fiddle with B factors
! so store the original ones to reset between trials
!
	vector do (BCOMP = B) (all)

!
! Document starting conditional
!
	eval ($KSOL = 0.0)
	eval ($BSOL = 0.0)
	eval ($KBSOLNAME = "original")
	@bulksol.inp
	vector do (B = BCOMP) (all)

!
! Try out bulk solvent model with several possible values of Bsol
!
	eval ($KSOL = $TRIAL_KSOL)

	eval ($BSOL = 25.)
	eval ($KBSOLNAME = "b25")
	@bulksol.inp
	vector do (B = BCOMP) (all)

	eval ($BSOL = 50.)
	eval ($KBSOLNAME = "b50")
	@bulksol.inp
	vector do (B = BCOMP) (all)

	eval ($BSOL = 75.)
	eval ($KBSOLNAME = "b75")
	@bulksol.inp
	vector do (B = BCOMP) (all)

	eval ($BSOL = 100.)
	eval ($KBSOLNAME = "b100")
	@bulksol.inp
	vector do (B = BCOMP) (all)

stop