text:
# XXX iteration over all technologies to create variables and constraints
function createTech!(tInt::Int,part::TechPart,prepTech_dic::Dict{Symbol,NamedTuple},parDef_dic::Dict{Symbol,NamedTuple},ts_dic::Dict{Tuple{Int64,Int64},Array{Int64,1}},r_dic::Dict{Tuple{Int64,Int64},Array{Int64,1}},anyM::anyModel)

    cns_dic = Dict{Symbol,cnsCont}()
    newHerit_dic = Dict(:lowest => (:Ts_dis => :avg_any, :R_dis => :avg_any),:reference => (:Ts_dis => :up, :R_dis => :up)) # inheritance rules after presetting
    ratioVar_dic = Dict(:StIn => ("StIn" => "Conv"), :StOut => ("StOut" => "StIn"), :StSize => ("StSize" => "StIn")) # assignment of tech types for ratios stuff

    tech_str = createFullString(tInt,anyM.sets[:Te])
    # presets all dispatch parameter and obtains mode-dependant variables
    modeDep_dic = presetDispatchParameter!(part,prepTech_dic,parDef_dic,newHerit_dic,ts_dic,r_dic,anyM)

    # creates expansion and capacity variables
    createExpCap!(part,prepTech_dic,anyM,ratioVar_dic)

    # connect capacity and expansion variables
    if part.type != :stock createCapaCns!(part,prepTech_dic,cns_dic) end

    # create and control operated capacity variables
    if anyM.options.decomm != :none createOprVarCns!(part,cns_dic,anyM) end
    produceMessage(anyM.options,anyM.report, 3," - Created all variables and prepared all constraints related to expansion and capacity for technology $(tech_str)")

    # create dispatch variables
    createDispVar!(part,modeDep_dic,ts_dic,r_dic,anyM)
    produceMessage(anyM.options,anyM.report, 3," - Created all dispatch variables for technology $(tech_str)")

    # create conversion balance for conversion technologies
    if keys(part.carrier) |> (x -> any(map(y -> y in x,(:use,:stIntOut))) && any(map(y -> y in x,(:gen,:stIntIn))))
        cns_dic[:convBal] = createConvBal(part,anyM)
        produceMessage(anyM.options,anyM.report, 3," - Prepared conversion balance for technology $(tech_str)")
    end

    # create storage balance for storage technologies
    if :stLvl in keys(part.var)
        cns_dic[:stBal] = createStBal(part,anyM)
        produceMessage(anyM.options,anyM.report, 3," - Prepared storage balance for technology $(tech_str)")
    end

    # create capacity restrictions
    createCapaRestr!(part,ts_dic,r_dic,cns_dic,anyM)
    produceMessage(anyM.options,anyM.report, 3," - Prepared capacity restrictions for technology $(tech_str)")

    # create ratio constraints
    if any(map(x -> occursin("ratioEner",string(x)), collectKeys(keys(part.par))))
        createRatioCns!(part,cns_dic,anyM)
        produceMessage(anyM.options,anyM.report, 3," - Prepared constraints controlling energy ratios for technology $(tech_str)")
    end

    # all constraints are scaled and then written into their respective array position
    foreach(x -> scaleCnsExpr!(x[2].data,anyM.options.coefRng,anyM.options.checkRng), collect(cns_dic))

    produceMessage(anyM.options,anyM.report, 2," - Created all variables and prepared constraints for technology $(tech_str)")

    return cns_dic
end

# <editor-fold desc= prepare to create expansion and capacity variables"

# XXX sets dimensions for expansion and capacity variables
function prepareTechs!(techSym_arr::Array{Symbol,1},prepVar_dic::Dict{Symbol,Dict{Symbol,NamedTuple}},tsYear_dic::Dict{Int,Int},anyM::anyModel)
	for tSym in techSym_arr
		prepTech_dic = Dict{Symbol,NamedTuple}()
		part = anyM.parts.tech[tSym]
        tInt = techInt(tSym,anyM.sets[:Te])

	    # dimension of expansion and corresponding capacity variables
	    if part.type != :stock
	        prepareExpansion!(prepTech_dic, tsYear_dic, part, tInt, anyM)

			for expan in collectKeys(keys(prepTech_dic))
				prepareCapacity!(part,prepTech_dic,vcat(map(x -> x[!,removeVal(x)],prepTech_dic[expan])...),Symbol(replace(string(expan),"exp" => "capa")),anyM, tech = tInt)
			end
		end

		# check for capacities variables that have to be created, because of residual capacities provided
		addResidualCapa!(prepTech_dic, part, tInt, anyM)

		# map required capacity constraints
		createCapaRestrMap!(tSym, anyM)

		# if any capacity variables or residuals were prepared, add these to overall dictionary
	    if collect(values(prepTech_dic)) |> (z -> any(map(x -> any(.!isempty.(getfield.(z,x))), (:var,:resi))))
			prepVar_dic[tSym] = prepTech_dic
			# write reporting, if not all 3 kind of storage capacities will be created
			if map(x -> x in keys(prepTech_dic),(:capaStIn,:capaStOut,:capaStSize)) |> (y -> any(y) && !all(y))
				push!(anyM.report,(3,"technology dimensions","storage","in case of $(string(tSym)) information for one storage capacity is missing (capaStIn, capaStOut or capaStSize)"))
			end
		end
	end
end

# XXX dimensions for expansion variables
function prepareExpansion!(prepTech_dic::Dict{Symbol,NamedTuple},tsYear_dic::Dict{Int,Int},part::AbstractModelPart,tInt::Int,anyM::anyModel)

	# extract tech info
	carGrp_ntup = part.carrier
	balLvl_ntup = part.balLvl
	defPar_tup = tuple(keys(part.par)...)

	tsExp_arr, rExp_arr   = [getfield.(getNodesLvl(anyM.sets[x[2]], balLvl_ntup.exp[x[1]]),:idx) for x in enumerate([:Ts,:R])]
	tsExpSup_arr = map(x -> getDescendants(x,anyM.sets[:Ts],false,anyM.supTs.lvl) |> (y -> typeof(y) == Array{Int,1} ? y : [y] ), tsExp_arr)
	if anyM.options.interCapa != :linear tsExp_arr = map(x -> [minimum(x)],tsExp_arr) end

	expDim_arr = vcat(collect(Iterators.product(Iterators.zip(tsExp_arr,tsExpSup_arr),rExp_arr))...)
	allMap_df =  getindex.(expDim_arr,1) |> (x -> DataFrame(Ts_exp = getindex.(x,1), Ts_expSup = getindex.(x,2), R_exp = getindex.(expDim_arr,2), Te = fill(tInt,length(expDim_arr))))

	# prepares expansion dimensions for conversion capacity
	if !isempty(intersect((:gen,:use),keys(carGrp_ntup)))
		# filters cases where expansion is fixed to zero
		convMap_df = removeEntries([filterZero(allMap_df,getLimPar(anyM.parts.lim,:expConvFix, anyM.sets[:Te], tech = tInt),anyM)],allMap_df)

		if !isempty(convMap_df) prepTech_dic[:expConv] =  (var = addSupTsToExp(convMap_df,part.par,:Conv,tsYear_dic,anyM), ratio = DataFrame(), resi = DataFrame()) end
	end

	stCar_arr::Array{Int,1} = unique(vcat(collect.(map(x -> getproperty(carGrp_ntup,x),intersect(keys(carGrp_ntup),(:stExtIn,:stExtOut,:stIntIn,:stIntOut))))...))

	# prepares expansion dimensions for storage capacity
	if !isempty(stCar_arr)
		stMap_df = combine(groupby(allMap_df,namesSym(allMap_df)), :Te => (x -> stCar_arr) => :C)

		for st in (:StIn, :StOut, :StSize)
			remove_arr = Array{DataFrame,1}()

			# filters cases where expansion is fixed to zero
			if Symbol(:exp,st,:Fix) in defPar_tup
				push!(remove_arr,filterZero(stMap_df,getLimPar(anyM.parts.lim,Symbol(:exp,st,:Fix), anyM.sets[:Te], tech = tInt),anyM))
			end

			remove_arr, ratioTab_df = findStorageRatio(tInt,stMap_df,st, remove_arr, part, :exp, anyM)
			specStMap_df = removeEntries(remove_arr,stMap_df)

			if !(isempty(specStMap_df) && isempty(ratioTab_df))
				prepTech_dic[Symbol(:exp,st)] =  (var = addSupTsToExp(specStMap_df,part.par,st,tsYear_dic,anyM), ratio = addSupTsToExp(ratioTab_df,part.par,st,tsYear_dic,anyM), resi = DataFrame())
			end
		end
	end
end

# XXX dimensions for capacity variables
function prepareCapacity!(part::AbstractModelPart,prep_dic::Dict{Symbol,NamedTuple},exp_df::DataFrame,capaVar::Symbol,anyM::anyModel; tech::Int = 0)

	# XXX initialize assignments and data
	defPar_tup = tuple(keys(part.par)...)
	techType_sym = :type in fieldnames(typeof(part)) ? part.type : :mature

	capaVar_df = expandExpToCapa(exp_df)

	# groups by expansion time steps in case of mature technologies
	if techType_sym == :mature
		select!(capaVar_df,Not(:Ts_expSup))
		capaVar_df = unique(capaVar_df)
		capaVar_df[!,:Ts_expSup] .= 0
	end

	# filters cases where capacity is fixed to zero
    varFix_sym = Symbol(capaVar,:Fix)

	if varFix_sym in defPar_tup
		capaVar_df = removeEntries([filterZero(capaVar_df,getLimPar(anyM.parts.lim,Symbol(capaVar,:Fix),anyM.sets[:Te], tech = tech),anyM)],capaVar_df)
	end

	# for exchange capacities add column to indicate these values are symmetric
	if capaVar == :capaExc capaVar_df[!,:dir] .= false; select!(capaVar_df,Not(:Ts_expSup)) end

	prep_dic[capaVar] =  (var = orderDf(capaVar_df), ratio = DataFrame(), resi = DataFrame())
end

# XXX checks if a storage capacity is defined via a ratio somewhere
function findStorageRatio(t_int::Int,find_df::DataFrame,st_sym::Symbol,remove_arr::Array{DataFrame,1},part::AbstractModelPart,kind_sym::Symbol,anyM::anyModel)

	expTypeRatio_dic = Dict(:StIn => :stInToConv, :StOut => :stOutToStIn, :StSize => :sizeToStIn)
	strRatio_dic = Dict(:StIn => "storage input", :StOut => "storage output", :StSize => "storage size")
	strLim_dic = Dict(:Up => "upper limits", :Low => "lower limits", :Fix => "fixed limits",:Resi => "residual capacities")

	ratio_sym = expTypeRatio_dic[st_sym]

	if ratio_sym in tuple(keys(part.par)...)
		ratioTab_df = filter(r -> r.ratio != 0, matchSetParameter(find_df,part.par[ratio_sym],anyM.sets, newCol = :ratio))
		push!(remove_arr, ratioTab_df[!,Not(:ratio)])
	else
		ratioTab_df = filter(x -> false,find_df)
	end

	# writes a report, if limits (upper/lower/fixed) on the storage expansion were ignored due to these ratios provided
	if !isempty(ratioTab_df) && (kind_sym == :exp || (kind_sym == :capa && part.type == :stock) )
		for limPar in intersect(keys(part.par),map(x -> Symbol(kind_sym,st_sym,x),kind_sym == :exp ? [:Up,:Low,:Fix] : [:Up,:Low,:Fix,:Resi]))
			lim_obj = getLimPar(anyM.parts.lim,limPar,anyM.sets[:Te], tech = t_int)
			if !(isdefined(lim_obj,:name)) continue end
			both_df =  innerjoin(ratioTab_df,lim_obj.data, on = intersect(intCol(ratioTab_df),intCol(lim_obj.data)))
			if !isempty(both_df)
				push!(anyM.report,(1,:variable,:expansion,"for $(join(part.name, " < ")) $(strLim_dic[Symbol(split(string(limPar),string(st_sym))[end])]) for $(strRatio_dic[st_sym]) were ignored since an conversion/storage input ratio was provided"))
			end
		end
	end

	return remove_arr, ratioTab_df
end

# </editor-fold>

# <editor-fold desc= create technology related variables"

# XXX create expansion and capacity variables
function createExpCap!(part::AbstractModelPart,prep_dic::Dict{Symbol,NamedTuple},anyM::anyModel,ratioVar_dic::Dict{Symbol,Pair{String,String}} = Dict{Symbol,Pair{String,String}}())
	for expVar in sort(collectKeys(keys(prep_dic)))

		varMap_tup = prep_dic[expVar]
		# create dataframe of capacity or expansion variables by creating the required capacity variables and join them with pure residual values
		var_df = createVar(varMap_tup.var,string(expVar),anyM.options.bound.capa,anyM.optModel,anyM.lock,anyM.sets, scaFac = anyM.options.scaFac.capa)
		if !isempty(varMap_tup.resi)
			if expVar == :capaExc # flips and repeats entries for directed exchange variabes before moving on
				var_df = filter(r -> r.dir,var_df) |> (x -> vcat(filter(r -> !r.dir,var_df),vcat(x,rename(x,replace(namesSym(x),:R_to => :R_from, :R_from => :R_to)))))
			end
			join_arr = intCol(var_df,:dir)
			var_df = combine(x -> (var = x.var + x.var_1,), groupby(joinMissing(var_df,varMap_tup.resi[!,vcat(:var,join_arr...)], join_arr, :outer, Dict(:var => AffExpr(),:var_1 => AffExpr()),true),intCol(var_df,:dir)))
		end

		# expands table of expansion variables to superordinate timesteps and modifies expansion variable accordingly
		if occursin("exp",string(expVar)) && !isempty(var_df)
			noExpCol_arr = intCol(var_df)
			allDf_arr = map(eachrow(var_df)) do x
				l_int = length(x.Ts_disSup)
				rem_df = repeat(DataFrame(x[noExpCol_arr]), inner = l_int, outer = 1)
				ext_df = DataFrame(Ts_expSup = x.Ts_expSup, Ts_disSup = x.Ts_disSup, var = x.var ./ fill(l_int,l_int) )
				return hcat(rem_df,ext_df)
			end
			var_df = vcat(allDf_arr...)
		end

		# check for ratios for expansion or (in case of stock technologies) capacities
		if !isempty(varMap_tup.ratio)
			ratioVar_sym = Symbol(replace(string(expVar),ratioVar_dic[Symbol(replace(replace(string(expVar),"exp" => ""),"capa" => ""))]))
			if occursin("exp",string(expVar)) # ratios controlling expansion
				noExpCol_arr = vcat(:ratio,intCol(var_df)...)

				ratio_arr = map(eachrow(varMap_tup.ratio)) do x
					l_int = length(x.Ts_expSup)
					rem_df = repeat(DataFrame(x[noExpCol_arr]), inner = l_int, outer = 1)
					ext_df = DataFrame(Ts_expSup = x.Ts_expSup)
					return hcat(rem_df,ext_df)
				end
				preRatio_df = vcat(ratio_arr...)
			else # ratios controlling stock capacities
				preRatio_df = varMap_tup.ratio
			end

			join_arr = intCol(part.var[ratioVar_sym])
			# join ratios and corresponding
			ratio_df = select(innerjoin(preRatio_df,part.var[ratioVar_sym]; on = join_arr),unique(vcat(join_arr,[:var,:ratio,:Ts_disSup])))
			ratio_df[!,:var] = ratio_df[!,:var] .* ratio_df[!,:ratio]

			var_df = ratio_df[!,Not(:ratio)] |> (x -> isempty(var_df) ? x : vcat(x,antijoin(var_df,x, on = join_arr)))
		end

		if !isempty(var_df)	part.var[expVar] = orderDf(var_df) end
	end
end

# XXX create all dispatch variables
function createDispVar!(part::TechPart,modeDep_dic::Dict{Symbol,DataFrame},ts_dic::Dict{Tuple{Int64,Int64},Array{Int64,1}},r_dic::Dict{Tuple{Int64,Int64},Array{Int64,1}},anyM::anyModel)
	# assign relevant availability parameters to each type of variable
	relAva_dic = Dict(:gen => (:avaConv,), :use => (:avaConv,), :stIntIn => (:avaConv, :avaStIn), :stIntOut => (:avaConv, :avaStOut), :stExtIn => (:avaStIn,), :stExtOut => (:avaStOut,), :stLvl => (:avaStSize,))

	for va in collectKeys(keys(part.carrier)) |> (x -> :capaStIn in keys(part.var) ? [:stLvl,x...]  : x) # loop over all relevant kind of variables
		conv_boo = va in (:gen,:use)
		# obtains relevant capacity variable
		if conv_boo
			basis_df = copy(part.var[:capaConv])[!,Not(:var)]
			basis_df[!,:C] .= [collect(getfield(part.carrier,va))]
			basis_df = orderDf(flatten(basis_df,:C))
		else
            lock(anyM.lock)
			basis_df = orderDf(copy(part.var[:capaStIn])[!,Not(:var)])
			unlock(anyM.lock)
			# filter carriers that are can be actively stored, although they got descendants
			intC_arr = union(collect(part.actSt),map(y -> part.carrier[y],filter(x -> x in keys(part.carrier),[:stIntIn,:stIntOut])) |> (y -> isempty(y) ? Int[] : union(y...)))
			basis_df = replCarLeafs(basis_df,anyM.sets[:C],noLeaf = intC_arr)
			# filter entries that are already descendants of carrier being actively stored
			unique(vcat(map(x -> filter(y -> x != y,getDescendants(x,anyM.sets[:C],true)),unique(basis_df[!,:C]))...)) |> (z -> filter!(x -> !(x.C in z) || x.C in intC_arr,basis_df))
		end

		# adds temporal and spatial level to dataframe
		cToLvl_dic = Dict(x => (anyM.cInfo[x].tsDis, part.disAgg ? part.balLvl.exp[2] : anyM.cInfo[x].rDis) for x in unique(basis_df[!,:C]))
		basis_df[!,:lvlTs] = map(x -> cToLvl_dic[x][1],basis_df[!,:C])
		basis_df[!,:lvlR] = map(x -> cToLvl_dic[x][2],basis_df[!,:C])
		allVar_df = orderDf(expandExpToDisp(basis_df,ts_dic,r_dic,true))

		# add mode dependencies
		modeDep_df = copy(modeDep_dic[va])
		modeDep_df[!,:M] .= isempty(modeDep_df) ? [0] : [collect(part.modes)]
		modeDep_df = flatten(modeDep_df,:M)

		allVar_df = joinMissing(allVar_df,modeDep_df,namesSym(modeDep_dic[va]),:left,Dict(:M => 0))

		# filter entries where availability is zero
		for avaPar in relAva_dic[va]
			if !isempty(part.par[avaPar].data) && 0.0 in part.par[avaPar].data[!,:val]
				allVar_df = filter(x -> x.val != 0.0,  matchSetParameter(allVar_df,part.par[avaPar],anyM.sets))[!,Not(:val)]
			end
		end

		# computes value to scale up the global limit on dispatch variable that is provied per hour and create variable
        if conv_boo
			scaFac_fl = anyM.options.scaFac.dispConv
		else
			scaFac_fl = anyM.options.scaFac.dispSt
		end
		part.var[va] = orderDf(createVar(allVar_df,string(va), getUpBound(allVar_df,anyM.options.bound.disp / scaFac_fl,anyM.supTs,anyM.sets[:Ts]),anyM.optModel,anyM.lock,anyM.sets, scaFac = scaFac_fl))
	end
end

# XXX create variables and constraints regarding operated variables
function createOprVarCns!(part::AbstractModelPart,cns_dic::Dict{Symbol,cnsCont},anyM::anyModel)
	for capaVar in filter(x -> occursin("capa",string(x)),keys(part.var))
		oprVar_sym = string(capaVar) |> (x -> Symbol(:opr,uppercase(x[1]),x[2:end]))
		# XXX create operated variable
		var_df = copy(part.var[capaVar])[!,Not(:var)]
		var_df = createVar(var_df,string(oprVar_sym),NaN,anyM.optModel,anyM.lock,anyM.sets, scaFac = anyM.options.scaFac.oprCapa)

		part.var[oprVar_sym] = orderDf(var_df)

		# XXX create constraint to connect operated and installed capacity
		var_df[!,:cnsExpr] = map(x -> x[1] - x[2],zip(var_df[!,:var],part.var[capaVar][!,:var]))
		cns_dic[oprVar_sym] = cnsCont(select(var_df,Not(:var)),:smaller)

		# XXX create constraint to prevent re-commissioning of capacity once decommissioned
		if anyM.options.decomm == :decomm
			# add previous period and its capacity variable to table
			prevTs_dic = Dict(x => anyM.supTs.step[findall(x .== anyM.supTs.step)[1]]-1 for x in anyM.supTs.step[2:end])
			select!(var_df, Not(:cnsExpr))
			cns_df = rename(filter(r -> r.Ts_disSup != anyM.supTs.step[1],var_df),:var => :oprNow)
			cns_df[!,:Ts_disSupPrev] = map(x -> prevTs_dic[x] ,cns_df[!,:Ts_disSup])
			cns_df = rename(innerjoin(cns_df,var_df; on = intCol(var_df,:dir) |> (x -> Pair.(replace(x,:Ts_disSup => :Ts_disSupPrev),x))),:var => :oprPrev)

			# add expansion variable to dataframe
			if (!(:type in fieldnames(typeof(part))) || part.type != :stock) && Symbol(replace(string(capaVar),"capa" => "exp")) in keys(part.var)
				exp_df = part.var[Symbol(replace(string(capaVar),"capa" => "exp"))][!,Not(:Ts_disSup)]
				join_arr = filter(x -> x != :Ts_expSup,intCol(var_df))

				cns_df = joinMissing(cns_df,exp_df, Pair.(join_arr,replace(join_arr,:Ts_disSup => :Ts_expSup)),:left,Dict(:var => AffExpr(),:Ts_exp => 0))
				cns_df = rename(cns_df[!,Not(:Ts_exp)],:var => :expNow)
			else
				cns_df[!,:expNow] .= AffExpr()
			end

			# add residual capacities of current and previous period
			joinResi_arr = filter(x -> x != :Ts_disSupPrev, intCol(cns_df,:dir))
			cns_df = rename(innerjoin(cns_df,part.var[capaVar],on = joinResi_arr),:var => :resiNow)
			cns_df[!,:resiNow] = getfield.(cns_df[!,:resiNow],:constant)
			cns_df = rename(joinMissing(cns_df, part.var[capaVar], Pair.(replace(joinResi_arr,:Ts_disSup => :Ts_disSupPrev),joinResi_arr),:left, Dict(:resiNow => AffExpr(),:var => AffExpr())),:var => :resiPrev)
			cns_df[!,:resiPrev] = getfield.(cns_df[!,:resiPrev],:constant)

			# create actual constraint information
			cns_df[!,:cnsExpr]  = map(x -> - x.oprNow + x.oprPrev + x.expNow + (x.resiNow - x.resiPrev |> (l -> l > 0.0 ? l : 0.0)),eachrow(cns_df))
			select!(cns_df,Not([:Ts_disSupPrev,:oprNow,:oprPrev,:expNow,:resiNow,:resiPrev]))
			cns_dic[string(oprVar_sym) |> (x -> Symbol(:re,uppercase(x[1]),x[2:end]))] = cnsCont(orderDf(cns_df),:greater)
		end
	end
end

# XXX capacity values for stock technologies
function addResidualCapa!(prepTech_dic::Dict{Symbol,NamedTuple},part::TechPart,tInt::Int,anyM::anyModel)

	carGrp_ntup = part.carrier
	stCar_arr = unique(vcat(collect.(map(x -> getproperty(carGrp_ntup,x),intersect(keys(carGrp_ntup),(:stExtIn,:stExtOut,:stIntIn,:stIntOut))))...))
	defPar_tup = tuple(keys(part.par)...)

	# checks conversion capacities for stock data
	if !isempty(intersect((:gen,:use),keys(carGrp_ntup)))
		permutDim_arr = [getindex.(vcat(collect(Iterators.product(getfield.(getNodesLvl(anyM.sets[:R], part.balLvl.exp[2]),:idx), anyM.supTs.step))...),i) for i in (1,2)]
		potCapa_df = DataFrame(Ts_expSup = fill(0,length(permutDim_arr[1])), Ts_disSup = permutDim_arr[2], R_exp = permutDim_arr[1], Te = fill(tInt,length(permutDim_arr[1])))
		# filters stock technologies with residual values provided and saves these in dictionary
		capaResi_df = checkResiCapa(:capaConv,potCapa_df, part, anyM)

		if !isempty(capaResi_df)
			mergePrepDic!(:capaConv,prepTech_dic,capaResi_df)
		end
	end

	# checks storage capacities for stock data
	if !isempty(stCar_arr)
		permutDim_arr = [getindex.(vcat(collect(Iterators.product(getfield.(getNodesLvl(anyM.sets[:R], part.balLvl.exp[2]),:idx), anyM.supTs.step,stCar_arr))...),i) for i in (1,2,3)]
		potCapa_df = DataFrame(Ts_expSup = fill(0,length(permutDim_arr[1])), Ts_disSup = permutDim_arr[2], R_exp = permutDim_arr[1], C = permutDim_arr[3], Te = fill(tInt,length(permutDim_arr[1])))

		# loops over different type of storage
		for st in (:StIn, :StOut, :StSize)
			capaResi_df = checkResiCapa(Symbol(:capa,st),potCapa_df, part, anyM)
			capa_sym = Symbol(:capa,st)
			remove_arr, ratioTab_df = findStorageRatio(tInt,potCapa_df,st, Array{DataFrame,1}(), part, :capa, anyM)
			if !isempty(capaResi_df) capaResi_df = removeEntries(remove_arr,capaResi_df) end

			if !(isempty(capaResi_df) && isempty(ratioTab_df))
				mergePrepDic!(capa_sym,prepTech_dic,capaResi_df,ratioTab_df)
			end
		end

	end
end

# </editor-fold>

# <editor-fold desc= create technology related constraints"

# XXX create conversion balance
function createConvBal(part::TechPart,anyM::anyModel)

	cns_df = rename(part.par[:effConv].data,:val => :eff)
	sort!(cns_df,sort(intCol(cns_df)))
	agg_arr = filter(x -> !(x in (:M, :Te)) && (part.type == :emerging || x != :Ts_expSup), intCol(cns_df))

	# defines tuple specificing dimension of aggregation later
	if part.type == :emerging
		srcRes_ntup = part.balLvl |> (x -> (Ts_expSup = anyM.supTs.lvl, Ts_dis = x.ref[1], R_dis = x.ref[2]))
	else
		srcRes_ntup = part.balLvl |> (x -> (Ts_dis = x.ref[1], R_dis = x.ref[2]))
	end

	# if modes are specified, gets rows of conversion dataframe where they are relevant and creates different tuples to define grouping dimensions
	if :M in namesSym(cns_df)
		srcResM_ntup = (; zip(tuple(:M,keys(srcRes_ntup)...),tuple(1,values(srcRes_ntup)...))...)
		srcResNoM_ntup = (; zip(tuple(:M,keys(srcRes_ntup)...),tuple(0,values(srcRes_ntup)...))...)
		m_arr = findall(0 .!= cns_df[!,:M])
		noM_arr = setdiff(1:size(cns_df,1),m_arr)
	end

	# add variables via aggregation
	in_arr = intersect(keys(part.carrier),(:use,:stIntOut))
	out_arr = intersect(keys(part.carrier),(:gen,:stIntIn))

	for va in union(in_arr,out_arr)
		if :M in namesSym(cns_df) # aggregated dispatch variables, if a mode is specified somewhere, mode dependant and non-mode dependant balances have to be aggregated separately
			cns_df[!,va] .= AffExpr()
			cns_df[m_arr,va] = aggUniVar(part.var[va], select(cns_df[m_arr,:],intCol(cns_df)), [:M,agg_arr...], srcResM_ntup, anyM.sets)
			cns_df[noM_arr,va] = aggUniVar(part.var[va], select(cns_df[noM_arr,:],intCol(cns_df)), [:M,agg_arr...], srcResNoM_ntup, anyM.sets)
		else
			cns_df[!,va] = aggUniVar(part.var[va], select(cns_df,intCol(cns_df)), agg_arr, srcRes_ntup, anyM.sets)
		end
	end

	# aggregate in and out variables respectively
	cns_df[!,:in] = map(x -> sum(x),eachrow(cns_df[!,in_arr]))
	cns_df[!,:out] = map(x -> sum(x),eachrow(cns_df[!,out_arr]))
	select(cns_df,Not(vcat(in_arr,out_arr)))

	# create actual constraint
	cns_df[!,:cnsExpr] = map(x -> x.in*x.eff - x.out,eachrow(cns_df))
	return cnsCont(orderDf(cns_df[!,[intCol(cns_df)...,:cnsExpr]]),:equal)
end

# XXX create storage balance
function createStBal(part::TechPart,anyM::anyModel)

	# XXX get variables for storage level
	# get variables for current storage level
	cns_df = rename(part.var[:stLvl],:var => :stLvl)
	cnsDim_arr = filter(x -> x != :Ts_disSup, intCol(cns_df))

	# join variables for previous storage level
	tsChildren_dic = Dict((x,y) => getDescendants(x,anyM.sets[:Ts],false,y) for x in anyM.supTs.step, y in unique(map(x -> getfield(anyM.sets[:Ts].nodes[x],:lvl), cns_df[!,:Ts_dis])))
	firstLastTs_dic = Dict(minimum(tsChildren_dic[z]) => maximum(tsChildren_dic[z]) for z in keys(tsChildren_dic))
	firstTs_arr = collect(keys(firstLastTs_dic))

	cns_df[!,:Ts_disPrev] = map(x -> x in firstTs_arr ? firstLastTs_dic[x] : x - 1, cns_df[!,:Ts_dis])
	cns_df = rename(joinMissing(cns_df,part.var[:stLvl], intCol(part.var[:stLvl]) |> (x -> Pair.(replace(x,:Ts_dis => :Ts_disPrev),x)), :left, Dict(:var => AffExpr())),:var => :stLvlPrev)

	# determines dimensions for aggregating dispatch variables
	agg_arr = filter(x -> !(x in (:M, :Te)) && (part.type == :emerging || x != :Ts_expSup), cnsDim_arr)

	# obtain all different carriers of level variable and create array to store the respective level constraint data
	uniC_arr = unique(cns_df[!,:C])
	cCns_arr = Array{DataFrame}(undef,length(uniC_arr))

	for (idx,c) in enumerate(uniC_arr)

		# get constraints relevant for carrier and find rows where mode is specified
		cnsC_df = filter(x -> x.C == c,cns_df)
		sort!(cnsC_df,sort(intCol(cnsC_df)))

		m_arr = findall(0 .!= cnsC_df[!,:M])
		noM_arr = setdiff(1:size(cnsC_df,1),m_arr)

		if part.type == :emerging
			srcRes_ntup = anyM.cInfo[c] |> (x -> (Ts_expSup = anyM.supTs.lvl, Ts_dis = x.tsDis, R_dis = x.rDis, C = anyM.sets[:C].nodes[c].lvl, M = 1))
		else
			srcRes_ntup = anyM.cInfo[c] |> (x -> (Ts_dis = x.tsDis, R_dis = x.rDis, C = anyM.sets[:C].nodes[c].lvl, M = 1))
		end

		# XXX join in and out dispatch variables and adds efficiency to them (hence efficiency can be specific for different carriers that are stored in and out)
		for typ in (:in,:out)
			typVar_df = copy(cns_df[!,cnsDim_arr])
			# create array of all dispatch variables
			allType_arr = intersect(keys(part.carrier),typ == :in ? (:stExtIn,:stIntIn) : (:stExtOut,:stIntOut))
			effPar_sym = typ == :in ? :effStIn : :effStOut
			# adds dispatch variables
			typExpr_arr = map(allType_arr) do va
				typVar_df = filter(x -> x.C == c,part.par[effPar_sym].data) |> (x -> innerjoin(part.var[va],x; on = intCol(x)))
				if typ == :in
					typVar_df[!,:var] = typVar_df[!,:var] .* typVar_df[!,:val]
				else
					typVar_df[!,:var] = typVar_df[!,:var] ./ typVar_df[!,:val]
				end
				return typVar_df[!,Not(:val)]
			end

			# adds dispatch variable to constraint dataframe, mode dependant and non-mode dependant balances have to be aggregated separately
			dispVar_df = vcat(typExpr_arr...)
			cnsC_df[!,typ] .= AffExpr()
			if isempty(dispVar_df) continue end
			cnsC_df[m_arr,typ] = aggUniVar(dispVar_df, select(cnsC_df[m_arr,:],intCol(cnsC_df)), [:M,agg_arr...], (M = 1,), anyM.sets)
			cnsC_df[noM_arr,typ] = aggUniVar(dispVar_df, select(cnsC_df[noM_arr,:],intCol(cnsC_df)), [:M,agg_arr...], (M = 0,), anyM.sets)
		end

		# XXX adds further parameters that depend on the carrier specified in storage level (superordinate or the same as dispatch carriers)
		sca_arr = getResize(cnsC_df,anyM.sets[:Ts],anyM.supTs)

		# add discharge parameter, if defined
		if :stDis in keys(part.par)
			part.par[:stDis].defVal = 0.0
			cnsC_df = matchSetParameter(cnsC_df,part.par[:stDis],anyM.sets)
			cnsC_df[!,:stDis] =   (1 .- cnsC_df[!,:val]) .^ sca_arr
			select!(cnsC_df,Not(:val))
		else
			cnsC_df[!,:stDis] .= 1.0
		end

		# add inflow parameter, if defined
		if :stInflow in keys(part.par)
			part.par[:stInflow].defVal = 0.0
			cnsC_df = matchSetParameter(cnsC_df,part.par[:stInflow],anyM.sets, newCol = :stInflow)
			if !isempty(part.modes)
            	cnsC_df[!,:stInflow] = cnsC_df[!,:stInflow] ./ length(part.modes) .* sca_arr
			end
		else
			cnsC_df[!,:stInflow] .= 0.0
		end

		# XXX create final equation
		cnsC_df[!,:cnsExpr] = map(x -> x.stLvlPrev * x.stDis + x.stInflow + x.in - x.out - x.stLvl,eachrow(cnsC_df))
		cCns_arr[idx] = cnsC_df
	end

	cns_df =  vcat(cCns_arr...)
	return cnsCont(orderDf(cns_df[!,[cnsDim_arr...,:cnsExpr]]),:equal)
end

# XXX create all capacity restrictions for technology
function createCapaRestr!(part::TechPart,ts_dic::Dict{Tuple{Int64,Int64},Array{Int64,1}},r_dic::Dict{Tuple{Int64,Int64},Array{Int64,1}},cns_dic::Dict{Symbol,cnsCont},anyM::anyModel)

	cnstrType_dic = Dict(:out => (dis = (:gen, :stIntIn), capa = :Conv), :in => (dis = (:use,:stIntOut), capa = :Conv),
							:stIn => (dis = (:stExtIn, :stIntIn), capa = :StIn), :stOut => (dis = (:stExtOut, :stIntOut), capa = :StOut), :stSize => (dis = (:stLvl,), capa = :StSize))

	capa_sym = anyM.options.decomm != :none ? :oprCapa : :capa
	capaRestr_gdf = groupby(part.capaRestr,:cnstrType)

	# loop over groups of capacity restrictions (like out, stIn, ...)
	for restrGrp in capaRestr_gdf
		# relevant capacity variables
		type_sym = Symbol(restrGrp.cnstrType[1])
		info_ntup = cnstrType_dic[type_sym]

		allCns_arr = Array{DataFrame}(undef,size(restrGrp,1))

		# loop over indiviudal constraints
		for (idx,restr) in enumerate(eachrow(restrGrp))
			allCns_arr[idx] = createRestr(part,copy(part.var[Symbol(capa_sym,info_ntup.capa)]),restr,type_sym,info_ntup,ts_dic,r_dic,anyM.sets,anyM.supTs)
		end

		allCns_df = vcat(allCns_arr...)

		# add all constraints to part
		allCns_df[!,:cnsExpr] = map(x -> x.disp - x.capa,eachrow(allCns_df))
		cns_dic[Symbol(type_sym,:Restr)] = cnsCont(orderDf(allCns_df[!,[intCol(allCns_df)...,:cnsExpr]]),:smaller)
	end
end

# XXX sub-function to create restriction
function createRestr(part::TechPart, capaVar_df::DataFrame, restr::DataFrameRow, type_sym::Symbol, info_ntup::NamedTuple,
															ts_dic::Dict{Tuple{Int64,Int64},Array{Int64,1}}, r_dic::Dict{Tuple{Int64,Int64},Array{Int64,1}}, sets_dic::Dict{Symbol,Tree}, supTs_ntup::NamedTuple)

	conv_boo = type_sym in (:out,:in)
	dim_arr = conv_boo ? [:Ts_expSup,:Ts_dis,:R_dis,:Te] : [:Ts_expSup,:Ts_dis,:R_dis,:C,:Te]
	agg_arr = [:Ts_expSup,:Ts_dis,:R_dis] |> (x -> filter(x -> part.type == :emerging || x != :Ts_expSup,x))

	# get relevant carriers for conversion and storage variables
	relConv_arr = restr.car
	intC_arr = union(collect(part.actSt),map(y -> part.carrier[y],filter(x -> x in keys(part.carrier),[:stIntIn,:stIntOut])) |> (y -> isempty(y) ? Int[] : union(y...)))
	relSt_arr = filter(y -> isempty(sets_dic[:C].nodes[y].down) || y in intC_arr, [restr.car[1],getDescendants(restr.car[1],sets_dic[:C],true)...])

	# determines dimensions for aggregating dispatch variables
	capaVar_df[!,:lvlTs] .= restr.lvlTs
	capaVar_df[!,:lvlR] .= restr.lvlR

	# resize capacity variables (expect for stSize since these are already provided in energy units)
	if type_sym != :stSize
		capaVar_df[!,:var]  = capaVar_df[!,:var] .* map(x -> supTs_ntup.sca[(x.Ts_disSup,x.lvlTs)],	eachrow(capaVar_df[!,[:Ts_disSup,:lvlTs]]))
	end

	# replaces expansion with dispatch regions and aggregates capacity variables accordingy if required
	grpCapaVar_df = copy(select(capaVar_df,Not(:var))) |> (y -> unique(combine(x -> (R_dis = r_dic[(x.R_exp[1],x.lvlR[1])],),groupby(y,namesSym(y)))[!,Not([:R_exp,:lvlR])]))
	resExp_ntup = :Ts_expSup in agg_arr ? (Ts_expSup = part.balLvl.exp[1], Ts_disSup = supTs_ntup.lvl, R_dis = restr.lvlR, scr = 1) : (Ts_disSup = supTs_ntup.lvl, R_dis = restr.lvlR, scr = 1)
	sort!(grpCapaVar_df,sort(intCol(grpCapaVar_df)))
	grpCapaVar_df[!,:var] = aggUniVar(rename(capaVar_df,:R_exp => :R_dis),grpCapaVar_df,replace(agg_arr,:Ts_dis => :Ts_disSup),resExp_ntup,sets_dic)

	# expand capacity to dimension of dispatch
	capaDim_df = combine(x -> (Ts_dis = ts_dic[(x.Ts_disSup[1],x.lvlTs[1])],), groupby(grpCapaVar_df[!,Not(:var)],namesSym(grpCapaVar_df[!,Not(:var)])))[!,Not(:lvlTs)]
	sort!(capaDim_df,sort(intCol(capaDim_df)))
	select!(grpCapaVar_df,Not(:lvlTs))

	# obtain all relevant dispatch variables
	dispVar_arr = type_sym != :stSize ? intersect(keys(part.carrier),info_ntup.dis) : collect(info_ntup.dis)
	resDis_ntup = :Ts_expSup in agg_arr ? (Ts_expSup = part.balLvl.exp[1], Ts_dis = restr.lvlTs, R_dis = restr.lvlR) : (Ts_dis = restr.lvlTs, R_dis = restr.lvlR)
	for va in dispVar_arr
		# filter dispatch variables not belonging to relevant carrier
		if va in (:gen,:use)
			relC_arr = relConv_arr
		else
			relC_arr = relSt_arr
		end
		allVar_df = filter(r -> r.C in relC_arr, part.var[va])[!,Not(:Ts_disSup)]

		# get availablity (and in case of paramter of type out also efficiency since capacities refer to input capacity) parameter and add to dispatch variable
		ava_arr = matchSetParameter(allVar_df,part.par[Symbol(:ava,info_ntup.capa)],sets_dic, newCol = :ava)[!,:ava]
		if type_sym in (:out,:stOut)
			ava_arr = matchSetParameter(allVar_df,part.par[type_sym == :out ? :effConv : :effStOut],sets_dic,newCol = :eff)[!,:eff] .* ava_arr
		end
		allVar_df[!,:var] = allVar_df[!,:var] .* 1 ./ ava_arr
		# aggregate dispatch variables
		capaDim_df[!,va] = aggUniVar(allVar_df, select(capaDim_df,intCol(capaDim_df)), agg_arr, resDis_ntup, sets_dic)
	end
	# sum dispatch variables and filter cases without any
	capaDim_df[!,:disp] = map(x -> sum(x),eachrow(capaDim_df[!,dispVar_arr]))
	select!(capaDim_df,Not(dispVar_arr))
	capaDim_df = filter(x -> !(x.disp == AffExpr()),capaDim_df)

	# join capacity and dispatch variables to create final constraint
	grpCapaVar_df = combine(groupby(grpCapaVar_df,replace(dim_arr,:Ts_dis => :Ts_disSup)), :var => (x -> sum(x)) => :capa)
	cns_df = innerjoin(capaDim_df,grpCapaVar_df,on = intCol(grpCapaVar_df))
	return cns_df
end

# XXX create ratio constraints (Fix, Low and Up for use and gen)
function createRatioCns!(part::TechPart,cns_dic::Dict{Symbol,cnsCont},anyM::anyModel)

	# collects all tables for equations
	for type in intersect(keys(part.carrier),(:use, :gen)), limit in (:Up, :Low, :Fix)

		typeUp_sym = string(type) |> (x -> Symbol(uppercase(x[1]),x[2:end]))
		ratioName_sym = Symbol(:ratioEner,typeUp_sym,limit)
		if !(ratioName_sym in keys(part.par)) continue end

		# obtain variable name and parameter data
		cns_df = rename(copy(part.par[ratioName_sym].data),:val => :ratio)
		sort!(cns_df,sort(intCol(cns_df)))

		# joins parameter data with ratio controlled variable and all variables
		agg_arr = filter(r -> r != :Te && (part.type == :emerging || r != :Ts_expSup), intCol(cns_df))

		if part.type == :emerging
			srcRes_ntup = (anyM.sets[:Ts].nodes[cns_df[1,:Ts_dis]].lvl, anyM.sets[:R].nodes[cns_df[1,:R_dis]].lvl) |> (x -> (Ts_expSup = anyM.supTs.lvl, Ts_dis = x[1], R_dis = x[2]))
		else
			srcRes_ntup = (anyM.sets[:Ts].nodes[cns_df[1,:Ts_dis]].lvl, anyM.sets[:R].nodes[cns_df[1,:R_dis]].lvl) |> (x -> (Ts_dis = x[1], R_dis = x[2]))
		end

		if :M in namesSym(cns_df) # aggregates dispatch variables, if a mode is specified somewhere, mode dependant and non-mode dependant balances have to be aggregated separately
			# find cases where ratio constraint is mode dependant
			srcResM_ntup = (; zip(tuple(:M,keys(srcRes_ntup)...),tuple(1,values(srcRes_ntup)...))...)
			srcResNoM_ntup = (; zip(tuple(:M,keys(srcRes_ntup)...),tuple(0,values(srcRes_ntup)...))...)
			m_arr = findall(0 .!= cns_df[!,:M])
			noM_arr = setdiff(1:size(cns_df,1),m_arr)
			# aggregate variables with defined ratio
			cns_df[!,:ratioVar] .= AffExpr()
			cns_df[m_arr,:ratioVar] = aggUniVar(part.var[type], select(cns_df[m_arr,:],intCol(cns_df)), agg_arr, srcResM_ntup, anyM.sets)
			cns_df[noM_arr,:ratioVar] = aggUniVar(part.var[type], select(cns_df[noM_arr,:],intCol(cns_df)), agg_arr, srcResNoM_ntup, anyM.sets)
			# aggregate all variables
			cns_df[!,:allVar] .= AffExpr()
			cns_df[m_arr,:allVar] =	aggUniVar(part.var[type], select(cns_df[m_arr,:],intCol(cns_df)), filter(x -> x != :C,agg_arr), srcResM_ntup, anyM.sets)
			cns_df[noM_arr,:allVar] =	aggUniVar(part.var[type], select(cns_df[noM_arr,:],intCol(cns_df)), filter(x -> x != :C,agg_arr), srcResNoM_ntup, anyM.sets)
		else
			cns_df[!,:ratioVar] = aggUniVar(part.var[type], select(cns_df,intCol(cns_df)), agg_arr, srcRes_ntup, anyM.sets)
			cns_df[!,:allVar] =	aggUniVar(part.var[type], select(cns_df,intCol(cns_df)), filter(x -> x != :C,agg_arr), srcRes_ntup, anyM.sets)
		end

		# create corresponding constraint
		if occursin("Fix",string(limit))
			sign_sym = :equal
		elseif occursin("Low",string(limit))
			sign_sym = :greater
		else
			sign_sym = :smaller
		end

		cns_df[!,:cnsExpr] = map(x -> x.ratioVar - x.allVar * x.ratio, eachrow(cns_df))

		cns_dic[ratioName_sym] = cnsCont(orderDf(cns_df[!,[intCol(cns_df)...,:cnsExpr]]),sign_sym)
	end
end

# </editor-fold>

max_stars_repo_name: leonardgoeke/anyMOD

max_issues_repo_name: leonardgoeke/anyMOD

max_forks_repo_name: leonardgoeke/anyMOD


text:
using TBT
using Test

@testset "TBT.jl" begin
    # Write your tests here.
end

ext: jl

max_stars_repo_name: xiaodaigh/TBT

max_issues_repo_name: xiaodaigh/TBT

max_forks_repo_name: xiaodaigh/TBT


text:
module DataFrameHelpers

export isin

"""Function for comparing if items in a vector appear in another vector"""
function isin(y::Vector, x::Vector)
    output = Vector{BitVector}() # generate a vector for the bitwise outputs
    for i in 1:length(x) # loop through each item in the second vector
        push!(output, y .== x[i]) # elementwise check if the item in the second vector matches the first vector
    end
    # fairly certain I can do this with count. but I haven't checked
    int_to_vec(x) = x > 0 ? true : false # define a function to turn integer to true/false
    return int_to_vec.(sum(hcat(output))) # returns a bitvector
end

end # module

max_stars_repo_name: killerfridge/DataFrameHelpers

max_issues_repo_name: killerfridge/DataFrameHelpers

max_forks_repo_name: killerfridge/DataFrameHelpers


text:
type Sprite <: Drawable
    ptr::Ptr{Void}
    _texture::Texture

    function Sprite(ptr::Ptr{Void})
        s = new(ptr)
        finalizer(s, destroy)
        s
    end
end

function Sprite()
    Sprite(ccall((:sfSprite_create, libcsfml_graphics), Ptr{Void}, ()))
end

function Sprite(texture::Texture)
    s = Sprite()
    set_texture(s, texture)
    s
end

function copy(sprite::Sprite)
    return Sprite(ccall((:sfSprite_copy, libcsfml_graphics), Ptr{Void}, (Ptr{Void},), sprite.ptr))
end

function destroy(sprite::Sprite)
    ccall((:sfSprite_destroy, libcsfml_graphics), Void, (Ptr{Void},), sprite.ptr)
end

function set_position(sprite::Sprite, position::Vector2f)
    ccall((:sfSprite_setPosition, libcsfml_graphics), Void, (Ptr{Void}, Vector2f,), sprite.ptr, position)
end

function set_rotation(sprite::Sprite, angle::Real)
    ccall((:sfSprite_setRotation, libcsfml_graphics), Void, (Ptr{Void}, Cfloat,), sprite.ptr, angle)
end

function set_scale(sprite::Sprite, scale::Vector2f)
    ccall((:sfSprite_setScale, libcsfml_graphics), Void, (Ptr{Void}, Vector2f,), sprite.ptr, scale)
end

function set_origin(sprite::Sprite, origin::Vector2f)
    ccall((:sfSprite_setOrigin, libcsfml_graphics), Void, (Ptr{Void}, Vector2f,), sprite.ptr, origin)
end

function get_position(sprite::Sprite)
    return ccall((:sfSprite_getPosition, libcsfml_graphics), Vector2f, (Ptr{Void},), sprite.ptr)
end

function get_rotation(sprite::Sprite)
    return Real(ccall((:sfSprite_getRotation, libcsfml_graphics), Cfloat, (Ptr{Void},), sprite.ptr))
end

function get_origin(sprite::Sprite)
    return ccall((:sfSprite_getOrigin, libcsfml_graphics), Vector2f, (Ptr{Void},), sprite.ptr)
end

function get_scale(sprite::Sprite)
    ccall((:sfSprite_getScale, libcsfml_graphics), Vector2f, (Ptr{Void},), sprite.ptr)
end

function move(sprite::Sprite, offset::Vector2f)
    ccall((:sfSprite_move, libcsfml_graphics), Void, (Ptr{Void}, Vector2f,), sprite.ptr, offset)
end

function rotate(sprite::Sprite, angle::Real)
    ccall((:sfSprite_rotate, libcsfml_graphics), Void, (Ptr{Void}, Cfloat,), sprite.ptr, angle)
end

function scale(sprite::Sprite, factors::Vector2f)
    ccall((:sfSprite_scale, libcsfml_graphics), Void, (Ptr{Void}, Vector2f,), sprite.ptr, factors)
end

function set_texture(sprite::Sprite, texture::Texture)
    sprite._texture = texture
    ccall((:sfSprite_setTexture, libcsfml_graphics), Void, (Ptr{Void}, Ptr{Void}, Int32,), sprite.ptr, texture.ptr, 1)
end

function set_color(sprite::Sprite, color::Color)
    ccall((:sfSprite_setColor, libcsfml_graphics), Void, (Ptr{Void}, Color,), sprite.ptr, color)
end

function get_texture(sprite::Sprite)
    return Texture(ccall((:sfSprite_getTexture, libcsfml_graphics), Ptr{Void}, (Ptr{Void},), sprite.ptr))
end

function get_color(sprite::Sprite)
    return ccall((:sfSprite_getColor, libcsfml_graphics), Color, (Ptr{Void},), sprite.ptr)
end

function get_localbounds(sprite::Sprite)
    return ccall((:sfSprite_getLocalBounds, libcsfml_graphics), FloatRect, (Ptr{Void},), sprite.ptr)
end

function get_globalbounds(sprite::Sprite)
    return ccall((:sfSprite_getGlobalBounds, libcsfml_graphics), FloatRect, (Ptr{Void},), sprite.ptr)
end


text:
#=
exp_xform.jl
timing test
=#

using MIRT: exp_xform
using BenchmarkTools: @btime

# timing test
function exp_xform_time( ; T::DataType = ComplexF32)

    modes = (:element, :row, :column)

    N = 10000
    M = 100
    D = 20
    L = 10
    X = randn(T, N, L)
    U = randn(T, D, N)
    V = randn(T, D, M)

    # timing tests: :matrix is fastest, with :row a close 2nd
    for mode in (:matrix, modes...)
        @info "Case :$mode"
        tmp = (x, u, v) -> exp_xform(x, u, v ; mode=mode)
        @btime $tmp($X, $U, $V)
    end

#=
ir28 with 1.4.2:

:matrix   25.248 ms (8 allocations: 22.90 MiB)
:element 353.762 ms (4000001 allocations: 762.95 MiB)
:row      40.402 ms (1001 allocations: 22.99 MiB)
:column   51.424 ms (70001 allocations: 107.12 MiB)
=#

end

exp_xform_time()


text:
using CSV, DataFrames, SparseArrays, MKLSparse

export latency_benchmark

function latency_benchmark()
    kernel = "src/latency/kernel.jl"
    nwait = 1
    nworkers = 3
    niterations = 10
    timeout = 0.1
    nbytes = 100
    nrows = 100
    ncols = 100
    ncomponents = 3
    density = 0.1
    outputfile = tempname()
    mpiexec(cmd -> run(```
        $cmd -n $(nworkers+1) julia --project $kernel $outputfile
            --niterations $niterations
            --nbytes $nbytes
            --nrows $nrows
            --ncols $ncols
            --ncomponents $ncomponents
            --density $density
            --nwait $nwait
            --timeout $timeout
    ```))
    return df_from_latency_file(outputfile)
end

"""
    mymul!(C::Matrix, A::SparseMatrixCSC, B::Matrix)

Compute the sparse-dense matrix-matrix multiplication `A*B` and store the result in `C`.
Reference implementation for comparing performance scaling of different sparse multiplication rountines.
"""
function mymul!(C::Matrix, A::SparseMatrixCSC, B::Matrix)
    @boundscheck size(C, 1) == size(A, 1) || throw(DimensionMismatch("C has dimensions $(size(C)), but A has dimensions $(size(A))"))
    @boundscheck size(C, 2) == size(B, 2) || throw(DimensionMismatch("C has dimensions $(size(C)), but B has dimensions $(size(B))"))
    @boundscheck size(A, 2) == size(B, 1) || throw(DimensionMismatch("A has dimensions $(size(A)), but B has dimensions $(size(B))"))
    m, n = size(A)
    k = size(B, 2)
    rows = rowvals(A)
    vals = nonzeros(A)
    C .= 0
    for col = 1:n # column of A
        for i in nzrange(A, col)
            @inbounds row = rows[i] # row of A
            @inbounds val = vals[i] # A[row, col]
            @simd for j = 1:k # columns indices of B
                @inbounds C[row, j] += val * B[col, j]
            end
        end
    end
    C
end

sample!(V, W, X) = mymul!(W, X, V)

# function sample!(V, W, X)
#     mul!(W, X, V)
#     # mul!(V, X', W)
# end

function samples(X, V; nsamples::Integer)
    W = similar(V, size(X, 1), size(V, 2))
    timestamps = zeros(nsamples)
    latencies = zeros(nsamples)

    # force JIT
    t0 = time_ns()
    timestamps[1] = (time_ns() - t0) / 1e9
    latencies[1] = @elapsed sample!(V, W, X)

    # take samples
    t0 = time_ns()
    for i in 1:nsamples
        timestamps[i] = (time_ns() - t0) / 1e9
        latencies[i] = @elapsed sample!(V, W, X)
    end
    timestamps, latencies
end

function latency_sweep()

    ncols = 1000
    ncomponents = 3
    density = 0.05
    nsamples = 100

    df = DataFrame()
    V = randn(ncols, ncomponents)

    # dense matrices
    # nrows_all = [1000, 2000, 4000, 8000, 16000, 32000]
    # for nrows in nrows_all
    #     dfi = DataFrame()
    #     Xd = randn(nrows, ncols)
    #     timestamps, latencies = samples(Xd, V; nsamples)
    #     dfi.timestamp = timestamps
    #     dfi.latency = latencies
    #     dfi.nrows = nrows
    #     dfi.ncols = ncols
    #     dfi.ncomponents = ncomponents
    #     dfi.density = 1.0
    #     df = vcat(df, dfi, cols=:union)
    # end

    # sparse matrices
    nrows_all = [1000, 2000, 4000, 8000, 16000, 32000, 64000, 128000, 256000, 512000]
    for nrows in nrows_all
        dfi = DataFrame()
        Xs = sprand(nrows, ncols, density)
        timestamps, latencies = samples(Xs, V; nsamples)
        dfi.timestamp = timestamps
        dfi.latency = latencies
        dfi.nrows = nrows
        dfi.ncols = ncols
        dfi.ncomponents = ncomponents
        dfi.density = density
        df = vcat(df, dfi, cols=:union)
    end

    df.nflops = 2 .* df.nrows .* df.ncols .* df.ncomponents .* df.density
    df
end


text:
using JuMP

function homogeneous_system(model::JuMP.Model)
    A = full(JuMP.prepConstrMatrix(model))
    c, lb, ub = JuMP.prepProblemBounds(model)
    l, u = model.colLower, model.colUpper

    m, n = size(A)
    @assert m == length(lb) == length(ub)
    @assert model.nlpdata == nothing
    @assert isempty(model.quadconstr)
    @assert isempty(model.sosconstr)

    C = Array(Float64, 0, n)
    b = Float64[]
    for i in 1:m
        if !isinf(lb[i])
            C = [C, A[i,:]]
            push!(b, -lb[i])
        end
        if !isinf(ub[i])
            C = [C; -A[i,:]]
            push!(b, ub[i])
        end
    end

    F = eye(n, n)
    B = Array(Float64, 0, n)
    d = Float64[]
    for i in 1:n
        if !isinf(l[i])
            B = [B; F[i,:]]
            push!(d, -l[i])
        end
        if !isinf(u[i])
            B = [B; -F[i,:]]
            push!(d, u[i])
        end
    end

    return [C;B], [b;d]
end

function write_ine(model::JuMP.Model, fname::String)
    A, b = homogeneous_system(model)
    try
        write_ine(fname, convert(Matrix{Int},A), convert(Vector{Int},b))
    catch InexactError
        try
            write_ine(fname, convert(Matrix{Int},2A), convert(Vector{Int},2b))
        catch InexactError
            try
                write_ine(fname, convert(Matrix{Int},4A), convert(Vector{Int},4b))
            catch InexactError
                write_ine(fname, A, b)
            end
        end
    end
end

function get_extrema(model::JuMP.Model)
    A, b = homogeneous_system(model)
    dd = double_description([b A])

    vertices, rays = Vector{Float64}[], Vector{Float64}[]
    for cv in dd.R
        r = vec(cv)
        # if abs(r[1]) > 10eps()
        if abs(r[1] - 1) < ε
            push!(vertices, r[2:end])
            # push!(vertices, r[2:end] ./ r[1])
        else
            @assert abs(r[1]) < ε
            push!(rays, r[2:end] ./ minimum(r[find(r)]))
        end
    end
    return vertices, rays
end


text:
using StringBuilders

abstract type JLPAbstValue end

abstract type JLPAbstAsmFunc end

struct JLPMacro
    name::String
    args::JLPAbstList
    returnlist::JLPAbstList
end

struct JLPFn <: JLPAbstValue
    returntype::JLPAbstValue
    args::Vector{JLPAbstValue}
end

struct JLPStruct <: JLPAbstValue
    name::String
    params::Vector{JLPAbstValue}
end

struct JLPAnyValue <: JLPAbstValue
    name::String
end

struct JLPGenericStruct <: JLPAbstValue
    name::String
    typeparams::Vector{JLPAnyValue}
    params::Vector{Pair{String,JLPAbstValue}}
end

struct JLPPrimitive <: JLPAbstValue
    name::String
    bitlen::UInt
    option::String
    JLPPrimitive(name, bitlen) = new(name, bitlen, "")
    JLPPrimitive(name, bitlen, option) = new(name, bitlen, option)
end

struct JLPVectorType <: JLPAbstValue
    name::String
    T::JLPAbstValue
end

struct JLPFloats <: JLPAbstValue
    name::String
    bitlen::UInt
end

struct JLPAsmVar
    name::String
    ty::JLPAbstValue
    localvar::Bool
    JLPAsmVar(name::String, ty::JLPAbstValue, localvar::Bool) = new(name, ty, localvar)
    JLPAsmVar(name::String, ty::JLPAbstValue) = new(name, ty, true)
    JLPAsmVar(name::Integer, ty::JLPAbstValue) = new(string(name), ty)
end

struct JLPVector <: JLPAbstAsmFunc
    name::String
    T::JLPAbstValue
    length::JLPAsmVar
end

struct JLPAssignment
    target::JLPAsmVar
    val::JLPAbstAsmFunc
end

struct JLPAsmLabel
    label::String
    JLPAsmLabel(label::String) = new(label)
    JLPAsmLabel(label::Integer) = new("label" * string(label))
end

struct JLPAsmBr <: JLPAbstAsmFunc
    condition::Union{JLPAsmVar,Nothing}
    labels
    JLPAsmBr(condition, labels::Vector) = new(condition, labels)
    JLPAsmBr(label::JLPAsmLabel) = new(nothing, label)
end

struct JLPAsmPhi <: JLPAbstAsmFunc
    returntype::JLPAbstValue
    vallabel
end

struct JLPGetelementptr <:JLPAbstValue
    name::String
    fn::JLPFn
    boname::String
end

struct JLPAsmStore <: JLPAbstAsmFunc
    name::String
    fn::JLPFn
    boname::String
end

struct JLPAsmLoad <: JLPAbstAsmFunc
    name::String
    fn::JLPFn
    boname::String
end

struct JLPBinaryOp <: JLPAbstAsmFunc
    name::String
    fn::JLPFn
    boname::String
    op::String
end

struct JLPPrimNew <: JLPAbstAsmFunc
    type::JLPPrimitive
    val
end

struct JLPAsmRet <: JLPAbstAsmFunc
    val::JLPAsmVar
end

struct JLPAlloca <: JLPAbstAsmFunc
    ty::JLPAbstValue
end

struct JLPAsmFunc <: JLPAbstAsmFunc
    name::String
    asmname::String
    fn::JLPFn
end

struct JLPAsmCallFunc <: JLPAbstAsmFunc
    fn::JLPAbstAsmFunc
    args::Vector
end

abstract type JLPLLVMIR end

struct JLPLLVMFunc <: JLPLLVMIR
    name::String
    returntype
    args::Vector
    asms::Vector
end

struct JLPLLVMDef <: JLPLLVMIR
    name::String
end

mutable struct Count
    c::Int
    labelc::Int
end


text:
# ---------------------------------------------------------------------------- #
#
#   luteos.jl
#
#   High-order FEM Library in Julia
#
#   λυτέος
#   Fall 2017
#
#   Max Opgenoord
#
# ---------------------------------------------------------------------------- #

__precompile__()

"""
    luteos

Julia library for high-order FEM problems.
`luteos` is the Latin transliteration of the ancient Greek word λυτέος, which
means "one must solve".

Max Opgenoord

Fall 2017
"""
module luteos

import IterativeSolvers
import ILU
import Dierckx # this is only used for getPressureFunc.jl, perhaps there is a better solution for that.

export CDR, Elas, Material

# problem setup
include("general/problem.jl")
include("general/setup.jl")

# integration
export Master2D, Master3D
export PG1, PG2, PG3, PG4, PG5, PG6, PG7, PG8, PG9, PG10, PG11
include("integration/pgauss.jl")
export P1, P2, P3, P4
include("integration/porder.jl")

# mesh
export Mesh2D, Mesh3D
include("mesh/mesh.jl")
include("mesh/master.jl")
include("mesh/compJacob.jl")

# links to meshers
export writeMsh, readBAMG, runMesher, MesherFEFLOA, MesherBAMG, writeMsh
export readFEFLOA_2D, readFEFLOA_3D, writeMetric
include("meshers/mesher.jl")

# io
export writeTecplot, writeTecplotP1, convertToNodal
include("io/writeTecplotElas.jl")
include("io/writeTecplotCDR.jl")
include("io/writeTecplotMesh.jl")
include("io/writeTecplotP1.jl")
include("io/readSU2.jl")

# solution methods
export hdgSolve
include("solve/hdgSolveElas.jl")
include("solve/hdgSolveCDR.jl")

# utilities
export pressureFunc
include("util/getPressureFunc.jl")
include("util/getPressureFunc2p5D.jl")
export compMises2D, compMises3D
include("util/compMises.jl")

function __init__()

    setup()

end

end


text:
@testset "Robots: particle" begin
    # TODO: add tests
    q0 = nominal_configuration(particle)
    # visualizer
    vis = RoboDojo.Visualizer();
    @test visualize!(vis, particle, [q0], Δt=0.1);
end


text:
using CUDA

include("LineOfResponse.jl")
include("Util.jl")
include("Recon3D.jl")
include("Visualisation.jl")

function gpu_recon_test()
    DIMX = 512
    DIMY = 512
    DIMZ = 512

    image = ones(Float32, DIMX, DIMY, DIMZ)
    c_image = CUDA.CuArray(image)

    events, = read_3D("Triple_line_source.lmdT")
    c_events = CUDA.CuArray(events)

    c_corr = CUDA.CuArray(zeros(Float32, DIMX, DIMY, DIMZ))
    c_tmp_total_values = CUDA.CuArray(zeros(Float32, length(c_events)))

    @cuda blocks=length(c_events) threads=DIMX gpu_kernel(c_events, c_image, c_corr, c_tmp_total_values, DIMX, DIMY, DIMZ)

    c_image = c_image .* c_corr
    c_image = c_image ./ maximum(c_image)

	a = Array(c_image);
	save_image(a)

    visualise(c_image)
end

max_stars_repo_name: michielvangendt/PET-Julia

max_issues_repo_name: michielvangendt/PET-Julia

max_forks_repo_name: michielvangendt/PET-Julia


text:
# Uniform distribution

export Uniform

@parameterized Uniform() ≪ Lebesgue(𝕀)
@kwstruct Uniform()

# Standard Uniform

distproxy(::Uniform{()}) = Dists.Uniform()

logdensity(d::Uniform{()}, x) = 0.0

TV.as(::Uniform) = as𝕀

Base.rand(rng::AbstractRNG, T::Type, μ::Uniform{()}) = rand(rng, T)


text:
# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # #
#
# Description
#
#   The exponential atmospheric model.
#
# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # #
#
# References
#
#   [1] Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications.
#       Microcosm Press, Hawthorn, CA, USA.
#
# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # ==#

export expatmosphere

"""
    expatmosphere(h::Number)

Compute the atmospheric density [kg/m³] at the altitude `h` \\[m] (above the
ellipsoid) using the exponential atmospheric model:

                    ┌            ┐
                    │    h - h₀  │
    ρ(h) = ρ₀ ⋅ exp │ - ──────── │ ,
                    │      H     │
                    └            ┘

in which `ρ₀`, `h₀`, and `H` are parameters obtained from tables that depend
only on `h`.

"""
function expatmosphere(h::Number)
    # Check the bounds.
    (h < 0) && error("h must be positive.")

    # Transform `h` to km.
    h /= 1000

    # Get the values for the exponential model.
    id = (h >= 1000) ? 28 : findfirst( (x)->x > 0 , _expatmosphere_h₀ .- h ) - 1

    @inbounds begin
        h₀ = _expatmosphere_h₀[id]
        ρ₀ = _expatmosphere_ρ₀[id]
        H  = _expatmosphere_H[id]
    end

    # Compute the density.
    ρ₀ * exp( - ( h - h₀ )/H )
end

################################################################################
#                              Private Variables
################################################################################

"""
    _expatmosphere_h₀

Base altitude for the exponential atmospheric model [km].

"""
const _expatmosphere_h₀ = [0, 25, 30, 40, 50, 60, 70, 80, 90, 100, 110, 120,
                           130, 140, 150, 180, 200, 250, 300, 350, 400, 450,
                           500, 600, 700, 800, 900, 1000]

"""
    _expatmosphere_ρ₀

Nominal density for the exponential atmospheric model [kg/m³].

"""
const _expatmosphere_ρ₀ = [1.225, 3.899e-2, 1.774e-2, 3.972e-3, 1.057e-3,
                           3.206e-4, 8.770e-5, 1.905e-5, 3.396e-6, 5.297e-7,
                           9.661e-8, 2.438e-8, 8.484e-9, 3.845e-9, 2.070e-9,
                           5.464e-10, 2.789e-10, 7.248e-11, 2.418e-11,
                           9.518e-12, 3.725e-12, 1.585e-12, 6.967e-13,
                           1.454e-13, 3.614e-14, 1.170e-14, 5.245e-15,
                           3.019e-15]

"""
    _expatmosphere_H

Scale height for the exponential atmospheric model [km].

"""
const _expatmosphere_H = [7.249, 6.349, 6.682, 7.554, 8.382, 7.714, 6.549,
                          5.799, 5.382, 5.877, 7.263, 9.473, 12.636, 16.149,
                          22.523, 29.740, 37.105, 45.546, 53.628, 53.298,
                          58.515, 60.828, 63.822, 71.835, 88.667, 124.64,
                          181.05, 268.00]


text:
export BrickletMultiTouchV2TouchStateCallbackConfiguration
struct BrickletMultiTouchV2TouchStateCallbackConfiguration
    period::Integer
    value_has_to_change::Bool
end

export BrickletMultiTouchV2SPITFPErrorCount
struct BrickletMultiTouchV2SPITFPErrorCount
    error_count_ack_checksum::Integer
    error_count_message_checksum::Integer
    error_count_frame::Integer
    error_count_overflow::Integer
end

export BrickletMultiTouchV2Identity
struct BrickletMultiTouchV2Identity
    uid::String
    connected_uid::String
    position::Char
    hardware_version::Vector{Integer}
    firmware_version::Vector{Integer}
    device_identifier::Integer
end

export BrickletMultiTouchV2
"""
Capacitive touch sensor for 12 electrodes
"""
mutable struct BrickletMultiTouchV2 <: TinkerforgeDevice
    ipcon::IPConnection
    deviceInternal::PyObject

    """
    Creates an object with the unique device ID *uid* and adds it to
    the IP Connection *ipcon*.
    """
    function BrickletMultiTouchV2(uid::String, ipcon::IPConnection)
        package = pyimport("tinkerforge.bricklet_multi_touch_v2")
        deviceInternal = package.BrickletMultiTouchV2(uid, ipcon.ipconInternal)

        return new(ipcon, deviceInternal)
    end
end

export get_touch_state
"""
    $(SIGNATURES)

Returns the current touch state. The state is given as a array of
bools.

Element 0 to 11 represent the 12 electrodes and element 12 represents
the proximity.

If an electrode is touched, the corresponding element is *true*. If
a hand or similar is in proximity to the electrodes, element 12 is
*true*.

The proximity is activated with a distance of 1-2cm. An electrode
is already counted as touched if a finger is nearly touching the
electrode. This means that you can put a piece of paper or foil
or similar on top of a electrode to build a touch panel with
a professional look.

If you want to get the value periodically, it is recommended to use the
:cb:`Touch State` callback. You can set the callback configuration
with :func:`Set Touch State Callback Configuration`.
"""
function get_touch_state(device::BrickletMultiTouchV2)
    return device.deviceInternal.get_touch_state()
end

export set_touch_state_callback_configuration
"""
    $(SIGNATURES)

The period is the period with which the :cb:`Touch State` callback
is triggered periodically. A value of 0 turns the callback off.

If the `value has to change`-parameter is set to true, the callback is only
triggered after the value has changed. If the value didn't change within the
period, the callback is triggered immediately on change.

If it is set to false, the callback is continuously triggered with the period,
independent of the value.
"""
function set_touch_state_callback_configuration(device::BrickletMultiTouchV2, period, value_has_to_change)
    device.deviceInternal.set_touch_state_callback_configuration(period, value_has_to_change)
end

export get_touch_state_callback_configuration
"""
    $(SIGNATURES)

Returns the callback configuration as set by
:func:`Set Touch State Callback Configuration`.
"""
function get_touch_state_callback_configuration(device::BrickletMultiTouchV2)
    return device.deviceInternal.get_touch_state_callback_configuration()
end

export recalibrate
"""
    $(SIGNATURES)

Recalibrates the electrodes. Call this function whenever you changed
or moved you electrodes.
"""
function recalibrate(device::BrickletMultiTouchV2)
    device.deviceInternal.recalibrate()
end

export set_electrode_config
"""
    $(SIGNATURES)

Enables/disables electrodes with a bool array (see :func:`Get Touch State`).

*True* enables the electrode, *false* disables the electrode. A
disabled electrode will always return *false* as its state. If you
don't need all electrodes you can disable the electrodes that are
not needed.

It is recommended that you disable the proximity electrode (element 12) if
the proximity feature is not needed. This will reduce the amount of
traffic that is produced by the :cb:`Touch State` callback.

Disabling electrodes will also reduce power consumption.
"""
function set_electrode_config(device::BrickletMultiTouchV2, enabled_electrodes)
    device.deviceInternal.set_electrode_config(enabled_electrodes)
end

export get_electrode_config
"""
    $(SIGNATURES)

Returns the electrode configuration, as set by :func:`Set Electrode Config`.
"""
function get_electrode_config(device::BrickletMultiTouchV2)
    return device.deviceInternal.get_electrode_config()
end

export set_electrode_sensitivity
"""
    $(SIGNATURES)

Sets the sensitivity of the electrodes. An electrode with a high sensitivity
will register a touch earlier then an electrode with a low sensitivity.

If you build a big electrode you might need to decrease the sensitivity, since
the area that can be charged will get bigger. If you want to be able to
activate an electrode from further away you need to increase the sensitivity.

After a new sensitivity is set, you likely want to call :func:`Recalibrate`
to calibrate the electrodes with the newly defined sensitivity.
"""
function set_electrode_sensitivity(device::BrickletMultiTouchV2, sensitivity)
    device.deviceInternal.set_electrode_sensitivity(sensitivity)
end

export get_electrode_sensitivity
"""
    $(SIGNATURES)

Returns the current sensitivity, as set by :func:`Set Electrode Sensitivity`.
"""
function get_electrode_sensitivity(device::BrickletMultiTouchV2)
    return device.deviceInternal.get_electrode_sensitivity()
end

export set_touch_led_config
"""
    $(SIGNATURES)

Configures the touch LED to be either turned off, turned on, blink in
heartbeat mode or show the touch state (electrode touched = LED on).
"""
function set_touch_led_config(device::BrickletMultiTouchV2, config)
    device.deviceInternal.set_touch_led_config(config)
end

export get_touch_led_config
"""
    $(SIGNATURES)

Returns the LED configuration as set by :func:`Set Touch LED Config`
"""
function get_touch_led_config(device::BrickletMultiTouchV2)
    return device.deviceInternal.get_touch_led_config()
end

export get_spitfp_error_count
"""
    $(SIGNATURES)

Returns the error count for the communication between Brick and Bricklet.

The errors are divided into

* ACK checksum errors,
* message checksum errors,
* framing errors and
* overflow errors.

The errors counts are for errors that occur on the Bricklet side. All
Bricks have a similar function that returns the errors on the Brick side.
"""
function get_spitfp_error_count(device::BrickletMultiTouchV2)
    return device.deviceInternal.get_spitfp_error_count()
end

export set_bootloader_mode
"""
    $(SIGNATURES)

Sets the bootloader mode and returns the status after the requested
mode change was instigated.

You can change from bootloader mode to firmware mode and vice versa. A change
from bootloader mode to firmware mode will only take place if the entry function,
device identifier and CRC are present and correct.

This function is used by Brick Viewer during flashing. It should not be
necessary to call it in a normal user program.
"""
function set_bootloader_mode(device::BrickletMultiTouchV2, mode)
    return device.deviceInternal.set_bootloader_mode(mode)
end

export get_bootloader_mode
"""
    $(SIGNATURES)

Returns the current bootloader mode, see :func:`Set Bootloader Mode`.
"""
function get_bootloader_mode(device::BrickletMultiTouchV2)
    return device.deviceInternal.get_bootloader_mode()
end

export set_write_firmware_pointer
"""
    $(SIGNATURES)

Sets the firmware pointer for :func:`Write Firmware`. The pointer has
to be increased by chunks of size 64. The data is written to flash
every 4 chunks (which equals to one page of size 256).

This function is used by Brick Viewer during flashing. It should not be
necessary to call it in a normal user program.
"""
function set_write_firmware_pointer(device::BrickletMultiTouchV2, pointer)
    device.deviceInternal.set_write_firmware_pointer(pointer)
end

export write_firmware
"""
    $(SIGNATURES)

Writes 64 Bytes of firmware at the position as written by
:func:`Set Write Firmware Pointer` before. The firmware is written
to flash every 4 chunks.

You can only write firmware in bootloader mode.

This function is used by Brick Viewer during flashing. It should not be
necessary to call it in a normal user program.
"""
function write_firmware(device::BrickletMultiTouchV2, data)
    return device.deviceInternal.write_firmware(data)
end

export set_status_led_config
"""
    $(SIGNATURES)

Sets the status LED configuration. By default the LED shows
communication traffic between Brick and Bricklet, it flickers once
for every 10 received data packets.

You can also turn the LED permanently on/off or show a heartbeat.

If the Bricklet is in bootloader mode, the LED is will show heartbeat by default.
"""
function set_status_led_config(device::BrickletMultiTouchV2, config)
    device.deviceInternal.set_status_led_config(config)
end

export get_status_led_config
"""
    $(SIGNATURES)

Returns the configuration as set by :func:`Set Status LED Config`
"""
function get_status_led_config(device::BrickletMultiTouchV2)
    return device.deviceInternal.get_status_led_config()
end

export get_chip_temperature
"""
    $(SIGNATURES)

Returns the temperature as measured inside the microcontroller. The
value returned is not the ambient temperature!

The temperature is only proportional to the real temperature and it has bad
accuracy. Practically it is only useful as an indicator for
temperature changes.
"""
function get_chip_temperature(device::BrickletMultiTouchV2)
    return device.deviceInternal.get_chip_temperature()
end

export reset
"""
    $(SIGNATURES)

Calling this function will reset the Bricklet. All configurations
will be lost.

After a reset you have to create new device objects,
calling functions on the existing ones will result in
undefined behavior!
"""
function reset(device::BrickletMultiTouchV2)
    device.deviceInternal.reset()
end

export write_uid
"""
    $(SIGNATURES)

Writes a new UID into flash. If you want to set a new UID
you have to decode the Base58 encoded UID string into an
integer first.

We recommend that you use Brick Viewer to change the UID.
"""
function write_uid(device::BrickletMultiTouchV2, uid)
    device.deviceInternal.write_uid(uid)
end

export read_uid
"""
    $(SIGNATURES)

Returns the current UID as an integer. Encode as
Base58 to get the usual string version.
"""
function read_uid(device::BrickletMultiTouchV2)
    return device.deviceInternal.read_uid()
end

export get_identity
"""
    $(SIGNATURES)

Returns the UID, the UID where the Bricklet is connected to,
the position, the hardware and firmware version as well as the
device identifier.

The position can be 'a', 'b', 'c', 'd', 'e', 'f', 'g' or 'h' (Bricklet Port).
A Bricklet connected to an :ref:`Isolator Bricklet <isolator_bricklet>` is always at
position 'z'.

The device identifier numbers can be found :ref:`here <device_identifier>`.
|device_identifier_constant|
"""
function get_identity(device::BrickletMultiTouchV2)
    return device.deviceInternal.get_identity()
end

export register_callback
"""
Registers the given *function* with the given *callback_id*.
"""
function register_callback(device::BrickletMultiTouchV2, callback_id, function_)
    if isnothing(function_)
        device.registered_callbacks.pop(callback_id, None)
    else
        device.registered_callbacks[callback_id] = function_
    end
end


text:
using Statistics
using Distributions
using ProgressMeter
using KernelDensity
# using NumericalIntegration
using Trapz
import StatsBase: ordinalrank

#=
References
----------
    [1] Borgonovo, E. (2007). "A new uncertainty importance measure."
       Reliability Engineering & System Safety, 92(6):771-784,
       doi:10.1016/j.ress.2006.04.015.

    [2] Plischke, E., E. Borgonovo, and C. L. Smith (2013). "Global
       sensitivity measures from given data." European Journal of
       Operational Research, 226(3):536-550, doi:10.1016/j.ejor.2012.11.047.
=#

"""
    analyze(data::DeltaData, model_input::AbstractArray{<:Number, S1}, model_output::AbstractArray{<:Number, S2}; num_resamples::Int = 1_000, conf_level::Number = 0.95, progress_meter::Bool = true, N_override::Union{Nothing, Integer}=nothing) where S1 where S2

Performs a Delta Moment-Independent Analysis on the `model_output` produced with
the problem defined by the information in `data` and `model_input` and returns
a dictionary of results with the sensitivity indices and respective confidence
intervals for each of the parameters defined using the `num_resamples` and `conf_level`
keyword args.  The `progress_meter` keyword argument indicates whether a progress meter
will be displayed and defaults to true. The `N_override` keyword argument allows users
to override the `N` used in a specific `analyze` call to analyze just a subset
(useful for convergence graphs).
"""
function analyze(data::DeltaData, model_input::AbstractArray{<:Number, S1}, model_output::AbstractArray{<:Number, S2}; num_resamples::Int = 1_000, conf_level::Number = 0.95, progress_meter::Bool = true, N_override::Union{Nothing, Integer}=nothing) where S1 where S2

    # this method requires a confidence interval and num_resamples, so we do not
    # currently allow these to be set to Nothing as we did with the Sobol method ...
    # if a user wants to limit run time for intial tries they can set a very low
    # number of resamples
    # conf_flag = _check_conf_flag(num_resamples, conf_level)

    if size(model_output, 2) != 1
        error("Model output for analyzing DeltaData has more than one col, not handled yet.")
    else
        model_output = vec(model_output)
    end

    # define constants
    D = length(data.params) # number of uncertain parameters in problem
    # deal with overriding N
    if N_override === nothing
        N = data.N # number of samples
    else
        N_override > data.N ? error("N_override ($N_override) cannot be greater than original N used in sampling ($(data.N))") : nothing
        N = N_override # number of samples

        # reduce the output to just what should be considered for this N
        model_output = model_output[1:N]
        model_input = model_input[1:N,:]
    end

    exp =  (2 / (7 + tanh((1500 - N) / 500)))
    M = Int(min(ceil(N ^ exp), 48))
    m = collect(LinRange(0, N, M + 1))
    model_output_grid = collect(LinRange(minimum(model_output), maximum(model_output), 100)) # note that SALIb defaults the YGridlength to 100, while GlobalSensitivity.jl defaults to 2048

    # preallocate arrays
    delta_biased = Array{Float64}(undef, D)
    delta = Array{Float64}(undef, D)
    delta_conf = Array{Float64}(undef, D)
    firstorder = Array{Float64}(undef, D)
    firstorder_conf = Array{Float64}(undef, D)

    # set up progress meter
    counter = 0
    progress_meter ? p = Progress(D, counter, "Calculating indices for $D parameters ...") : nothing

    for i in 1:D
        # increment progress meter
        counter += 1
        progress_meter ? ProgressMeter.update!(p, counter) : nothing
        delta_biased[i], delta[i], delta_conf[i] = bias_reduced_delta(model_output, model_output_grid, model_input[:,i], m, num_resamples, conf_level)
        firstorder[i] = sobol_first(model_output, model_input[:,i], m)
        firstorder_conf[i] = sobol_first_conf(model_output, model_input[:,i], m, num_resamples, conf_level)
    end

    results = Dict(
        :delta_biased       => delta_biased,
        :delta              => delta,
        :delta_conf         => delta_conf,
        :firstorder         => firstorder,
        :firstorder_conf    => firstorder_conf
    )
    return results
end

"""
    calc_delta(model_output::AbstractArray{<:Number, S1}, model_output_grid::AbstractArray{<:Number, 1},
        model_input::AbstractArray{<:Number, S2}, m::AbstractArray{<:Number, 1}) where S1 where S2

    Plischke et al. 2013 estimator (eqn 26) for d_hat
"""
function calc_delta(model_output::AbstractArray{<:Number, S1}, model_output_grid::AbstractArray{<:Number, 1}, model_input::AbstractArray{<:Number, S2}, m::AbstractArray{<:Number, 1}) where S1 where S2
    N = length(model_output)
    k = KernelDensity.kde(model_output) # defaults are kernel = normal and bandwidth = Silverman which match SALib
    fy = Distributions.pdf(k, model_output_grid) # eq 23.1
    model_input_ranks = ordinalrank(model_input)

    d_hat = 0
    for j = 1:length(m) - 1
        ix = findall((model_input_ranks .> m[j]) .& (model_input_ranks .<= m[j + 1]))
        nm = length(ix)

        # get the separation between the total y pdf and the condition y pdf induced by this class

        # Estimated conditional distribution of y (using normal kernel)
        k = KernelDensity.kde(model_output[ix]) # defaults are kernel = normal and bandwidth = Silverman which match SALib
        fyc = Distributions.pdf(k, model_output_grid) # eq 23.2
        pdf_diff = abs.(fy - fyc) # eq 24

        # Use trapezoidal rule to estimate the difference between the curves.
        class_separation = Trapz.trapz(model_output_grid, pdf_diff) # eq 25

        # Increment Estimator
        d_hat += (nm / (2 * N)) * class_separation # eq 26

    end
    return d_hat

end

"""
    bias_reduced_delta(model_output::AbstractArray{<:Number, S1}, model_output_grid::AbstractArray{<:Number, 1},
        model_input::AbstractArray{<:Number, S2}, m::AbstractArray{<:Number, 1}, num_resamples::Int,
        conf_level::Number) where S1 where S2

Plischke et al. 2013 bias reduction technique (eqn 30)
"""
function bias_reduced_delta(model_output::AbstractArray{<:Number, S1}, model_output_grid::AbstractArray{<:Number, 1}, model_input::AbstractArray{<:Number, S2}, m::AbstractArray{<:Number, 1}, num_resamples::Int, conf_level::Number) where S1 where S2
    d = zeros(num_resamples)
    d_hat = calc_delta(model_output, model_output_grid, model_input, m)

    for i = 1:num_resamples
        r = rand(1:length(model_output), length(model_output))
        d[i] = calc_delta(model_output[r], model_output_grid, model_input[r], m)
    end

    d = 2 * d_hat .- d
    Z = quantile(Normal(0.0, 1.0),1 - (1 - conf_level)/2)

    return (d_hat, mean(d), Z * std(d))
end

"""
    sobol_first(model_output::AbstractArray{<:Number, S1}, model_input::AbstractArray{<:Number, S2},
        m::AbstractArray{<:Number, 1}) where S1 where S2
"""
function sobol_first(model_output::AbstractArray{<:Number, S1}, model_input::AbstractArray{<:Number, S2}, m::AbstractArray{<:Number, 1}) where S1 where S2
    model_input_ranks = ordinalrank(model_input)
    Vi = 0
    N = length(model_output)
    for j in 1:length(m) - 1
        ix = findall((model_input_ranks .> m[j]) .& (model_input_ranks .<= m[j + 1]))
        nm = length(ix)
        Vi += (nm / N) * ((mean(model_output[ix]) - mean(model_output))^2)
    end

    return Vi / var(model_output)

end

"""
    sobol_first_conf(model_output::AbstractArray{<:Number, S2}, model_input::AbstractArray{<:Number, S1},
        m::AbstractArray{<:Number, 1}, num_resamples::Int, conf_level::Number) where S1 where S2
"""
function sobol_first_conf(model_output::AbstractArray{<:Number, S2}, model_input::AbstractArray{<:Number, S1}, m::AbstractArray{<:Number, 1}, num_resamples::Int, conf_level::Number) where S1 where S2
    s = zeros(num_resamples)
    Z = quantile(Normal(0.0, 1.0),1 - (1 - conf_level)/2)

    for i in 1:num_resamples
        r = rand(1:length(model_output), length(model_output))
        s[i] = sobol_first(model_output[r], model_input[r], m)
    end

    return Z * std(s)
end

max_stars_repo_name: lrennels/SAJulia

max_issues_repo_name: lrennels/SAJulia


text:
function tracem(x::Array{T, 2}) where T <: AbstractFloat

  # We require the number of rows to be greater than the number of columns, so that m is greater than one.

  trans = false

  if size(x,1) < size(x,2)
    x = x'
    trans = true
  end

  n = size(x,2)
  m = round(Int,size(x,1)/n)

  y = zeros(m,1)  # We want this to be a 2-d array for subsequent matrix multiplication

  for i = 1:m
    y[i,1] = tr(x[(n*(i-1)+1):i*n,1:n])
  end

  if trans == true
    y = y'
  end

  return y

end


text:
#Computes the total Leonard-Jones potential energy

  function u_total_naive(p,data)
    ut = 0.
    for i in 1:(data.N)-1, j in i+1:(data.N)
      dist = distpbc(p[i],p[j],data.side)
      if dist < data.cutoff
        ut += upair(dist,data)
      end
    end
  return ut
end


text:
using AgnosticBayesEnsemble
using Documenter

DocMeta.setdocmeta!(AgnosticBayesEnsemble, :DocTestSetup, :(using AgnosticBayesEnsemble); recursive=true)

makedocs(;
    modules=[AgnosticBayesEnsemble],
    authors="Jules Rasetaharison <Jules.rasetaharison@auticon.de>",
    repo="https://github.com/hondoRandale/AgnosticBayesEnsemble.jl/hondoRandale/AgnosticBayesEnsemble.jl/blob/{commit}{path}#{line}",
    sitename="AgnosticBayesEnsemble.jl",
    format=Documenter.HTML(;
        prettyurls=get(ENV, "CI", "false") == "true",
        canonical="https://hondoRandale.github.io/AgnosticBayesEnsemble.jl",
        assets=String[],
    ),
    pages=[
        "Home" => "index.md",
    ],
)

deploydocs(;
    repo="github.com/hondoRandale/AgnosticBayesEnsemble.jl/hondoRandale/AgnosticBayesEnsemble.jl",
)


text:
# DESCRIPTION
# Wraps a phase signal between 0 and 360 deg.
#
# ------------------------------------------------------------------------------

function wrapper(input)
    phase=zeros(length(input),1)
    phase[:]=input[:]
    for i=1:length(phase)
        if phase[i]>0
            n=floor(phase[i]/360)
            phase[i]=phase[i]-n*360
        else
            n=ceil(abs(phase[i])/360)
            phase[i]=phase[i]+n*360
        end
    end
    return phase
end

max_stars_repo_name: umbertosaetti/J-GenHel

max_issues_repo_name: umbertosaetti/J-SimpleHel

max_forks_repo_name: umbertosaetti/J-SimpleHel


text:
# so need to write a function that uses the helper functions to send a get request to tm1 base url

@api_default function ping(api::TM1API; options...)
  tm1_get_json(api)
end


text:
module PowerLogistic

using Distributions, Roots, Plots

export powerLogistic, mdpr, mdor, plotmdor

# power for simple logistic regression based on
# http://personal.health.usf.edu/ywu/logistic.pdf
# Hsieh et al. STATISTICS IN MEDICINE
# Statist. Med. 17, 1623-1634 (1998)
#

# created on Sept. 24, 2012

# OR = exp(β) odds ratio
# p1 = Pr(Y = 1)
function SSizeLogisticCon(p1, OR; alpha=0.05, power=0.8)
   β = log(OR)
   za = quantile(Normal(), 1 - alpha/2)
   zb = quantile(Normal(), power)
   return ceil(Int64,(za+zb)^2 / (p1*(1-p1)*β^2))
end

# OR = exp(β) odds ratio
# p1 = Pr(Y = 1)
function powerLogisticCon(n, p1, OR; alpha=0.05)
   β = log(OR)
   za = quantile(Normal(),1-alpha/2)
   return ccdf(Normal(),za-sqrt(n*β^2*p1*(1-p1)))
end

###########################
# logistic regression logit(p) = beta0+ beta1*X
# B=pr(X=1)
# p1 = Pr(D|X=0) # event rate at X=0
# p2 = Pr(D|X=1) # event rate at X=1
# alpha - type I error rate
function SSizeLogisticBin(p1, p2, B; alpha=0.05, power=0.8)
   za = quantile(Normal(),1-alpha/2)
   zb = quantile(Normal(),power)

   p=(1-B)*p1+B*p2
   part1 = za*sqrt(p*(1-p)/B)
   part2 = zb*sqrt( p1*(1-p1)+p2*(1-p2)*(1-B)/B)
   part3 = (p1-p2)^2*(1-B)
   return ceil(Int64,(part1+part2)^2/part3)
end

###########################
# logistic regression logit(p) = beta0+ beta1*X
# B=pr(X=1)
# p1 = Pr(D|X=0) # event rate at X=0
# p2 = Pr(D|X=1) # event rate at X=1
# alpha - type I error rate
function powerLogisticBin(n, p1, p2, B; alpha=0.05)
   za = quantile(Normal(),1-alpha/2)

   p=(1-B)*p1+B*p2
   a = za*sqrt(p*(1-p)/B)
   b = sqrt( p1*(1-p1)+p2*(1-p2)*(1-B)/B)
   myc = (p1-p2)^2*(1-B)
   return cdf(Normal(), (sqrt(n*myc) - a)/b )
end

#############################################################
# public functions
#############################################################

"""
   powerLogistic(n = 0, p1 = .0, p2 = .0, B = .0, OR = .0, alpha = 0.05, power = 0.8)

Estimates power or sample size at type I error = `alpha` and power ≥ `power`.
All options are keyword options. If `n` is not specified (i.e., n = 0), it produces sample size.
Otherwise, it produces power of the sample. Other options are:

- `p1` = Pr(Y = 1) for a continous variable. `p1` = Pr(D|X=0), event rate at X = 0.
- `p2` = Pr(D|X=1), event rate at X = 1.
- `B` = Pr(X=1)
- `alpha`: Type I error (default = 0.05)
- `power`: Power (default = 0.8)

### Example
```jldoctest
julia> powerLogistic(p1=.5,OR=exp(0.405),power=0.95)
317

julia> powerLogistic(n = 317, p1=.5,OR=exp(0.405))
0.950...

julia> powerLogistic(p1 = .4, p2 = .5, B = .5, OR=exp(0.405), alpha = 0.05, power = 0.95)
1281

julia> powerLogistic(n = 1281, p1 = .4, p2 = .5, B = .5, OR=exp(0.405), alpha = 0.05)
0.950...
```
"""
function powerLogistic(;
   n = 0,
   p1 = 0.0,
   p2 = 0.0,
   B = 0.0,
   OR = 0.0,
   alpha = 0.05,
   power = 0.8)

   # return sample size
   if n == 0
      if B == 0.0 # continuous
         if p1 == 0.0 || OR == 0.0
            error("`p1` and `OR` must be greater than zero.")
         end
         return SSizeLogisticCon(p1,OR,alpha = alpha, power = power)
      else
         if p1 == p2 && p1 > 0.0 && p2 > 0.0
            error("`p1` and `p2` must be different and must be greater than zero.")
         end
         return SSizeLogisticBin(p1,p2,B,alpha = alpha, power = power)
      end
   # return power
   else
      if B == 0.0 # continuous
         if p1 == 0.0 || OR == 0.0
            error("`p1` and `OR` must be greater than zero.")
         end
         return powerLogisticCon(n, p1, OR,alpha = alpha)
      else
         if p1 == p2 && p1 > 0.0 && p2 > 0.0
            error("`p1` and `p2` must be different and must be greater than zero.")
         end
         return powerLogisticBin(n, p1, p2, B,alpha = alpha)
      end
   end
end

function mdpr(;n = 0, p1 = 0.0, B = 0.0, alpha = 0.05, power = 0.8)

   if n == 0
      error("`n` is required. `n` is the total sample size.")
   end

   if p1 == 0.0
      error("`p1` is required. `p1` is the event rate at X = 0.")
   end

   if B == 0.0
      error("`B` is required. `B` is the proportion of the sample with X = 1.")
   end

   # find p2 that achieves power greater than the power specified
   return fzero(x->powerLogisticBin(n,p1,x,B,alpha=alpha) - power,p1+1e-10,1.0-1e-10)
end

function mdor(;n = 0, p1 = 0.0, B = 0.0, alpha = 0.05, power = 0.8)

   if n == 0
      error("`n` is required. `n` is the total sample size.")
   end

   if p1 == 0.0
      error("`p1` is required. `p1` is the event rate at X = 0.")
   end

   if B == 0.0
      error("`B` is required. `B` is the proportion of the sample with X = 1.")
   end

   # find p2 that achieves power greater than the power specified
   p2 = fzero(x->powerLogisticBin(n,p1,x,B,alpha=alpha) - power,p1 + 1e-9,1.0-1e-9)
   return oddsratio(p2,p1)
end

function oddsratio(p1,p2)
   return (p1/(1-p1))/(p2/(1-p2))
end

function plotmdor(; nvec::Vector = [], p1::Union{Float64,Vector} = 0.0, B::Union{Float64,Vector} = 0.0, alpha = 0.05, power = 0.8)

   len = length(nvec)
   if len < 2
      error("`nvec` must have 2 or more values for sample size")
   end

   orvec = Array{Float64,2}(length(p1),len)

   for (i,p) in enumerate(p1)
      orvec[i,:] = [mdor(n = x, p1 = p, B = B, alpha = alpha, power = power) for x in nvec]
   end

   # plotly
   plotly()

   xlim_lower = nvec[len] / 40
   orvec = orvec'

   # plot
   nlines = length(p1)

   # plot
   if nlines == 1
      plot(nvec,orvec[:,1],
         marker=(2,.5,:circle,:blue),
         left_margin = 8mm,
         title = "Minimum Detectable Odds Ratio",
         ylabel = "Odds Ratio",
         xlabel = "Sample Size",
         ylims=(1-0.1,maximum(orvec)+0.1),
         xlims=(-xlim_lower,nvec[len]+xlim_lower),
         label=false,
         legend = false
      )
   else
      plot(nvec,orvec[:,1],label = string("p1 = ",p1[1]),marker=(2,.5,:circle))
      for i=2:(nlines-1)
         plot!(nvec,orvec[:,i],label = string("p1 = ",p1[i]),marker=(2,.5,:circle))
      end
      plot!(nvec,orvec[:,nlines],
            label = string("p1 = ",p1[nlines]),
            marker=(2,.5,:circle),
            left_margin = 8mm,
            title = "Minimum Detectable Odds Ratio",
            ylabel = "Odds Ratio",
            xlabel = "Sample Size",
            ylims=(1-0.1,maximum(orvec)+0.1),
            xlims=(-xlim_lower,nvec[len]+xlim_lower)
         )
    end
end

end


text:
#=
Functionalities to take a string corresponding to Julia code and evaluate
that code in a given module while capturing stdout and redirecting it to
a file.
=#

"""
$SIGNATURES

Consumes a string with Julia code, returns a vector of expression(s).

Note: this function was adapted from the `parse_input` function from Weave.jl.
"""
function parse_code(code::AS)
    exs = Any[] # Expr, Symbol or Any Julia core value
    n   = sizeof(code)
    pos = 1
    while pos ≤ n
        ex, pos = Meta.parse(code, pos)
        isnothing(ex) && continue
        push!(exs, ex)
    end
    exs
end

"""
$SIGNATURES

Run some code in a given module while redirecting stdout to a given path.
Return the result of the evaluation or `nothing` if the code was empty or
the evaluation failed.
If the evaluation errors, the error is printed to output then a warning is
shown.

## Arguments

1. `mod`:      the module in which to evaluate the code,
1. `code`:     string corresponding to the code,
1. `out_path`: path where stdout should be redirected

## Keywords

* `warn_err=true`:  whether to show a warning in the REPL if there was an error
                    running the code.
* `strip=false`:    whether to strip the code, this may already have been done.
"""
function run_code(mod::Module, code::AS, out_path::AS;
                  warn_err::Bool=true, strip_code::Bool=true)
    isempty(code) && return nothing
    strip_code && (code = strip(code))
    exs  = parse_code(strip(code))
    ne   = length(exs)
    res  = nothing # to capture final result
    err  = nothing
    ispath(out_path) || mkpath(dirname(out_path))
    open(out_path, "w") do outf
        if !FD_ENV[:SILENT_MODE]
            rprint("→ evaluating code [$(out_path |> basename |> splitext |> first)] in ($(locvar("fd_rpath")))")
        end
        redirect_stdout(outf) do
            e = 1
            while e <= ne
                try
                    res = Core.eval(mod, exs[e])
                catch e
                    io = IOBuffer()
                    showerror(io, e)
                    println(String(take!(io)))
                    err = typeof(e)
                    break
                end
                e += 1
            end
        end
    end
    # if there was an error, return nothing and possibly show warning
    if !isnothing(err)
        # TODO: add more informative message, maybe show type of error
        # + parent path
        FD_ENV[:SILENT_MODE] || print("\n")
        warn_err && @warn "There was an error of type $err running the code."
        res = nothing
    end
    # if last bit ends with `;` return nothing (no display)
    code[end] == ';' && return nothing
    # if last line is a Julia value return
    isa(exs[end], Expr) || return res
    # if last line of the code is a `show`
    if length(exs[end].args) > 1 && exs[end].args[1] == Symbol("@show")
        return nothing
    end
    # otherwise return the result of the last expression
    return res
end


text:
let
    fontsize = 13
    commomm_params = (;dpi = 1000,
        thickness_scaling = 1.3,
        xguidefontsize = fontsize,
        yguidefontsize = fontsize
    )

    # let
    #     p = plot(;title = "ider_colors", commomm_params...)
    #     for (label, color) in ider_colors
    #         scatter!(p, [0.0], [0.0]; color, label)
    #     end
    #     mysavefig(p, "ider_color")
    # end

    # let
    #     p = plot(;title = "source_markers", commomm_params...)
    #     srcs = [Kd, Nd, Fd] # cronological order
    #     for src in srcs
    #         marker = source_markers[src]
    #         label = source_labels[src]
    #         scatter!(p, [0.0], [0.0];
    #             marker, label, color = :black
    #         )
    #     end
    #     mysavefig(p, "source_markers")
    # end

    # let
    #     p = plot(;title = "subs_markers", commomm_params...)
    #     for (sub, m) in subs_markers
    #         label = string(sub)
    #         scatter!(p, [0.0], [0.0];
    #             m, label, color = :black
    #         )
    #     end
    #     mysavefig(p, "subs_markers")
    # end

    # let
    #     p = plot(;title = "method_colors", commomm_params...)
    #     for (method, m) in method_colors
    #         label = string(method)
    #         scatter!(p, [0.0], [0.0];
    #             m, label, color = :black
    #         )
    #     end
    #     mysavefig(p, "method_colors")
    # end

    # let
    #     p = plot(;title = "method_colors", commomm_params...)
    #     for (method, m) in method_colors
    #         label = string(method)
    #         scatter!(p, [0.0], [0.0];
    #             m, label, color = :black
    #         )
    #     end
    #     mysavefig(p, "method_colors")
    # end

    let
        exp_label = "\$ \\textbf{exp.} \$"
        me_label = "\$\\textbf{ME}^{\\langle \\textbf{1} \\rangle}\$"
        fba_label = "\$ \\textbf{FBA} \$"
        p = plot(;title = "method_comparizon", commomm_params...)
        for (label, marker) in [
            (exp_label, (:star, :black)),
            (fba_label, (:square, :red)),
            (me_label, (:circle, :blue)),
        ]
        scatter!(p, [0.0], [0.0]; label, marker,
            grid = false,
        )
        vspan!(p, [-1.0, 1.0]; label = "",
            linecolor = :darkgray, fillcolor = :darkgray,
            fillalpha = 0.6, linealpha = 0.6,
            xlim = [-1, 1],
            ylim = [-1, 1],
        )
        end

        mysavefig(p, "method_comparizon")
    end
end


text:
# declare a return type for the flux variability calculation -
mutable struct VLFluxVariabilityResult

    optimal_flux_array::Array{Float64,1}
    flux_variability_array::Array{Float64,1}

    # constructor -
    function VLFluxVariabilityResult()
        this = new()
    end
end

# declare a return type -
mutable struct VLOptimalFluxResult

    objective_value::Float64
    flux_array::Array{Float64,1}
    dual_array::Array{Float64,1}
    uptake_array::Array{Float64,1}
    exit_flag::Int64
    status_flag::Int64
    flux_bounds_array::Array{Float64,2}

    # constructor -
    function VLOptimalFluxResult()
        this = new()
    end
end

mutable struct Gene

    # information about a Gene -
    symbol::AbstractString
    organism_code::AbstractString

    function Gene()
        this = new()
    end
end

mutable struct Promoter

    # information about a Gene -
    symbol::AbstractString
    organism_code::AbstractString

    function Promoter()
        this = new()
    end
end