# Julia S(q,w) Models

Similarly to the Python modules, a S(q,w) Julia module has to define the two functions "TakinInit" and "TakinSqw".

"TakinInit" is called after one or several parameters have changed (for example after each minimisation step in the convolution fitter). It can be used to check if e.g. pre-calculated tables, variables, etc. need to be recalculated.

The "TakinSqw" function receives four floating-point parameters h,k,l, and E and returns a single floating-point value S, the dynamical structure factor. The function is called for every Monte-Carlo point.

All global variables that are defined in an S(q,w) Julia module and that are prefixed with "g_" (for "global") are made available as settable parameters in the convolution dialog and as fit parameters for the convolution fitter.

The minimal interface to Takin is defined as follows (a full example can be found in the subdirectory "examples/sqw_py"): ``

    function TakinInit()
    	# reinitialise variables here
    end

    function TakinSqw(h::Float64, k::Float64, l::Float64, E::Float64)::Float64
    	S = 0.
    	# calculate S here
    	return Float64(S)
    end

A full example is given below (it can also be found in the subdirectory "examples/sqw_jl" of the source code distribution): ``

    #
    # Sample Julia S(q,w) module for ferromagnetic dispersions
    # @author Tobias Weber
    # @license GPLv2
    # @date dec-2016
    #

    # TODO: import constant
    kB = 1.23

    # example dispersion
    function disp_ferro(q, D, offs)
    	return D*q^2. + offs
    end

    # Gaussian peak
    function gauss(x, x0, sig, amp)
    	norm = (sqrt(2.*pi) * sig)
    	return amp * exp(-0.5*((x-x0)/sig)^2.) / norm
    end

    # Bose factor
    function bose(E, T)
    	n = 1./(exp(abs(E)/(kB*T)) - 1.)
    	if E >= 0.
    		n += 1.
    	end
    	return n
    end

    # Bose factor which is cut off below Ecut
    function bose_cutoff(E, T, Ecut=0.02)
    	Ecut = abs(Ecut)

    	b = 0.
    	if abs(E) < Ecut
    		b = bose(sign(E)*Ecut, T)
    	else
    		b = bose(E, T)
    	end

    	return b
    end

    # -----------------------------------------------------------------------------

    #
    # global variables which can be accessed / changed by Takin
    #
    g_G = vec([1., 1., 0.])		# Bragg peak
    g_D = 50.			# magnon stiffness
    g_offs = 0.			# energy gap
    g_sig = 0.02			# linewidth
    g_S0 = 10.			# intensity

    g_inc_sig = 0.02		# incoherent width
    g_inc_amp = 10.			# incoherent intensity

    g_T = 100.			# temperature
    g_bose_cut = 0.02		# Bose cutoff

    # -----------------------------------------------------------------------------

    #
    # the init function is called after Takin has changed a global variable (optional)
    #
    function TakinInit()
    	println("Calling TakinInit")
    end

    #
    # dispersion E(Q) and weight factor (optional)
    #
    function TakinDisp(h::Float64, k::Float64, l::Float64)
    	# momentum
    	Q = vec([h,k,l])
    	# reduced momentum
    	q = vecnorm(Q - g_G)

    	# energy
    	E_peak = disp_ferro(q, g_D, g_offs)
    	# weight
    	w_peak = 1.
    	return [[E_peak, -E_peak], [w_peak, w_peak]]
    end

    #
    # called for every Monte-Carlo point
    #
    function TakinSqw(h::Float64, k::Float64, l::Float64, E::Float64)::Float64
    	#println("Calling TakinSqw(", h, ", ", k, ", ", l, ", ", E, ") -> ", S)
    	Es, ws = TakinDisp(h,k,l)

    	S_p = gauss(E, Es[1], g_sig, g_S0*ws[1])
    	S_m = gauss(E, Es[2], g_sig, g_S0*ws[2])
    	incoh = gauss(E, 0., g_inc_sig, g_inc_amp)

    	b = 1.
    	#b = bose_cutoff(E, g_T, g_bose_cut)
    	S = (S_p + S_m)*b + incoh
    	return Float64(S)
    end

    # -----------------------------------------------------------------------------
    # test
    #
    #println(TakinSqw(1., 1., 0., 0.))
    #