program define epimodels_util
	version 16.0
	`0'
end

program define check_total_population
    version 16.0
	syntax anything

	if (`anything' < 1.00+`c(epsdouble)') {
		display ""
	    display as result "Warning! Total population size (`anything') is less than 1.00"
		display as result "Make sure the initial conditions at t0 represent population groups counts in persons."
	}
end

program define check_steps
    version 16.0
	syntax anything

	capture confirm integer number `anything'
	if (_rc) {
		display as error "Number of simulation steps in each simulation day must be an integer number."
		error 198
	}

	capture assert (`anything' > 0)
	if (_rc) {
	    display as error "Number of simulation steps in each simulation day must be 1 or more."
		error 198
	}

	capture assert (`anything' <= 1e4)
	if (_rc) {
	    display as error "Number of simulation steps in each simulation day must be no more than 10000."
		error 198
	}
end


program define check_days
	version 16.0
	syntax anything

	capture confirm integer number `anything'
	if (_rc) {
		display as error "Number of simulation steps (days) must an integer number."
		error 198
	}

	capture assert (`anything' > 0)
	if (_rc) {
	    display as error "Number of simulation steps (days) must be 1 or more."
		error 198
	}
end

program define check_day0_date
	version 16.0
	syntax [anything]

	if (`"`anything'"'!="") {
		if (date(`"`anything'"',"YMD")==.) {
			display as error "Option day0() is specified incorrectly."
			display as error "The date must be specified in the YYYY-MM-DD format, for example: 2020-02-29"
			error 111
		}
	}
end

program define makedatevar

	version 16.0
	syntax varname, [day0(string)] datefmt(string)

	if (`"`day0'"'!="") {
	  quietly replace t = t + date(`"`day0'"',"YMD")  // only affects the data, not the matrix
	  format t `datefmt'
	  label variable `varlist' "Date"
	}

end

program define ditable, rclass

    version 16.0

	syntax varname, ///
	    days(real) [day0(string)] datefmt(string) ///
	    [modeltitle(string)] mcolnames(string) ivar(string) ///
		varlabels(string) ylabel(string) [digits(real 2) comma(string) percent]

	quietly maxinfect `ivar' `varlist', ///
	  day0(`"`day0'"') datefmt("`datefmt'") `percent'

	local max=r(maxinfect)
	local dstar=r(d_maxinfect)
	local tstar=r(t_maxinfect)
	local ostar=r(o_maxinfect)
	local X=floor(t[`ostar'])
	local last=_N

	local title_t0 = "t0"
	local title_tX = "t`X'"
	local title_t1 = "t`days'"

	if (`"`day0'"' != "") {
	    local title_t0 = string(`varlist'[1], "`datefmt'")
		local title_tX = string(`varlist'[`=`dstar'+1'], "`datefmt'")
		local title_t1 = string(`varlist'[`last'], "`datefmt'")
	}

	local twid = 17
	local cwid = 14
	local lmarg = 2

	forvalues i = 2/`:word count `varlabels'' {
	    local vl `"`:word `i' of `varlabels''"'
		local twid=max(strlen(`"`vl'"')+1,`twid')
	}

	local twidth = `twid' + 3*(`cwid'+1)
	local titlepos = `lmarg' + floor((`twidth'-strlen(`"`modeltitle'"'))/2)
	if (`titlepos' < `lmarg') local titlepos=`lmarg'
	local titleoffset=`"_col(`=`titlepos'+1')"'

	display ""
    display `titleoffset' "`modeltitle'"
    tempname tab
	.`tab' = ._tab.new, col(4) lmarg(`lmarg') commas
	.`tab'.width    `twid' | `cwid' `cwid' `cwid'
	.`tab'.titlefmt .   %`cwid's   %`cwid's %`cwid's
	.`tab'.numfmt   .   %`cwid'.`digits'f`comma'  %`cwid'.`digits'f`comma' %`cwid'.`digits'f`comma'
	.`tab'.sep, top
	.`tab'.titles `"`ylabel'"' "`title_t0'" "`title_tX'" "`title_t1'"
	.`tab'.sep, mid

	local total_t0 = 0
	local total_tX = 0
	local total_t1 = 0

	forvalues i = 2/`:word count `varlabels'' {
	    local vl `"`:word `i' of `varlabels''"'
		local vn `:word `i' of `mcolnames''
	    .`tab'.row `"`vl'"' `=`vn'[1]' `=`vn'[`=`dstar'+1']' `=`vn'[`last']'
		local total_t0 = `total_t0' + `=`vn'[1]'
		local total_tX = `total_tX' + `=`vn'[`=`dstar'+1']'
		local total_t1 = `total_t1' + `=`vn'[`last']'
	}
	.`tab'.sep, middle
	.`tab'.row `"Total"' `total_t0' `total_tX' `total_t1'
	.`tab'.sep, bottom

	if (`"`day0'"'!="") {
	  local specdate `"(`=string(`=`varlist'[`ostar']',`"`datefmt'"')') "'
	  local specday =floor(t[`dstar']-t[1]+1)
	}
	else {
	  local specday =floor(t[`ostar'])
	}

	local maxtext `"`=string(`=`max'',"%20.4gc")'"'
	if (`"`percent'"'!="") local maxtext `"`maxtext'%"'
	display as text "The maximum size of the infected group {result:`maxtext'} is reached on day {result:`specday' `specdate'}of the simulation."

	if (`ostar'==`last') {
	    display as result "Warning! The peak is detected at the end of the simulation, and may change if you extend it."
	}
	display ""

	return scalar maxinfect= `max'
	return scalar t_maxinfect=`tstar'
	return scalar d_maxinfect=`dstar'
	return scalar o_maxinfect=`ostar'
end

program define maxinfect, rclass

    version 16.0

	syntax varlist, [day0(string) datefmt(string) percent]

	local ivar `"`: word 1 of `varlist''"'
	local tvar `"`: word 2 of `varlist''"'
	if (`"`datefmt'"'=="") local datefmt="%td"

	tempname tmax max dmax

	summarize `ivar' if (!missing(`ivar')), meanonly
	scalar `max' = r(max)

	summarize `tvar' if (abs(`ivar'-`max') < c(epsdouble)), meanonly
	scalar `tmax' = r(min) // in case of multiple identical values pick the first one
	quietly count if (`tvar' < `tmax')
    scalar `dmax' = `=r(N)'

	return scalar o_maxinfect=`dmax'+1
	return scalar t_maxinfect=`tmax'
	return scalar d_maxinfect=`dmax'
	return scalar maxinfect=`max'
end

program define pdfreport

	version 16.0

	syntax [varlist(default=none)], modelname(string) modelparams(string) ///
	                [modelgraph(string)] [appendixgraphs(string)] ///
					save(string)

	if (strlen(`"`modelparams'"')>40) local br `"`=char(10)'"'

	mata epimodels_about()

	putpdf begin
	putpdf paragraph
	putpdf text ("EPIMODELS Report"), bold font("Helvetica",28,steelblue)
	putpdf paragraph
	putpdf text ("This report was generated on `c(current_date)'"), font("Helvetica",16,steelblue)
	putpdf paragraph
	putpdf text ("`modelname' model with parameters:`br'`modelparams'."), bold

	capture findfile "epimodels_eq_`modelname'.png"
	if !_rc {
	    // model with equations
		putpdf paragraph
		putpdf text ("          ")
		putpdf image "`=r(fn)'" , width(4)    // todo: instream this??
    }
	if (`"`modelgraph'"'!="") {
		putpdf paragraph
		putpdf image `"`modelgraph'"'
	}

	putpdf paragraph
	putpdf text ("Generated with EPIMODELS version `epimodels_version' from `compile_date' built for Stata v`compile_version'"), font("Helvetica",12,steelblue)
	putpdf paragraph
	putpdf text ("For more information, visit EPIMODELS' homepage: http://www.radyakin.org/stata/epimodels/"), font("Helvetica",12,steelblue)
	putpdf pagebreak
	putpdf paragraph , halign(center)
	putpdf text ("Simulation results"), font("Helvetica",28,steelblue)

	if (`"`varlist'"'!="") {
		putpdf table t=data(`varlist'), varnames
		putpdf table t(.,.), halign(center) valign(center) nformat(%12.0fc) font("Consolas",8)
		putpdf table t(1,.), bgcolor("aliceblue")
		putpdf paragraph
	}

	foreach v in `varlist' {
		putpdf text ("`v': `:variable label `v''`=char(10)'")
	}

	// optional appendix
	if (`"`appendixgraphs'"'!="") {
	    putpdf sectionbreak, landscape
		local first=1
	    foreach f in `appendixgraphs' {
		    if (!`first') {
			    putpdf pagebreak
				local first=0
		    }
			putpdf paragraph
			putpdf image `"`f'"'
		}
	}

	putpdf save `"`save'"' , replace

end

program define popmatrix, rclass

	version 16.0

	syntax , agevar(varname) sexvar(varname) at(numlist) ///
	         [malecode(integer 1) femalecode(integer 2)]

	assert !missing(`agevar')
	assert !missing(`sexvar')
	assert (`malecode'!=`femalecode')
	local cc=subinstr(`"`at'"',","," ",.)
	local NG `: word count `cc''

	tempvar agegrp sexgrp
	quietly egen `agegrp'=cut(`agevar'), at(`at')
	quietly recode `sexvar' (`malecode'=1) (`femalecode'=2) , generate(`sexgrp')

	local N=_N
	local rn ""
	tempname M
	matrix `M' = J(`=`NG'-1',2,.)
	forval i=1/`=`NG'-1' {
	  local v1 =`: word `i' of `cc''
      local v2 =`: word `=`i'+1' of `cc''-1
	  local rn `"`rn' "`v1'-`v2'""'
      forval s=1/2 {
		quietly count if `agegrp'==`v1' & `sexgrp'==`s'
		local n=r(N)
		matrix `M'[`i',`s']=`n'/`N'
	  }
	}
	matrix colnames `M' = "Males" "Females"
	matrix rownames `M' = `rn'
	return matrix F=`M'
	return local N=_N
end

// END OF FILE