*! version 0.5.0 2014-08-21 | scott long | refactor with mata

//  matmath: math operations on matrices

// mata function for book keeping
mata

function mm_select(matsrc, matrows, matcols)
    {
        return(matsrc[matrows,matcols])
    }

function mm_sortmat(matsrc, sortcol)
    {
        return(sort(matsrc,sortcol))
    }

function mm_sortorder(matsrc, sortcol)
    {
        return(order(matsrc,sortcol))
    }
end

capture program drop _matmath
program _matmath, sclass

    version 11.2
    syntax, matout(string) matleft(string) ///
        type(string) [ matright(string) op(string)  ///
        list cols(string) rows(string) format(string) ///
        value(string) name(string) transpose ]

    tempname matrows matcols origreturns
    tempname newvec matinsert sortcol ordervec mat1 mat2
    _return hold `origreturns'
    _return restore `origreturns', hold

if "$debug"=="1" di "matout `matout'"

//  put left matrix into mata and get characteristics of the matrix

    local matleftrownms : rowfullnames `matleft', quoted
    local matleftcolnms : colfullnames `matleft', quoted
    local matleftrowsN = rowsof(`matleft')
    local matleftcolsN = colsof(`matleft')

    if "`type'"=="bi1" {

        mata: st_matrix(st_local("matout"), ///
              st_matrix(st_local("matleft")) ///
        :`op' st_matrix(st_local("matright")) )

        matrix colnames `matout' = `matleftcolnms'
        matrix rownames `matout' = `matleftrownms'
    }

    else if "`type'"=="bi2" {

        if "`op'"=="%of" { // what % 1st is of 2nd
            mata: st_matrix(st_local("matout"), ///
             100*(st_matrix(st_local("matleft")) ///
               :/ st_matrix(st_local("matright")) ) )
        }

        else if "`op'"=="%dif" { // % change 1st to 2nd
            mata:  st_matrix(st_local("matout"), ///
             100*( st_matrix(st_local("matleft")) ///
                :- st_matrix(st_local("matright")) ) ///
                :/ st_matrix(st_local("matleft")) ///
               )
        }

        matrix colnames `matout' = `matleftcolnms'
        matrix rownames `matout' = `matleftrownms'

    }

    else if "`type'"=="bi3" {

        if "`op'"=="before" ///
            mat_capp `matout' : `matleft' `matright', miss(.)
        else if "`op'"=="above" ///
            mat_rapp `matout' : `matleft' `matright', miss(.)

    }

   else if "`type'"=="un1" {

        mata: st_matrix(st_local("matout"), ///
         `op'(st_matrix(st_local("matleft")) ) )

        matrix rownames `matout' = `matleftrownms'
        matrix colnames `matout' = `matleftcolnms'
   }

   else if "`type'"=="un2" {

        if "`op'"=="sq" {
            mata: st_matrix(st_local("matout"), ///
                  st_matrix(st_local("matleft")) ///
               :* st_matrix(st_local("matleft")) )
        }
        else if "`op'"=="inv" {
            mata: st_matrix(st_local("matout"), ///
             1  :/ st_matrix(st_local("matleft")) )
        }
        else if "`op'"=="logit" {
            mata:         st_matrix(st_local("matout"), ///
               log (      st_matrix(st_local("matleft")) ///
                 :/ ( 1:- st_matrix(st_local("matleft")) ) ///
               )  )
        }
        else if "`op'"=="odds" {
            mata:         st_matrix(st_local("matout"), ///
                          st_matrix(st_local("matleft")) ///
                 :/ ( 1:- st_matrix(st_local("matleft")) ) ///
               )
        }

    }

    else if "`type'"=="un3" {

        if ("`op'"=="diag") matrix `matout' = `matleft'

        if "`op'"=="fill" {
            matrix `matout' = J(`matleftrowsN',`matleftcolsN',`value')
            matrix colnames `matout' = `matleftcolnms'
            matrix rownames `matout' = `matleftrownms'
        }
    }

   else if "`type'"=="book" {

        mata: `matleft' = st_matrix(st_local("matleft"))

        //  select rows/cols to keep

        if "`op'"=="keep" {

            matrix `matcols' = J(1,1,.)
            matrix `matrows' = J(1,1,.)

            if "`rows'"!="" {
                foreach r in `rows' {
                    if (`r'>`matleftrowsN') | (`r'<1) {
                        di as error "row `r' in keep() is invalid"
                        local iserror = 1
                    }
                }
                local tmp = subinstr("`rows'"," ","\",.)
                matrix `matrows' = (`tmp')
            }
            if "`cols'"!="" {
                foreach c in `cols' {
                    if (`c'>`matleftcolsN') | (`c'<1) {
                        di as error "column `c' in keep() is invalid"
                        local iserror = 1
                    }
                }
                local tmp = subinstr("`cols'"," ","\",.)
                matrix `matcols' = (`tmp')
            }
            if "`iserror'"=="1" {
                _return restore `origreturns'
                exit
            }

            mata: st_matrix(st_local("matout"), ///
                mm_select(st_matrix(st_local("matleft")), ///
                       st_matrix(st_local("matrows")), ///
                       st_matrix(st_local("matcols")) ) )
        } // keep

        //  select rows/cols to drop

        else if "`op'"=="drop" {

            local tmpindices ""
            forvalues row = 1/`matleftrowsN' {
                local ipos : list posof "`row'" in rows
                if (`ipos'==0) local tmpindices "`tmpindices' \ `row'"
            }
            local tmpindices = subinstr("`tmpindices'"," \ ","",1)
            matrix `matrows' = (`tmpindices')
            local rows = subinstr("`tmpindices'","\"," ",.) // used later

            local tmpindices ""
            forvalues col = 1/`matleftcolsN' {
                local ipos : list posof "`col'" in cols
                if (`ipos'==0) local tmpindices "`tmpindices' , `col'"
            }
            local tmpindices = subinstr("`tmpindices'"," , ","",1)
            matrix `matcols' = (`tmpindices')
            local cols = subinstr("`tmpindices'",","," ",.)

            mata: st_matrix(st_local("matout"), ///
                mm_select(st_matrix(st_local("matleft")), ///
                       st_matrix(st_local("matrows")), ///
                       st_matrix(st_local("matcols")) ) )

        } // drop

        else if "`op'"=="insert" {

            if "`rows'"!="" { // insert row

                matrix `newvec' = J(1,`matleftcolsN',`value')
                if ("`name'"!="") matrix rowname `newvec' = `name'

                if (`rows'== 0) matrix `matout' = `newvec' \ `matleft'
                else if `rows'==`matleftrowsN' ///
                    matrix `matout' = `matleft' \ `newvec'
                else {
                    local st1 = `rows' + 1
                    matrix `mat1' = `matleft'[1..`rows',1...]
                    matrix `mat2' = `matleft'[`st1'...,1...]
                    matrix `matout' = `mat1' \ `newvec' \ `mat2'
                }

            } // insert row

            if "`cols'"!="" { // insert column

                if "`rows'"!="" {
                    matrix `matinsert' = `matout' // in case rows were inserted
                    local matleftrowsN = rowsof(`matinsert')
                }
                else matrix `matinsert' = `matleft'
                matrix `newvec' = J(`matleftrowsN',1,`value')
                if ("`name'"!="") matrix colname `newvec' = `name'

                if (`cols'==0) matrix `matout' = `newvec' , `matinsert'
                else if `cols'==`matleftcolsN' ///
                    matrix `matout' = `matinsert' , `newvec'
                else {
                    local st1 = `cols' + 1
                    matrix `mat1' = `matinsert'[1...,1..`cols']
                    matrix `mat2' = `matinsert'[1...,`st1'...]
                    matrix `matout' = `mat1' , `newvec' , `mat2'
                }

            } // insert col

        } // insert

        else if "`op'"=="sort" {

            scalar `sortcol' = `rows'`cols'

            if "`rows'"!="" local sortbyrow = 1
            else local sortbyrow = 0

            matrix `matout' = `matleft'
            if `sortbyrow' matrix `matout' = `matout''

            mata: st_matrix(st_local("ordervec"), ///
                  mm_sortorder(st_matrix(st_local("matout")), ///
                               st_numscalar(st_local("sortcol")) ) )

            mata: st_matrix(st_local("matout"), ///
                  mm_sortmat(st_matrix(st_local("matout")), ///
                             st_numscalar(st_local("sortcol")) ) )

            * rearrange labels
            local orderrowsN = rowsof(`ordervec')
            forvalues i = 1/`orderrowsN' {
                local irow = `ordervec'[`i',1]
                if `sortbyrow' local rnm : word `irow' of `matleftcolnms'
                else local rnm : word `irow' of `matleftrownms'
                local newnms `" `newnms' `"`rnm'"' "'
            }
            if `sortbyrow' {
                matrix `matout' = `matout''
                matrix rownames `matout' = `matleftrownms'
                matrix colnames `matout' = `newnms'
            }
            else {
                matrix colnames `matout' = `matleftcolnms'
                matrix rownames `matout' = `newnms '
            }

        } // sort

        //  reattach labels for drop and keep

        if "`op'"=="keep" | "`op'"=="drop" {

            if "`rows'"!="" {
                foreach row in `rows' {
                    local rnm : word `row' of `matleftrownms'
                    * `"` '"' needed since names have ""s in them
                    local newrownms `" `newrownms' `"`rnm'"' "'
                }
                matrix rownames `matout' = `newrownms'
            }
            else matrix rownames `matout' = `matleftrownms'

            if "`cols'"!="" { // selected columns
                foreach col in `cols' {
                    local cnm : word `col' of `matleftcolnms'
                    local newcolnms `" `newcolnms' `"`cnm'"' "'
                }
                matrix colnames `matout' = `newcolnms'
            }
            else matrix colnames `matout' = `matleftcolnms'
        }

    } // book

// close up

    if ("`transpose'"=="transpose") matrix `matout' = `matout''

    if ("`list'"=="list") ///
        matlist `matout', title(matrix `matout') format(`format')

end
exit
* version 0.4.0 2014-08-18 | scott long | sort