r2rs-stats 0.1.1

use crate::{
    optim::getListElement,
    sexprec::{SEXP, SEXPREC, SEXPTYPE},
};
use ::libc;
extern "C" {
    #[no_mangle]
    fn memset(_: *mut libc::c_void, _: libc::c_int, _: libc::c_ulong) -> *mut libc::c_void;

    #[no_mangle]
    fn sqrt(_: libc::c_double) -> libc::c_double;
    /* IEEE -Inf */
    #[no_mangle]
    static mut R_NaReal: libc::c_double;
    /*
     *  R : A Computer Language for Statistical Data Analysis
     *  Copyright (C) 1998-2005   The R Core Team
     *
     *  This header file is free software; you can redistribute it and/or modify
     *  it under the terms of the GNU Lesser General Public License as published by
     *  the Free Software Foundation; either version 2.1 of the License, or
     *  (at your option) any later version.
     *
     *  This file is part of R. R is distributed under the terms of the
     *  GNU General Public License, either Version 2, June 1991 or Version 3,
     *  June 2007. See doc/COPYRIGHTS for details of the copyright status of R.
     *
     *  This program is distributed in the hope that it will be useful,
     *  but WITHOUT ANY WARRANTY; without even the implied warranty of
     *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
     *  GNU Lesser General Public License for more details.
     *
     *  You should have received a copy of the GNU Lesser General Public License
     *  along with this program; if not, a copy is available at
     *  https://www.R-project.org/Licenses/
     */
    /* Included by R.h: API */
    #[no_mangle]
    fn error(_: *const libc::c_char, _: ...) -> !;
    #[no_mangle]
    fn R_alloc(_: size_t, _: libc::c_int) -> *mut libc::c_char;
    #[no_mangle]
    fn TYPEOF(x: SEXP) -> libc::c_int;
    /* Vector Access Functions */
    #[no_mangle]
    fn LENGTH(x: SEXP) -> libc::c_int;
    #[no_mangle]
    fn INTEGER(x: SEXP) -> *mut libc::c_int;
    #[no_mangle]
    fn REAL(x: SEXP) -> *mut libc::c_double;
    #[no_mangle]
    fn SET_STRING_ELT(x: SEXP, i: R_xlen_t, v: SEXP);
    #[no_mangle]
    fn SET_VECTOR_ELT(x: SEXP, i: R_xlen_t, v: SEXP) -> SEXP;
    /* An empty environment at the root of the
    environment tree */
    #[no_mangle]
    static mut R_BaseEnv: SEXP;
    /* Current srcref, for debuggers */
    /* Special Values */
    #[no_mangle]
    static mut R_NilValue: SEXP;
    /* "name" */
    #[no_mangle]
    static mut R_NamesSymbol: SEXP;
    #[no_mangle]
    fn coerceVector(_: SEXP, _: SEXPTYPE) -> SEXP;
    #[no_mangle]
    fn asLogical(x: SEXP) -> libc::c_int;
    #[no_mangle]
    fn asInteger(x: SEXP) -> libc::c_int;
    #[no_mangle]
    fn asReal(x: SEXP) -> libc::c_double;
    #[no_mangle]
    fn allocMatrix(_: SEXPTYPE, _: libc::c_int, _: libc::c_int) -> SEXP;
    #[no_mangle]
    fn duplicate(_: SEXP) -> SEXP;
    #[no_mangle]
    fn eval(_: SEXP, _: SEXP) -> SEXP;
    #[no_mangle]
    fn install(_: *const libc::c_char) -> SEXP;
    #[no_mangle]
    fn mkChar(_: *const libc::c_char) -> SEXP;
    #[no_mangle]
    fn setAttrib(_: SEXP, _: SEXP, _: SEXP) -> SEXP;
    /* Defining NO_RINLINEDFUNS disables use to simulate platforms where
    this is not available */
    /* need remapped names here for use with R_NO_REMAP */
    /*
       These are the inlinable functions that are provided in Rinlinedfuns.h
       It is *essential* that these do not appear in any other header file,
       with or without the Rf_ prefix.
    */
    #[no_mangle]
    fn allocVector(_: SEXPTYPE, _: R_xlen_t) -> SEXP;
    #[no_mangle]
    fn lang4(_: SEXP, _: SEXP, _: SEXP, _: SEXP) -> SEXP;
    #[no_mangle]
    fn ScalarReal(_: libc::c_double) -> SEXP;
    #[no_mangle]
    fn protect(_: SEXP) -> SEXP;
    #[no_mangle]
    fn unprotect(_: libc::c_int);
    #[no_mangle]
    fn dcgettext(
        __domainname: *const libc::c_char,
        __msgid: *const libc::c_char,
        __category: libc::c_int,
    ) -> *mut libc::c_char;
    /*
     *  R : A Computer Language for Statistical Data Analysis
     *  Copyright (C) 2012   The R Core Team.
     *
     *  This program is free software; you can redistribute it and/or modify
     *  it under the terms of the GNU General Public License as published by
     *  the Free Software Foundation; either version 2 of the License, or
     *  (at your option) any later version.
     *
     *  This program is distributed in the hope that it will be useful,
     *  but WITHOUT ANY WARRANTY; without even the implied warranty of
     *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
     *  GNU General Public License for more details.
     *
     *  You should have received a copy of the GNU General Public License
     *  along with this program; if not, a copy is available at
     *  https://www.R-project.org/Licenses/
     */
    /* auxiliary */

}
pub type size_t = libc::c_ulong;
/*
 *  R : A Computer Language for Statistical Data Analysis
 *  Copyright (C) 2000, 2001 The R Core Team.
 *
 *  This header file is free software; you can redistribute it and/or modify
 *  it under the terms of the GNU Lesser General Public License as published by
 *  the Free Software Foundation; either version 2.1 of the License, or
 *  (at your option) any later version.
 *
 *  This file is part of R. R is distributed under the terms of the
 *  GNU General Public License, either Version 2, June 1991 or Version 3,
 *  June 2007. See doc/COPYRIGHTS for details of the copyright status of R.
 *
 *  This program is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU Lesser General Public License for more details.
 *
 *  You should have received a copy of the GNU Lesser General Public License
 *  along with this program; if not, a copy is available at
 *  https://www.R-project.org/Licenses/
 */
/* Included by R.h: API */
pub type Rboolean = libc::c_uint;
pub const TRUE: Rboolean = 1;
pub const FALSE: Rboolean = 0;
pub type ptrdiff_t = libc::c_long;
/*, MAYBE */
/* both config.h and Rconfig.h set SIZEOF_SIZE_T, but Rconfig.h is
skipped if config.h has already been included. */
pub type R_xlen_t = ptrdiff_t;
/* Fundamental Data Types:  These are largely Lisp
 * influenced structures, with the exception of LGLSXP,
 * INTSXP, REALSXP, CPLXSXP and STRSXP which are the
 * element types for S-like data objects.
 *
 *   --> TypeTable[] in ../main/util.c for  typeof()
 */
/* UUID identifying the internals version -- packages using compiled
code should be re-installed when this changes */
/*  These exact numeric values are seldom used, but they are, e.g., in
 *  ../main/subassign.c, and they are serialized.
*/
/* NOT YET using enum:
 *  1) The SEXPREC struct below has 'SEXPTYPE type : 5'
 * (making FUNSXP and CLOSXP equivalent in there),
 * giving (-Wall only ?) warnings all over the place
 * 2) Many switch(type) { case ... } statements need a final `default:'
 * added in order to avoid warnings like [e.g. l.170 of ../main/util.c]
 *   "enumeration value `FUNSXP' not handled in switch"
 */
/*
 *  R : A Computer Language for Statistical Data Analysis
 *  Copyright (C) 2002-2016   The R Core Team.
 *
 *  This program is free software; you can redistribute it and/or modify
 *  it under the terms of the GNU General Public License as published by
 *  the Free Software Foundation; either version 2 of the License, or
 *  (at your option) any later version.
 *
 *  This program is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU General Public License for more details.
 *
 *  You should have received a copy of the GNU General Public License
 *  along with this program; if not, a copy is available at
 *  https://www.R-project.org/Licenses/
 */
// for abs
/*
  KalmanLike, internal to StructTS:
  .Call(C_KalmanLike, y, mod$Z, mod$a, mod$P, mod$T, mod$V, mod$h, mod$Pn,
        nit, FALSE, update)
  KalmanRun:
  .Call(C_KalmanLike, y, mod$Z, mod$a, mod$P, mod$T, mod$V, mod$h, mod$Pn,
        nit, TRUE, update)
*/
/* y vector length n of observations
  Z vector length p for observation equation y_t = Za_t +  eps_t
  a vector length p of initial state
  P p x p matrix for initial state uncertainty (contemparaneous)
  T  p x p transition matrix
  V  p x p = RQR'
  h = var(eps_t)
  anew used for a[t|t-1]
  Pnew used for P[t|t -1]
  M used for M = P[t|t -1]Z

  op is FALSE for KalmanLike, TRUE for KalmanRun.
  The latter computes residuals and states and has
  a more elaborate return value.

  Almost no checking here!
*/
#[no_mangle]
pub unsafe extern "C" fn KalmanLike(
    mut sy: SEXP,
    mut mod_0: SEXP,
    mut sUP: SEXP,
    mut op: SEXP,
    mut update: SEXP,
) -> SEXP {
    let mut lop: libc::c_int = asLogical(op);
    mod_0 = protect(duplicate(mod_0));
    let mut sZ: SEXP = getListElement(
        mod_0,
        b"Z\x00" as *const u8 as *const libc::c_char as *mut libc::c_char,
    );
    let mut sa: SEXP = getListElement(
        mod_0,
        b"a\x00" as *const u8 as *const libc::c_char as *mut libc::c_char,
    );
    let mut sP: SEXP = getListElement(
        mod_0,
        b"P\x00" as *const u8 as *const libc::c_char as *mut libc::c_char,
    );
    let mut sT: SEXP = getListElement(
        mod_0,
        b"T\x00" as *const u8 as *const libc::c_char as *mut libc::c_char,
    );
    let mut sV: SEXP = getListElement(
        mod_0,
        b"V\x00" as *const u8 as *const libc::c_char as *mut libc::c_char,
    );
    let mut sh: SEXP = getListElement(
        mod_0,
        b"h\x00" as *const u8 as *const libc::c_char as *mut libc::c_char,
    );
    let mut sPn: SEXP = getListElement(
        mod_0,
        b"Pn\x00" as *const u8 as *const libc::c_char as *mut libc::c_char,
    );
    if TYPEOF(sy) != 14 as libc::c_int
        || TYPEOF(sZ) != 14 as libc::c_int
        || TYPEOF(sa) != 14 as libc::c_int
        || TYPEOF(sP) != 14 as libc::c_int
        || TYPEOF(sPn) != 14 as libc::c_int
        || TYPEOF(sT) != 14 as libc::c_int
        || TYPEOF(sV) != 14 as libc::c_int
    {
        error(dcgettext(
            b"stats\x00" as *const u8 as *const libc::c_char,
            b"invalid argument type\x00" as *const u8 as *const libc::c_char,
            5 as libc::c_int,
        ));
    }
    let mut n: libc::c_int = LENGTH(sy);
    let mut p: libc::c_int = LENGTH(sa);
    let mut y: *mut libc::c_double = REAL(sy);
    let mut Z: *mut libc::c_double = REAL(sZ);
    let mut T: *mut libc::c_double = REAL(sT);
    let mut V: *mut libc::c_double = REAL(sV);
    let mut P: *mut libc::c_double = REAL(sP);
    let mut a: *mut libc::c_double = REAL(sa);
    let mut Pnew: *mut libc::c_double = REAL(sPn);
    let mut h: libc::c_double = asReal(sh);
    let mut anew: *mut libc::c_double = R_alloc(
        p as size_t,
        ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
    ) as *mut libc::c_double;
    let mut M: *mut libc::c_double = R_alloc(
        p as size_t,
        ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
    ) as *mut libc::c_double;
    let mut mm: *mut libc::c_double = R_alloc(
        (p * p) as size_t,
        ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
    ) as *mut libc::c_double;
    // These are only used if(lop), but avoid -Wall trouble
    let mut ans: SEXP = R_NilValue;
    let mut resid: SEXP = R_NilValue;
    let mut states: SEXP = R_NilValue;
    if lop != 0 {
        ans = allocVector(19 as libc::c_int as SEXPTYPE, 3 as libc::c_int as R_xlen_t);
        protect(ans);
        resid = allocVector(14 as libc::c_int as SEXPTYPE, n as R_xlen_t);
        SET_VECTOR_ELT(ans, 1 as libc::c_int as R_xlen_t, resid);
        states = allocMatrix(14 as libc::c_int as SEXPTYPE, n, p);
        SET_VECTOR_ELT(ans, 2 as libc::c_int as R_xlen_t, states);
        let mut nm: SEXP = protect(allocVector(
            16 as libc::c_int as SEXPTYPE,
            3 as libc::c_int as R_xlen_t,
        ));
        SET_STRING_ELT(
            nm,
            0 as libc::c_int as R_xlen_t,
            mkChar(b"values\x00" as *const u8 as *const libc::c_char),
        );
        SET_STRING_ELT(
            nm,
            1 as libc::c_int as R_xlen_t,
            mkChar(b"resid\x00" as *const u8 as *const libc::c_char),
        );
        SET_STRING_ELT(
            nm,
            2 as libc::c_int as R_xlen_t,
            mkChar(b"states\x00" as *const u8 as *const libc::c_char),
        );
        setAttrib(ans, R_NamesSymbol, nm);
        unprotect(1 as libc::c_int);
    }
    let mut sumlog: libc::c_double = 0.0f64;
    let mut ssq: libc::c_double = 0.0f64;
    let mut nu: libc::c_int = 0 as libc::c_int;
    let mut l: libc::c_int = 0 as libc::c_int;
    while l < n {
        let mut i: libc::c_int = 0 as libc::c_int;
        while i < p {
            let mut tmp: libc::c_double = 0.0f64;
            let mut k: libc::c_int = 0 as libc::c_int;
            while k < p {
                tmp += *T.offset((i + p * k) as isize) * *a.offset(k as isize);
                k += 1
            }
            *anew.offset(i as isize) = tmp;
            i += 1
        }
        if l > asInteger(sUP) {
            let mut i_0: libc::c_int = 0 as libc::c_int;
            while i_0 < p {
                let mut j: libc::c_int = 0 as libc::c_int;
                while j < p {
                    let mut tmp_0: libc::c_double = 0.0f64;
                    let mut k_0: libc::c_int = 0 as libc::c_int;
                    while k_0 < p {
                        tmp_0 +=
                            *T.offset((i_0 + p * k_0) as isize) * *P.offset((k_0 + p * j) as isize);
                        k_0 += 1
                    }
                    *mm.offset((i_0 + p * j) as isize) = tmp_0;
                    j += 1
                }
                i_0 += 1
            }
            let mut i_1: libc::c_int = 0 as libc::c_int;
            while i_1 < p {
                let mut j_0: libc::c_int = 0 as libc::c_int;
                while j_0 < p {
                    let mut tmp_1: libc::c_double = *V.offset((i_1 + p * j_0) as isize);
                    let mut k_1: libc::c_int = 0 as libc::c_int;
                    while k_1 < p {
                        tmp_1 += *mm.offset((i_1 + p * k_1) as isize)
                            * *T.offset((j_0 + p * k_1) as isize);
                        k_1 += 1
                    }
                    *Pnew.offset((i_1 + p * j_0) as isize) = tmp_1;
                    j_0 += 1
                }
                i_1 += 1
            }
        }
        if !((*y.offset(l as isize)).is_nan() as i32 != 0 as libc::c_int) {
            nu += 1;
            let mut rr: *mut libc::c_double = 0 as *mut libc::c_double;
            if lop != 0 {
                rr = REAL(resid)
            }
            let mut resid0: libc::c_double = *y.offset(l as isize);
            let mut i_2: libc::c_int = 0 as libc::c_int;
            while i_2 < p {
                resid0 -= *Z.offset(i_2 as isize) * *anew.offset(i_2 as isize);
                i_2 += 1
            }
            let mut gain: libc::c_double = h;
            let mut i_3: libc::c_int = 0 as libc::c_int;
            while i_3 < p {
                let mut tmp_2: libc::c_double = 0.0f64;
                let mut j_1: libc::c_int = 0 as libc::c_int;
                while j_1 < p {
                    tmp_2 += *Pnew.offset((i_3 + j_1 * p) as isize) * *Z.offset(j_1 as isize);
                    j_1 += 1
                }
                *M.offset(i_3 as isize) = tmp_2;
                gain += *Z.offset(i_3 as isize) * *M.offset(i_3 as isize);
                i_3 += 1
            }
            ssq += resid0 * resid0 / gain;
            if lop != 0 {
                *rr.offset(l as isize) = resid0 / sqrt(gain)
            }
            sumlog += gain.ln();
            let mut i_4: libc::c_int = 0 as libc::c_int;
            while i_4 < p {
                *a.offset(i_4 as isize) =
                    *anew.offset(i_4 as isize) + *M.offset(i_4 as isize) * resid0 / gain;
                i_4 += 1
            }
            let mut i_5: libc::c_int = 0 as libc::c_int;
            while i_5 < p {
                let mut j_2: libc::c_int = 0 as libc::c_int;
                while j_2 < p {
                    *P.offset((i_5 + j_2 * p) as isize) = *Pnew.offset((i_5 + j_2 * p) as isize)
                        - *M.offset(i_5 as isize) * *M.offset(j_2 as isize) / gain;
                    j_2 += 1
                }
                i_5 += 1
            }
        /* -Wall */
        } else {
            let mut rr_0: *mut libc::c_double = 0 as *mut libc::c_double;
            if lop != 0 {
                rr_0 = REAL(resid)
            }
            let mut i_6: libc::c_int = 0 as libc::c_int;
            while i_6 < p {
                *a.offset(i_6 as isize) = *anew.offset(i_6 as isize);
                i_6 += 1
            }
            let mut i_7: libc::c_int = 0 as libc::c_int;
            while i_7 < p * p {
                *P.offset(i_7 as isize) = *Pnew.offset(i_7 as isize);
                i_7 += 1
            }
            if lop != 0 {
                *rr_0.offset(l as isize) = R_NaReal
            }
            /* -Wall */
        }
        if lop != 0 {
            let mut rs: *mut libc::c_double = REAL(states);
            let mut j_3: libc::c_int = 0 as libc::c_int;
            while j_3 < p {
                *rs.offset((l + n * j_3) as isize) = *a.offset(j_3 as isize);
                j_3 += 1
            }
        }
        l += 1
    }
    let mut res: SEXP = protect(allocVector(
        14 as libc::c_int as SEXPTYPE,
        2 as libc::c_int as R_xlen_t,
    ));
    *REAL(res).offset(0 as libc::c_int as isize) = ssq / nu as libc::c_double;
    *REAL(res).offset(1 as libc::c_int as isize) = sumlog / nu as libc::c_double;
    if lop != 0 {
        SET_VECTOR_ELT(ans, 0 as libc::c_int as R_xlen_t, res);
        if asLogical(update) != 0 {
            setAttrib(
                ans,
                install(b"mod\x00" as *const u8 as *const libc::c_char),
                mod_0,
            );
        }
        unprotect(3 as libc::c_int);
        return ans;
    } else {
        if asLogical(update) != 0 {
            setAttrib(
                res,
                install(b"mod\x00" as *const u8 as *const libc::c_char),
                mod_0,
            );
        }
        unprotect(2 as libc::c_int);
        return res;
    };
}
#[no_mangle]
pub unsafe extern "C" fn KalmanSmooth(mut sy: SEXP, mut mod_0: SEXP, mut sUP: SEXP) -> SEXP {
    let mut sZ: SEXP = getListElement(
        mod_0,
        b"Z\x00" as *const u8 as *const libc::c_char as *mut libc::c_char,
    );
    let mut sa: SEXP = getListElement(
        mod_0,
        b"a\x00" as *const u8 as *const libc::c_char as *mut libc::c_char,
    );
    let mut sP: SEXP = getListElement(
        mod_0,
        b"P\x00" as *const u8 as *const libc::c_char as *mut libc::c_char,
    );
    let mut sT: SEXP = getListElement(
        mod_0,
        b"T\x00" as *const u8 as *const libc::c_char as *mut libc::c_char,
    );
    let mut sV: SEXP = getListElement(
        mod_0,
        b"V\x00" as *const u8 as *const libc::c_char as *mut libc::c_char,
    );
    let mut sh: SEXP = getListElement(
        mod_0,
        b"h\x00" as *const u8 as *const libc::c_char as *mut libc::c_char,
    );
    let mut sPn: SEXP = getListElement(
        mod_0,
        b"Pn\x00" as *const u8 as *const libc::c_char as *mut libc::c_char,
    );
    if TYPEOF(sy) != 14 as libc::c_int
        || TYPEOF(sZ) != 14 as libc::c_int
        || TYPEOF(sa) != 14 as libc::c_int
        || TYPEOF(sP) != 14 as libc::c_int
        || TYPEOF(sT) != 14 as libc::c_int
        || TYPEOF(sV) != 14 as libc::c_int
    {
        error(dcgettext(
            b"stats\x00" as *const u8 as *const libc::c_char,
            b"invalid argument type\x00" as *const u8 as *const libc::c_char,
            5 as libc::c_int,
        ));
    }
    let mut ssa: SEXP = 0 as *mut SEXPREC;
    let mut ssP: SEXP = 0 as *mut SEXPREC;
    let mut ssPn: SEXP = 0 as *mut SEXPREC;
    let mut res: SEXP = 0 as *mut SEXPREC;
    let mut states: SEXP = R_NilValue;
    let mut sN: SEXP = 0 as *mut SEXPREC;
    let mut n: libc::c_int = LENGTH(sy);
    let mut p: libc::c_int = LENGTH(sa);
    let mut y: *mut libc::c_double = REAL(sy);
    let mut Z: *mut libc::c_double = REAL(sZ);
    let mut a: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut P: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut T: *mut libc::c_double = REAL(sT);
    let mut V: *mut libc::c_double = REAL(sV);
    let mut h: libc::c_double = asReal(sh);
    let mut Pnew: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut at: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut rt: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut Pt: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut gains: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut resids: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut Mt: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut L: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut gn: libc::c_double = 0.;
    let mut Nt: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut var: Rboolean = TRUE;
    ssa = duplicate(sa);
    protect(ssa);
    a = REAL(ssa);
    ssP = duplicate(sP);
    protect(ssP);
    P = REAL(ssP);
    ssPn = duplicate(sPn);
    protect(ssPn);
    Pnew = REAL(ssPn);
    res = allocVector(19 as libc::c_int as SEXPTYPE, 2 as libc::c_int as R_xlen_t);
    protect(res);
    let mut nm: SEXP = protect(allocVector(
        16 as libc::c_int as SEXPTYPE,
        2 as libc::c_int as R_xlen_t,
    ));
    SET_STRING_ELT(
        nm,
        0 as libc::c_int as R_xlen_t,
        mkChar(b"smooth\x00" as *const u8 as *const libc::c_char),
    );
    SET_STRING_ELT(
        nm,
        1 as libc::c_int as R_xlen_t,
        mkChar(b"var\x00" as *const u8 as *const libc::c_char),
    );
    setAttrib(res, R_NamesSymbol, nm);
    unprotect(1 as libc::c_int);
    states = allocMatrix(14 as libc::c_int as SEXPTYPE, n, p);
    SET_VECTOR_ELT(res, 0 as libc::c_int as R_xlen_t, states);
    at = REAL(states);
    sN = allocVector(14 as libc::c_int as SEXPTYPE, (n * p * p) as R_xlen_t);
    SET_VECTOR_ELT(res, 1 as libc::c_int as R_xlen_t, sN);
    Nt = REAL(sN);
    let mut anew: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut mm: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut M: *mut libc::c_double = 0 as *mut libc::c_double;
    anew = R_alloc(
        p as size_t,
        ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
    ) as *mut libc::c_double;
    M = R_alloc(
        p as size_t,
        ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
    ) as *mut libc::c_double;
    mm = R_alloc(
        (p * p) as size_t,
        ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
    ) as *mut libc::c_double;
    Pt = R_alloc(
        (n * p * p) as size_t,
        ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
    ) as *mut libc::c_double;
    gains = R_alloc(
        n as size_t,
        ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
    ) as *mut libc::c_double;
    resids = R_alloc(
        n as size_t,
        ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
    ) as *mut libc::c_double;
    Mt = R_alloc(
        (n * p) as size_t,
        ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
    ) as *mut libc::c_double;
    L = R_alloc(
        (p * p) as size_t,
        ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
    ) as *mut libc::c_double;
    let mut l: libc::c_int = 0 as libc::c_int;
    while l < n {
        let mut i: libc::c_int = 0 as libc::c_int;
        while i < p {
            let mut tmp: libc::c_double = 0.0f64;
            let mut k: libc::c_int = 0 as libc::c_int;
            while k < p {
                tmp += *T.offset((i + p * k) as isize) * *a.offset(k as isize);
                k += 1
            }
            *anew.offset(i as isize) = tmp;
            i += 1
        }
        if l > asInteger(sUP) {
            let mut i_0: libc::c_int = 0 as libc::c_int;
            while i_0 < p {
                let mut j: libc::c_int = 0 as libc::c_int;
                while j < p {
                    let mut tmp_0: libc::c_double = 0.0f64;
                    let mut k_0: libc::c_int = 0 as libc::c_int;
                    while k_0 < p {
                        tmp_0 +=
                            *T.offset((i_0 + p * k_0) as isize) * *P.offset((k_0 + p * j) as isize);
                        k_0 += 1
                    }
                    *mm.offset((i_0 + p * j) as isize) = tmp_0;
                    j += 1
                }
                i_0 += 1
            }
            let mut i_1: libc::c_int = 0 as libc::c_int;
            while i_1 < p {
                let mut j_0: libc::c_int = 0 as libc::c_int;
                while j_0 < p {
                    let mut tmp_1: libc::c_double = *V.offset((i_1 + p * j_0) as isize);
                    let mut k_1: libc::c_int = 0 as libc::c_int;
                    while k_1 < p {
                        tmp_1 += *mm.offset((i_1 + p * k_1) as isize)
                            * *T.offset((j_0 + p * k_1) as isize);
                        k_1 += 1
                    }
                    *Pnew.offset((i_1 + p * j_0) as isize) = tmp_1;
                    j_0 += 1
                }
                i_1 += 1
            }
        }
        let mut i_2: libc::c_int = 0 as libc::c_int;
        while i_2 < p {
            *at.offset((l + n * i_2) as isize) = *anew.offset(i_2 as isize);
            i_2 += 1
        }
        let mut i_3: libc::c_int = 0 as libc::c_int;
        while i_3 < p * p {
            *Pt.offset((l + n * i_3) as isize) = *Pnew.offset(i_3 as isize);
            i_3 += 1
        }
        if !((*y.offset(l as isize)).is_nan() as i32 != 0 as libc::c_int) {
            let mut resid0: libc::c_double = *y.offset(l as isize);
            let mut i_4: libc::c_int = 0 as libc::c_int;
            while i_4 < p {
                resid0 -= *Z.offset(i_4 as isize) * *anew.offset(i_4 as isize);
                i_4 += 1
            }
            let mut gain: libc::c_double = h;
            let mut i_5: libc::c_int = 0 as libc::c_int;
            while i_5 < p {
                let mut tmp_2: libc::c_double = 0.0f64;
                let mut j_1: libc::c_int = 0 as libc::c_int;
                while j_1 < p {
                    tmp_2 += *Pnew.offset((i_5 + j_1 * p) as isize) * *Z.offset(j_1 as isize);
                    j_1 += 1
                }
                let ref mut fresh0 = *M.offset(i_5 as isize);
                *fresh0 = tmp_2;
                *Mt.offset((l + n * i_5) as isize) = *fresh0;
                gain += *Z.offset(i_5 as isize) * *M.offset(i_5 as isize);
                i_5 += 1
            }
            *gains.offset(l as isize) = gain;
            *resids.offset(l as isize) = resid0;
            let mut i_6: libc::c_int = 0 as libc::c_int;
            while i_6 < p {
                *a.offset(i_6 as isize) =
                    *anew.offset(i_6 as isize) + *M.offset(i_6 as isize) * resid0 / gain;
                i_6 += 1
            }
            let mut i_7: libc::c_int = 0 as libc::c_int;
            while i_7 < p {
                let mut j_2: libc::c_int = 0 as libc::c_int;
                while j_2 < p {
                    *P.offset((i_7 + j_2 * p) as isize) = *Pnew.offset((i_7 + j_2 * p) as isize)
                        - *M.offset(i_7 as isize) * *M.offset(j_2 as isize) / gain;
                    j_2 += 1
                }
                i_7 += 1
            }
        } else {
            let mut i_8: libc::c_int = 0 as libc::c_int;
            while i_8 < p {
                *a.offset(i_8 as isize) = *anew.offset(i_8 as isize);
                *Mt.offset((l + n * i_8) as isize) = 0.0f64;
                i_8 += 1
            }
            let mut i_9: libc::c_int = 0 as libc::c_int;
            while i_9 < p * p {
                *P.offset(i_9 as isize) = *Pnew.offset(i_9 as isize);
                i_9 += 1
            }
            *gains.offset(l as isize) = R_NaReal;
            *resids.offset(l as isize) = R_NaReal
        }
        l += 1
    }
    /* rt stores r_{t-1} */
    rt = R_alloc(
        (n * p) as size_t,
        ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
    ) as *mut libc::c_double;
    let mut l_0: libc::c_int = n - 1 as libc::c_int;
    while l_0 >= 0 as libc::c_int {
        if !((*gains.offset(l_0 as isize)).is_nan() as i32 != 0 as libc::c_int) {
            gn = 1 as libc::c_int as libc::c_double / *gains.offset(l_0 as isize);
            let mut i_10: libc::c_int = 0 as libc::c_int;
            while i_10 < p {
                *rt.offset((l_0 + n * i_10) as isize) =
                    *Z.offset(i_10 as isize) * *resids.offset(l_0 as isize) * gn;
                i_10 += 1
            }
        } else {
            let mut i_11: libc::c_int = 0 as libc::c_int;
            while i_11 < p {
                *rt.offset((l_0 + n * i_11) as isize) = 0.0f64;
                i_11 += 1
            }
            gn = 0.0f64
        }
        if var as u64 != 0 {
            let mut i_12: libc::c_int = 0 as libc::c_int;
            while i_12 < p {
                let mut j_3: libc::c_int = 0 as libc::c_int;
                while j_3 < p {
                    *Nt.offset((l_0 + n * i_12 + n * p * j_3) as isize) =
                        *Z.offset(i_12 as isize) * *Z.offset(j_3 as isize) * gn;
                    j_3 += 1
                }
                i_12 += 1
            }
        }
        if l_0 < n - 1 as libc::c_int {
            /* compute r_{t-1} */
            let mut i_13: libc::c_int = 0 as libc::c_int;
            while i_13 < p {
                let mut j_4: libc::c_int = 0 as libc::c_int;
                while j_4 < p {
                    *mm.offset((i_13 + p * j_4) as isize) = (if i_13 == j_4 {
                        1 as libc::c_int
                    } else {
                        0 as libc::c_int
                    })
                        as libc::c_double
                        - *Mt.offset((l_0 + n * i_13) as isize) * *Z.offset(j_4 as isize) * gn;
                    j_4 += 1
                }
                i_13 += 1
            }
            let mut i_14: libc::c_int = 0 as libc::c_int;
            while i_14 < p {
                let mut j_5: libc::c_int = 0 as libc::c_int;
                while j_5 < p {
                    let mut tmp_3: libc::c_double = 0.0f64;
                    let mut k_2: libc::c_int = 0 as libc::c_int;
                    while k_2 < p {
                        tmp_3 += *T.offset((i_14 + p * k_2) as isize)
                            * *mm.offset((k_2 + p * j_5) as isize);
                        k_2 += 1
                    }
                    *L.offset((i_14 + p * j_5) as isize) = tmp_3;
                    j_5 += 1
                }
                i_14 += 1
            }
            let mut i_15: libc::c_int = 0 as libc::c_int;
            while i_15 < p {
                let mut tmp_4: libc::c_double = 0.0f64;
                let mut j_6: libc::c_int = 0 as libc::c_int;
                while j_6 < p {
                    tmp_4 += *L.offset((j_6 + p * i_15) as isize)
                        * *rt.offset((l_0 + 1 as libc::c_int + n * j_6) as isize);
                    j_6 += 1
                }
                *rt.offset((l_0 + n * i_15) as isize) += tmp_4;
                i_15 += 1
            }
            if var as u64 != 0 {
                /* compute N_{t-1} */
                let mut i_16: libc::c_int = 0 as libc::c_int;
                while i_16 < p {
                    let mut j_7: libc::c_int = 0 as libc::c_int;
                    while j_7 < p {
                        let mut tmp_5: libc::c_double = 0.0f64;
                        let mut k_3: libc::c_int = 0 as libc::c_int;
                        while k_3 < p {
                            tmp_5 += *L.offset((k_3 + p * i_16) as isize)
                                * *Nt.offset(
                                    (l_0 + 1 as libc::c_int + n * k_3 + n * p * j_7) as isize,
                                );
                            k_3 += 1
                        }
                        *mm.offset((i_16 + p * j_7) as isize) = tmp_5;
                        j_7 += 1
                    }
                    i_16 += 1
                }
                let mut i_17: libc::c_int = 0 as libc::c_int;
                while i_17 < p {
                    let mut j_8: libc::c_int = 0 as libc::c_int;
                    while j_8 < p {
                        let mut tmp_6: libc::c_double = 0.0f64;
                        let mut k_4: libc::c_int = 0 as libc::c_int;
                        while k_4 < p {
                            tmp_6 += *mm.offset((i_17 + p * k_4) as isize)
                                * *L.offset((k_4 + p * j_8) as isize);
                            k_4 += 1
                        }
                        *Nt.offset((l_0 + n * i_17 + n * p * j_8) as isize) += tmp_6;
                        j_8 += 1
                    }
                    i_17 += 1
                }
            }
        }
        let mut i_18: libc::c_int = 0 as libc::c_int;
        while i_18 < p {
            let mut tmp_7: libc::c_double = 0.0f64;
            let mut j_9: libc::c_int = 0 as libc::c_int;
            while j_9 < p {
                tmp_7 += *Pt.offset((l_0 + n * i_18 + n * p * j_9) as isize)
                    * *rt.offset((l_0 + n * j_9) as isize);
                j_9 += 1
            }
            *at.offset((l_0 + n * i_18) as isize) += tmp_7;
            i_18 += 1
        }
        l_0 -= 1
    }
    if var as u64 != 0 {
        let mut l_1: libc::c_int = 0 as libc::c_int;
        while l_1 < n {
            let mut i_19: libc::c_int = 0 as libc::c_int;
            while i_19 < p {
                let mut j_10: libc::c_int = 0 as libc::c_int;
                while j_10 < p {
                    let mut tmp_8: libc::c_double = 0.0f64;
                    let mut k_5: libc::c_int = 0 as libc::c_int;
                    while k_5 < p {
                        tmp_8 += *Pt.offset((l_1 + n * i_19 + n * p * k_5) as isize)
                            * *Nt.offset((l_1 + n * k_5 + n * p * j_10) as isize);
                        k_5 += 1
                    }
                    *mm.offset((i_19 + p * j_10) as isize) = tmp_8;
                    j_10 += 1
                }
                i_19 += 1
            }
            let mut i_20: libc::c_int = 0 as libc::c_int;
            while i_20 < p {
                let mut j_11: libc::c_int = 0 as libc::c_int;
                while j_11 < p {
                    let mut tmp_9: libc::c_double =
                        *Pt.offset((l_1 + n * i_20 + n * p * j_11) as isize);
                    let mut k_6: libc::c_int = 0 as libc::c_int;
                    while k_6 < p {
                        tmp_9 -= *mm.offset((i_20 + p * k_6) as isize)
                            * *Pt.offset((l_1 + n * k_6 + n * p * j_11) as isize);
                        k_6 += 1
                    }
                    *Nt.offset((l_1 + n * i_20 + n * p * j_11) as isize) = tmp_9;
                    j_11 += 1
                }
                i_20 += 1
            }
            l_1 += 1
        }
    }
    unprotect(4 as libc::c_int);
    return res;
}
#[no_mangle]
pub unsafe extern "C" fn KalmanFore(mut nahead: SEXP, mut mod_0: SEXP, mut update: SEXP) -> SEXP {
    mod_0 = protect(duplicate(mod_0));
    let mut sZ: SEXP = getListElement(
        mod_0,
        b"Z\x00" as *const u8 as *const libc::c_char as *mut libc::c_char,
    );
    let mut sa: SEXP = getListElement(
        mod_0,
        b"a\x00" as *const u8 as *const libc::c_char as *mut libc::c_char,
    );
    let mut sP: SEXP = getListElement(
        mod_0,
        b"P\x00" as *const u8 as *const libc::c_char as *mut libc::c_char,
    );
    let mut sT: SEXP = getListElement(
        mod_0,
        b"T\x00" as *const u8 as *const libc::c_char as *mut libc::c_char,
    );
    let mut sV: SEXP = getListElement(
        mod_0,
        b"V\x00" as *const u8 as *const libc::c_char as *mut libc::c_char,
    );
    let mut sh: SEXP = getListElement(
        mod_0,
        b"h\x00" as *const u8 as *const libc::c_char as *mut libc::c_char,
    );
    if TYPEOF(sZ) != 14 as libc::c_int
        || TYPEOF(sa) != 14 as libc::c_int
        || TYPEOF(sP) != 14 as libc::c_int
        || TYPEOF(sT) != 14 as libc::c_int
        || TYPEOF(sV) != 14 as libc::c_int
    {
        error(dcgettext(
            b"stats\x00" as *const u8 as *const libc::c_char,
            b"invalid argument type\x00" as *const u8 as *const libc::c_char,
            5 as libc::c_int,
        ));
    }
    let mut n: libc::c_int = asInteger(nahead);
    let mut p: libc::c_int = LENGTH(sa);
    let mut Z: *mut libc::c_double = REAL(sZ);
    let mut a: *mut libc::c_double = REAL(sa);
    let mut P: *mut libc::c_double = REAL(sP);
    let mut T: *mut libc::c_double = REAL(sT);
    let mut V: *mut libc::c_double = REAL(sV);
    let mut h: libc::c_double = asReal(sh);
    let mut mm: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut anew: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut Pnew: *mut libc::c_double = 0 as *mut libc::c_double;
    anew = R_alloc(
        p as size_t,
        ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
    ) as *mut libc::c_double;
    Pnew = R_alloc(
        (p * p) as size_t,
        ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
    ) as *mut libc::c_double;
    mm = R_alloc(
        (p * p) as size_t,
        ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
    ) as *mut libc::c_double;
    let mut res: SEXP = 0 as *mut SEXPREC;
    let mut forecasts: SEXP = 0 as *mut SEXPREC;
    let mut se: SEXP = 0 as *mut SEXPREC;
    res = allocVector(19 as libc::c_int as SEXPTYPE, 2 as libc::c_int as R_xlen_t);
    protect(res);
    forecasts = allocVector(14 as libc::c_int as SEXPTYPE, n as R_xlen_t);
    SET_VECTOR_ELT(res, 0 as libc::c_int as R_xlen_t, forecasts);
    se = allocVector(14 as libc::c_int as SEXPTYPE, n as R_xlen_t);
    SET_VECTOR_ELT(res, 1 as libc::c_int as R_xlen_t, se);
    let mut nm: SEXP = protect(allocVector(
        16 as libc::c_int as SEXPTYPE,
        2 as libc::c_int as R_xlen_t,
    ));
    SET_STRING_ELT(
        nm,
        0 as libc::c_int as R_xlen_t,
        mkChar(b"pred\x00" as *const u8 as *const libc::c_char),
    );
    SET_STRING_ELT(
        nm,
        1 as libc::c_int as R_xlen_t,
        mkChar(b"var\x00" as *const u8 as *const libc::c_char),
    );
    setAttrib(res, R_NamesSymbol, nm);
    unprotect(1 as libc::c_int);
    let mut l: libc::c_int = 0 as libc::c_int;
    while l < n {
        let mut fc: libc::c_double = 0.0f64;
        let mut i: libc::c_int = 0 as libc::c_int;
        while i < p {
            let mut tmp: libc::c_double = 0.0f64;
            let mut k: libc::c_int = 0 as libc::c_int;
            while k < p {
                tmp += *T.offset((i + p * k) as isize) * *a.offset(k as isize);
                k += 1
            }
            *anew.offset(i as isize) = tmp;
            fc += tmp * *Z.offset(i as isize);
            i += 1
        }
        let mut i_0: libc::c_int = 0 as libc::c_int;
        while i_0 < p {
            *a.offset(i_0 as isize) = *anew.offset(i_0 as isize);
            i_0 += 1
        }
        *REAL(forecasts).offset(l as isize) = fc;
        let mut i_1: libc::c_int = 0 as libc::c_int;
        while i_1 < p {
            let mut j: libc::c_int = 0 as libc::c_int;
            while j < p {
                let mut tmp_0: libc::c_double = 0.0f64;
                let mut k_0: libc::c_int = 0 as libc::c_int;
                while k_0 < p {
                    tmp_0 +=
                        *T.offset((i_1 + p * k_0) as isize) * *P.offset((k_0 + p * j) as isize);
                    k_0 += 1
                }
                *mm.offset((i_1 + p * j) as isize) = tmp_0;
                j += 1
            }
            i_1 += 1
        }
        let mut i_2: libc::c_int = 0 as libc::c_int;
        while i_2 < p {
            let mut j_0: libc::c_int = 0 as libc::c_int;
            while j_0 < p {
                let mut tmp_1: libc::c_double = *V.offset((i_2 + p * j_0) as isize);
                let mut k_1: libc::c_int = 0 as libc::c_int;
                while k_1 < p {
                    tmp_1 +=
                        *mm.offset((i_2 + p * k_1) as isize) * *T.offset((j_0 + p * k_1) as isize);
                    k_1 += 1
                }
                *Pnew.offset((i_2 + p * j_0) as isize) = tmp_1;
                j_0 += 1
            }
            i_2 += 1
        }
        let mut tmp_2: libc::c_double = h;
        let mut i_3: libc::c_int = 0 as libc::c_int;
        while i_3 < p {
            let mut j_1: libc::c_int = 0 as libc::c_int;
            while j_1 < p {
                *P.offset((i_3 + j_1 * p) as isize) = *Pnew.offset((i_3 + j_1 * p) as isize);
                tmp_2 += *Z.offset(i_3 as isize)
                    * *Z.offset(j_1 as isize)
                    * *P.offset((i_3 + j_1 * p) as isize);
                j_1 += 1
            }
            i_3 += 1
        }
        *REAL(se).offset(l as isize) = tmp_2;
        l += 1
    }
    if asLogical(update) != 0 {
        setAttrib(
            res,
            install(b"mod\x00" as *const u8 as *const libc::c_char),
            mod_0,
        );
    }
    unprotect(2 as libc::c_int);
    return res;
}
unsafe extern "C" fn partrans(
    mut p: libc::c_int,
    mut raw: *mut libc::c_double,
    mut new: *mut libc::c_double,
) {
    let mut j: libc::c_int = 0;
    let mut k: libc::c_int = 0;
    let mut a: libc::c_double = 0.;
    let mut work: [libc::c_double; 100] = [0.; 100];
    if p > 100 as libc::c_int {
        error(dcgettext(
            b"stats\x00" as *const u8 as *const libc::c_char,
            b"can only transform 100 pars in arima0\x00" as *const u8 as *const libc::c_char,
            5 as libc::c_int,
        ));
    }
    /* Step one: map (-Inf, Inf) to (-1, 1) via tanh
    The parameters are now the pacf phi_{kk} */
    j = 0 as libc::c_int;
    while j < p {
        let ref mut fresh1 = *new.offset(j as isize);
        *fresh1 = (*raw.offset(j as isize)).tanh();
        work[j as usize] = *fresh1;
        j += 1
    }
    /* Step two: run the Durbin-Levinson recursions to find phi_{j.},
    j = 2, ..., p and phi_{p.} are the autoregression coefficients */
    j = 1 as libc::c_int;
    while j < p {
        a = *new.offset(j as isize);
        k = 0 as libc::c_int;
        while k < j {
            work[k as usize] -= a * *new.offset((j - k - 1 as libc::c_int) as isize);
            k += 1
        }
        k = 0 as libc::c_int;
        while k < j {
            *new.offset(k as isize) = work[k as usize];
            k += 1
        }
        j += 1
    }
}
#[no_mangle]
pub unsafe extern "C" fn ARIMA_undoPars(mut sin: SEXP, mut sarma: SEXP) -> SEXP {
    let mut arma: *mut libc::c_int = INTEGER(sarma);
    let mut mp: libc::c_int = *arma.offset(0 as libc::c_int as isize);
    let mut mq: libc::c_int = *arma.offset(1 as libc::c_int as isize);
    let mut msp: libc::c_int = *arma.offset(2 as libc::c_int as isize);
    let mut v: libc::c_int = 0;
    let mut n: libc::c_int = LENGTH(sin);
    let mut params: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut in_0: *mut libc::c_double = REAL(sin);
    let mut res: SEXP = allocVector(14 as libc::c_int as SEXPTYPE, n as R_xlen_t);
    params = REAL(res);
    let mut i: libc::c_int = 0 as libc::c_int;
    while i < n {
        *params.offset(i as isize) = *in_0.offset(i as isize);
        i += 1
    }
    if mp > 0 as libc::c_int {
        partrans(mp, in_0, params);
    }
    v = mp + mq;
    if msp > 0 as libc::c_int {
        partrans(msp, in_0.offset(v as isize), params.offset(v as isize));
    }
    return res;
}
#[no_mangle]
pub unsafe extern "C" fn ARIMA_transPars(mut sin: SEXP, mut sarma: SEXP, mut strans: SEXP) -> SEXP {
    let mut arma: *mut libc::c_int = INTEGER(sarma);
    let mut trans: libc::c_int = asLogical(strans);
    let mut mp: libc::c_int = *arma.offset(0 as libc::c_int as isize);
    let mut mq: libc::c_int = *arma.offset(1 as libc::c_int as isize);
    let mut msp: libc::c_int = *arma.offset(2 as libc::c_int as isize);
    let mut msq: libc::c_int = *arma.offset(3 as libc::c_int as isize);
    let mut ns: libc::c_int = *arma.offset(4 as libc::c_int as isize);
    let mut i: libc::c_int = 0;
    let mut j: libc::c_int = 0;
    let mut p: libc::c_int = mp + ns * msp;
    let mut q: libc::c_int = mq + ns * msq;
    let mut v: libc::c_int = 0;
    let mut in_0: *mut libc::c_double = REAL(sin);
    let mut params: *mut libc::c_double = REAL(sin);
    let mut phi: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut theta: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut res: SEXP = 0 as *mut SEXPREC;
    let mut sPhi: SEXP = 0 as *mut SEXPREC;
    let mut sTheta: SEXP = 0 as *mut SEXPREC;
    res = allocVector(19 as libc::c_int as SEXPTYPE, 2 as libc::c_int as R_xlen_t);
    protect(res);
    sPhi = allocVector(14 as libc::c_int as SEXPTYPE, p as R_xlen_t);
    SET_VECTOR_ELT(res, 0 as libc::c_int as R_xlen_t, sPhi);
    sTheta = allocVector(14 as libc::c_int as SEXPTYPE, q as R_xlen_t);
    SET_VECTOR_ELT(res, 1 as libc::c_int as R_xlen_t, sTheta);
    phi = REAL(sPhi);
    theta = REAL(sTheta);
    if trans != 0 {
        let mut n: libc::c_int = mp + mq + msp + msq;
        params = R_alloc(
            n as size_t,
            ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
        ) as *mut libc::c_double;
        i = 0 as libc::c_int;
        while i < n {
            *params.offset(i as isize) = *in_0.offset(i as isize);
            i += 1
        }
        if mp > 0 as libc::c_int {
            partrans(mp, in_0, params);
        }
        v = mp + mq;
        if msp > 0 as libc::c_int {
            partrans(msp, in_0.offset(v as isize), params.offset(v as isize));
        }
    }
    if ns > 0 as libc::c_int {
        /* expand out seasonal ARMA models */
        i = 0 as libc::c_int;
        while i < mp {
            *phi.offset(i as isize) = *params.offset(i as isize);
            i += 1
        }
        i = 0 as libc::c_int;
        while i < mq {
            *theta.offset(i as isize) = *params.offset((i + mp) as isize);
            i += 1
        }
        i = mp;
        while i < p {
            *phi.offset(i as isize) = 0.0f64;
            i += 1
        }
        i = mq;
        while i < q {
            *theta.offset(i as isize) = 0.0f64;
            i += 1
        }
        j = 0 as libc::c_int;
        while j < msp {
            *phi.offset(((j + 1 as libc::c_int) * ns - 1 as libc::c_int) as isize) +=
                *params.offset((j + mp + mq) as isize);
            i = 0 as libc::c_int;
            while i < mp {
                *phi.offset(((j + 1 as libc::c_int) * ns + i) as isize) -=
                    *params.offset(i as isize) * *params.offset((j + mp + mq) as isize);
                i += 1
            }
            j += 1
        }
        j = 0 as libc::c_int;
        while j < msq {
            *theta.offset(((j + 1 as libc::c_int) * ns - 1 as libc::c_int) as isize) +=
                *params.offset((j + mp + mq + msp) as isize);
            i = 0 as libc::c_int;
            while i < mq {
                *theta.offset(((j + 1 as libc::c_int) * ns + i) as isize) += *params
                    .offset((i + mp) as isize)
                    * *params.offset((j + mp + mq + msp) as isize);
                i += 1
            }
            j += 1
        }
    } else {
        i = 0 as libc::c_int;
        while i < mp {
            *phi.offset(i as isize) = *params.offset(i as isize);
            i += 1
        }
        i = 0 as libc::c_int;
        while i < mq {
            *theta.offset(i as isize) = *params.offset((i + mp) as isize);
            i += 1
        }
    }
    unprotect(1 as libc::c_int);
    return res;
}
unsafe extern "C" fn invpartrans(
    mut p: libc::c_int,
    mut phi: *mut libc::c_double,
    mut new: *mut libc::c_double,
) {
    let mut j: libc::c_int = 0;
    let mut k: libc::c_int = 0;
    let mut a: libc::c_double = 0.;
    let mut work: [libc::c_double; 100] = [0.; 100];
    if p > 100 as libc::c_int {
        error(dcgettext(
            b"stats\x00" as *const u8 as *const libc::c_char,
            b"can only transform 100 pars in arima0\x00" as *const u8 as *const libc::c_char,
            5 as libc::c_int,
        ));
    }
    j = 0 as libc::c_int;
    while j < p {
        let ref mut fresh2 = *new.offset(j as isize);
        *fresh2 = *phi.offset(j as isize);
        work[j as usize] = *fresh2;
        j += 1
    }
    /* Run the Durbin-Levinson recursions backwards
    to find the PACF phi_{j.} from the autoregression coefficients */
    j = p - 1 as libc::c_int;
    while j > 0 as libc::c_int {
        a = *new.offset(j as isize);
        k = 0 as libc::c_int;
        while k < j {
            work[k as usize] = (*new.offset(k as isize)
                + a * *new.offset((j - k - 1 as libc::c_int) as isize))
                / (1 as libc::c_int as libc::c_double - a * a);
            k += 1
        }
        k = 0 as libc::c_int;
        while k < j {
            *new.offset(k as isize) = work[k as usize];
            k += 1
        }
        j -= 1
    }
    j = 0 as libc::c_int;
    while j < p {
        *new.offset(j as isize) = (*new.offset(j as isize)).atanh();
        j += 1
    }
}
#[no_mangle]
pub unsafe extern "C" fn ARIMA_Invtrans(mut in_0: SEXP, mut sarma: SEXP) -> SEXP {
    let mut arma: *mut libc::c_int = INTEGER(sarma);
    let mut mp: libc::c_int = *arma.offset(0 as libc::c_int as isize);
    let mut mq: libc::c_int = *arma.offset(1 as libc::c_int as isize);
    let mut msp: libc::c_int = *arma.offset(2 as libc::c_int as isize);
    let mut i: libc::c_int = 0;
    let mut v: libc::c_int = 0;
    let mut n: libc::c_int = LENGTH(in_0);
    let mut y: SEXP = allocVector(14 as libc::c_int as SEXPTYPE, n as R_xlen_t);
    let mut raw: *mut libc::c_double = REAL(in_0);
    let mut new: *mut libc::c_double = REAL(y);
    i = 0 as libc::c_int;
    while i < n {
        *new.offset(i as isize) = *raw.offset(i as isize);
        i += 1
    }
    if mp > 0 as libc::c_int {
        invpartrans(mp, raw, new);
    }
    v = mp + mq;
    if msp > 0 as libc::c_int {
        invpartrans(msp, raw.offset(v as isize), new.offset(v as isize));
    }
    return y;
}
#[no_mangle]
pub unsafe extern "C" fn ARIMA_Gradtrans(mut in_0: SEXP, mut sarma: SEXP) -> SEXP {
    let mut arma: *mut libc::c_int = INTEGER(sarma);
    let mut mp: libc::c_int = *arma.offset(0 as libc::c_int as isize);
    let mut mq: libc::c_int = *arma.offset(1 as libc::c_int as isize);
    let mut msp: libc::c_int = *arma.offset(2 as libc::c_int as isize);
    let mut n: libc::c_int = LENGTH(in_0);
    let mut y: SEXP = allocMatrix(14 as libc::c_int as SEXPTYPE, n, n);
    let mut raw: *mut libc::c_double = REAL(in_0);
    let mut A: *mut libc::c_double = REAL(y);
    let mut w1: [libc::c_double; 100] = [0.; 100];
    let mut w2: [libc::c_double; 100] = [0.; 100];
    let mut w3: [libc::c_double; 100] = [0.; 100];
    let mut i: libc::c_int = 0 as libc::c_int;
    while i < n {
        let mut j: libc::c_int = 0 as libc::c_int;
        while j < n {
            *A.offset((i + j * n) as isize) = (i == j) as libc::c_int as libc::c_double;
            j += 1
        }
        i += 1
    }
    if mp > 0 as libc::c_int {
        let mut i_0: libc::c_int = 0 as libc::c_int;
        while i_0 < mp {
            w1[i_0 as usize] = *raw.offset(i_0 as isize);
            i_0 += 1
        }
        partrans(mp, w1.as_mut_ptr(), w2.as_mut_ptr());
        let mut i_1: libc::c_int = 0 as libc::c_int;
        while i_1 < mp {
            w1[i_1 as usize] += 1e-3f64;
            partrans(mp, w1.as_mut_ptr(), w3.as_mut_ptr());
            let mut j_0: libc::c_int = 0 as libc::c_int;
            while j_0 < mp {
                *A.offset((i_1 + j_0 * n) as isize) =
                    (w3[j_0 as usize] - w2[j_0 as usize]) / 1e-3f64;
                j_0 += 1
            }
            w1[i_1 as usize] -= 1e-3f64;
            i_1 += 1
        }
    }
    if msp > 0 as libc::c_int {
        let mut v: libc::c_int = mp + mq;
        let mut i_2: libc::c_int = 0 as libc::c_int;
        while i_2 < msp {
            w1[i_2 as usize] = *raw.offset((i_2 + v) as isize);
            i_2 += 1
        }
        partrans(msp, w1.as_mut_ptr(), w2.as_mut_ptr());
        let mut i_3: libc::c_int = 0 as libc::c_int;
        while i_3 < msp {
            w1[i_3 as usize] += 1e-3f64;
            partrans(msp, w1.as_mut_ptr(), w3.as_mut_ptr());
            let mut j_1: libc::c_int = 0 as libc::c_int;
            while j_1 < msp {
                *A.offset((i_3 + v + (j_1 + v) * n) as isize) =
                    (w3[j_1 as usize] - w2[j_1 as usize]) / 1e-3f64;
                j_1 += 1
            }
            w1[i_3 as usize] -= 1e-3f64;
            i_3 += 1
        }
    }
    return y;
}
#[no_mangle]
pub unsafe extern "C" fn ARIMA_Like(
    mut sy: SEXP,
    mut mod_0: SEXP,
    mut sUP: SEXP,
    mut giveResid: SEXP,
) -> SEXP {
    let mut sPhi: SEXP = getListElement(
        mod_0,
        b"phi\x00" as *const u8 as *const libc::c_char as *mut libc::c_char,
    );
    let mut sTheta: SEXP = getListElement(
        mod_0,
        b"theta\x00" as *const u8 as *const libc::c_char as *mut libc::c_char,
    );
    let mut sDelta: SEXP = getListElement(
        mod_0,
        b"Delta\x00" as *const u8 as *const libc::c_char as *mut libc::c_char,
    );
    let mut sa: SEXP = getListElement(
        mod_0,
        b"a\x00" as *const u8 as *const libc::c_char as *mut libc::c_char,
    );
    let mut sP: SEXP = getListElement(
        mod_0,
        b"P\x00" as *const u8 as *const libc::c_char as *mut libc::c_char,
    );
    let mut sPn: SEXP = getListElement(
        mod_0,
        b"Pn\x00" as *const u8 as *const libc::c_char as *mut libc::c_char,
    );
    if TYPEOF(sPhi) != 14 as libc::c_int
        || TYPEOF(sTheta) != 14 as libc::c_int
        || TYPEOF(sDelta) != 14 as libc::c_int
        || TYPEOF(sa) != 14 as libc::c_int
        || TYPEOF(sP) != 14 as libc::c_int
        || TYPEOF(sPn) != 14 as libc::c_int
    {
        error(dcgettext(
            b"stats\x00" as *const u8 as *const libc::c_char,
            b"invalid argument type\x00" as *const u8 as *const libc::c_char,
            5 as libc::c_int,
        ));
    }
    let mut res: SEXP = 0 as *mut SEXPREC;
    let mut nres: SEXP = 0 as *mut SEXPREC;
    let mut sResid: SEXP = R_NilValue;
    let mut n: libc::c_int = LENGTH(sy);
    let mut rd: libc::c_int = LENGTH(sa);
    let mut p: libc::c_int = LENGTH(sPhi);
    let mut q: libc::c_int = LENGTH(sTheta);
    let mut d: libc::c_int = LENGTH(sDelta);
    let mut r: libc::c_int = rd - d;
    let mut y: *mut libc::c_double = REAL(sy);
    let mut a: *mut libc::c_double = REAL(sa);
    let mut P: *mut libc::c_double = REAL(sP);
    let mut Pnew: *mut libc::c_double = REAL(sPn);
    let mut phi: *mut libc::c_double = REAL(sPhi);
    let mut theta: *mut libc::c_double = REAL(sTheta);
    let mut delta: *mut libc::c_double = REAL(sDelta);
    let mut sumlog: libc::c_double = 0.0f64;
    let mut ssq: libc::c_double = 0 as libc::c_int as libc::c_double;
    let mut anew: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut mm: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut M: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut nu: libc::c_int = 0 as libc::c_int;
    let mut useResid: Rboolean = asLogical(giveResid) as Rboolean;
    let mut rsResid: *mut libc::c_double = 0 as *mut libc::c_double;
    anew = R_alloc(
        rd as size_t,
        ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
    ) as *mut libc::c_double;
    M = R_alloc(
        rd as size_t,
        ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
    ) as *mut libc::c_double;
    if d > 0 as libc::c_int {
        mm = R_alloc(
            (rd * rd) as size_t,
            ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
        ) as *mut libc::c_double
    }
    if useResid as u64 != 0 {
        sResid = allocVector(14 as libc::c_int as SEXPTYPE, n as R_xlen_t);
        protect(sResid);
        rsResid = REAL(sResid)
    }
    let mut l: libc::c_int = 0 as libc::c_int;
    while l < n {
        let mut i: libc::c_int = 0 as libc::c_int;
        while i < r {
            let mut tmp: libc::c_double = if i < r - 1 as libc::c_int {
                *a.offset((i + 1 as libc::c_int) as isize)
            } else {
                0.0f64
            };
            if i < p {
                tmp += *phi.offset(i as isize) * *a.offset(0 as libc::c_int as isize)
            }
            *anew.offset(i as isize) = tmp;
            i += 1
        }
        if d > 0 as libc::c_int {
            let mut i_0: libc::c_int = r + 1 as libc::c_int;
            while i_0 < rd {
                *anew.offset(i_0 as isize) = *a.offset((i_0 - 1 as libc::c_int) as isize);
                i_0 += 1
            }
            let mut tmp_0: libc::c_double = *a.offset(0 as libc::c_int as isize);
            let mut i_1: libc::c_int = 0 as libc::c_int;
            while i_1 < d {
                tmp_0 += *delta.offset(i_1 as isize) * *a.offset((r + i_1) as isize);
                i_1 += 1
            }
            *anew.offset(r as isize) = tmp_0
        }
        if l > asInteger(sUP) {
            if d == 0 as libc::c_int {
                let mut i_2: libc::c_int = 0 as libc::c_int;
                while i_2 < r {
                    let mut vi: libc::c_double = 0.0f64;
                    if i_2 == 0 as libc::c_int {
                        vi = 1.0f64
                    } else if (i_2 - 1 as libc::c_int) < q {
                        vi = *theta.offset((i_2 - 1 as libc::c_int) as isize)
                    }
                    let mut j: libc::c_int = 0 as libc::c_int;
                    while j < r {
                        let mut tmp_1: libc::c_double = 0.0f64;
                        if j == 0 as libc::c_int {
                            tmp_1 = vi
                        } else if (j - 1 as libc::c_int) < q {
                            tmp_1 = vi * *theta.offset((j - 1 as libc::c_int) as isize)
                        }
                        if i_2 < p && j < p {
                            tmp_1 += *phi.offset(i_2 as isize)
                                * *phi.offset(j as isize)
                                * *P.offset(0 as libc::c_int as isize)
                        }
                        if i_2 < r - 1 as libc::c_int && j < r - 1 as libc::c_int {
                            tmp_1 += *P.offset(
                                (i_2 + 1 as libc::c_int + r * (j + 1 as libc::c_int)) as isize,
                            )
                        }
                        if i_2 < p && j < r - 1 as libc::c_int {
                            tmp_1 += *phi.offset(i_2 as isize)
                                * *P.offset((j + 1 as libc::c_int) as isize)
                        }
                        if j < p && i_2 < r - 1 as libc::c_int {
                            tmp_1 += *phi.offset(j as isize)
                                * *P.offset((i_2 + 1 as libc::c_int) as isize)
                        }
                        *Pnew.offset((i_2 + r * j) as isize) = tmp_1;
                        j += 1
                    }
                    i_2 += 1
                }
            } else {
                /* -Wall */
                /* mm = TP */
                let mut i_3: libc::c_int = 0 as libc::c_int;
                while i_3 < r {
                    let mut j_0: libc::c_int = 0 as libc::c_int;
                    while j_0 < rd {
                        let mut tmp_2: libc::c_double = 0.0f64;
                        if i_3 < p {
                            tmp_2 += *phi.offset(i_3 as isize) * *P.offset((rd * j_0) as isize)
                        }
                        if i_3 < r - 1 as libc::c_int {
                            tmp_2 += *P.offset((i_3 + 1 as libc::c_int + rd * j_0) as isize)
                        }
                        *mm.offset((i_3 + rd * j_0) as isize) = tmp_2;
                        j_0 += 1
                    }
                    i_3 += 1
                }
                let mut j_1: libc::c_int = 0 as libc::c_int;
                while j_1 < rd {
                    let mut tmp_3: libc::c_double = *P.offset((rd * j_1) as isize);
                    let mut k: libc::c_int = 0 as libc::c_int;
                    while k < d {
                        tmp_3 += *delta.offset(k as isize) * *P.offset((r + k + rd * j_1) as isize);
                        k += 1
                    }
                    *mm.offset((r + rd * j_1) as isize) = tmp_3;
                    j_1 += 1
                }
                let mut i_4: libc::c_int = 1 as libc::c_int;
                while i_4 < d {
                    let mut j_2: libc::c_int = 0 as libc::c_int;
                    while j_2 < rd {
                        *mm.offset((r + i_4 + rd * j_2) as isize) =
                            *P.offset((r + i_4 - 1 as libc::c_int + rd * j_2) as isize);
                        j_2 += 1
                    }
                    i_4 += 1
                }
                /* Pnew = mmT' */
                let mut i_5: libc::c_int = 0 as libc::c_int;
                while i_5 < r {
                    let mut j_3: libc::c_int = 0 as libc::c_int;
                    while j_3 < rd {
                        let mut tmp_4: libc::c_double = 0.0f64;
                        if i_5 < p {
                            tmp_4 += *phi.offset(i_5 as isize) * *mm.offset(j_3 as isize)
                        }
                        if i_5 < r - 1 as libc::c_int {
                            tmp_4 += *mm.offset((rd * (i_5 + 1 as libc::c_int) + j_3) as isize)
                        }
                        *Pnew.offset((j_3 + rd * i_5) as isize) = tmp_4;
                        j_3 += 1
                    }
                    i_5 += 1
                }
                let mut j_4: libc::c_int = 0 as libc::c_int;
                while j_4 < rd {
                    let mut tmp_5: libc::c_double = *mm.offset(j_4 as isize);
                    let mut k_0: libc::c_int = 0 as libc::c_int;
                    while k_0 < d {
                        tmp_5 += *delta.offset(k_0 as isize)
                            * *mm.offset((rd * (r + k_0) + j_4) as isize);
                        k_0 += 1
                    }
                    *Pnew.offset((rd * r + j_4) as isize) = tmp_5;
                    j_4 += 1
                }
                let mut i_6: libc::c_int = 1 as libc::c_int;
                while i_6 < d {
                    let mut j_5: libc::c_int = 0 as libc::c_int;
                    while j_5 < rd {
                        *Pnew.offset((rd * (r + i_6) + j_5) as isize) =
                            *mm.offset((rd * (r + i_6 - 1 as libc::c_int) + j_5) as isize);
                        j_5 += 1
                    }
                    i_6 += 1
                }
                /* Pnew <- Pnew + (1 theta) %o% (1 theta) */
                let mut i_7: libc::c_int = 0 as libc::c_int;
                while i_7 <= q {
                    let mut vi_0: libc::c_double = if i_7 == 0 as libc::c_int {
                        1.0f64
                    } else {
                        *theta.offset((i_7 - 1 as libc::c_int) as isize)
                    };
                    let mut j_6: libc::c_int = 0 as libc::c_int;
                    while j_6 <= q {
                        *Pnew.offset((i_7 + rd * j_6) as isize) += vi_0
                            * (if j_6 == 0 as libc::c_int {
                                1.0f64
                            } else {
                                *theta.offset((j_6 - 1 as libc::c_int) as isize)
                            });
                        j_6 += 1
                    }
                    i_7 += 1
                }
            }
        }
        if !((*y.offset(l as isize)).is_nan() as i32 != 0 as libc::c_int) {
            let mut resid: libc::c_double =
                *y.offset(l as isize) - *anew.offset(0 as libc::c_int as isize);
            let mut i_8: libc::c_int = 0 as libc::c_int;
            while i_8 < d {
                resid -= *delta.offset(i_8 as isize) * *anew.offset((r + i_8) as isize);
                i_8 += 1
            }
            let mut i_9: libc::c_int = 0 as libc::c_int;
            while i_9 < rd {
                let mut tmp_6: libc::c_double = *Pnew.offset(i_9 as isize);
                let mut j_7: libc::c_int = 0 as libc::c_int;
                while j_7 < d {
                    tmp_6 +=
                        *Pnew.offset((i_9 + (r + j_7) * rd) as isize) * *delta.offset(j_7 as isize);
                    j_7 += 1
                }
                *M.offset(i_9 as isize) = tmp_6;
                i_9 += 1
            }
            let mut gain: libc::c_double = *M.offset(0 as libc::c_int as isize);
            let mut j_8: libc::c_int = 0 as libc::c_int;
            while j_8 < d {
                gain += *delta.offset(j_8 as isize) * *M.offset((r + j_8) as isize);
                j_8 += 1
            }
            if gain < 1e4f64 {
                nu += 1;
                ssq += resid * resid / gain;
                sumlog += gain.ln()
            }
            if useResid as u64 != 0 {
                *rsResid.offset(l as isize) = resid / sqrt(gain)
            }
            let mut i_10: libc::c_int = 0 as libc::c_int;
            while i_10 < rd {
                *a.offset(i_10 as isize) =
                    *anew.offset(i_10 as isize) + *M.offset(i_10 as isize) * resid / gain;
                i_10 += 1
            }
            let mut i_11: libc::c_int = 0 as libc::c_int;
            while i_11 < rd {
                let mut j_9: libc::c_int = 0 as libc::c_int;
                while j_9 < rd {
                    *P.offset((i_11 + j_9 * rd) as isize) = *Pnew
                        .offset((i_11 + j_9 * rd) as isize)
                        - *M.offset(i_11 as isize) * *M.offset(j_9 as isize) / gain;
                    j_9 += 1
                }
                i_11 += 1
            }
        } else {
            let mut i_12: libc::c_int = 0 as libc::c_int;
            while i_12 < rd {
                *a.offset(i_12 as isize) = *anew.offset(i_12 as isize);
                i_12 += 1
            }
            let mut i_13: libc::c_int = 0 as libc::c_int;
            while i_13 < rd * rd {
                *P.offset(i_13 as isize) = *Pnew.offset(i_13 as isize);
                i_13 += 1
            }
            if useResid as u64 != 0 {
                *rsResid.offset(l as isize) = R_NaReal
            }
        }
        l += 1
    }
    if useResid as u64 != 0 {
        res = allocVector(19 as libc::c_int as SEXPTYPE, 3 as libc::c_int as R_xlen_t);
        protect(res);
        nres = allocVector(14 as libc::c_int as SEXPTYPE, 3 as libc::c_int as R_xlen_t);
        SET_VECTOR_ELT(res, 0 as libc::c_int as R_xlen_t, nres);
        *REAL(nres).offset(0 as libc::c_int as isize) = ssq;
        *REAL(nres).offset(1 as libc::c_int as isize) = sumlog;
        *REAL(nres).offset(2 as libc::c_int as isize) = nu as libc::c_double;
        SET_VECTOR_ELT(res, 1 as libc::c_int as R_xlen_t, sResid);
        unprotect(2 as libc::c_int);
        return res;
    } else {
        nres = allocVector(14 as libc::c_int as SEXPTYPE, 3 as libc::c_int as R_xlen_t);
        *REAL(nres).offset(0 as libc::c_int as isize) = ssq;
        *REAL(nres).offset(1 as libc::c_int as isize) = sumlog;
        *REAL(nres).offset(2 as libc::c_int as isize) = nu as libc::c_double;
        return nres;
    };
}
/* do differencing here */
/* arma is p, q, sp, sq, ns, d, sd */
#[no_mangle]
pub unsafe extern "C" fn ARIMA_CSS(
    mut sy: SEXP,
    mut sarma: SEXP,
    mut sPhi: SEXP,
    mut sTheta: SEXP,
    mut sncond: SEXP,
    mut giveResid: SEXP,
) -> SEXP {
    let mut res: SEXP = 0 as *mut SEXPREC;
    let mut sResid: SEXP = R_NilValue;
    let mut ssq: libc::c_double = 0.0f64;
    let mut y: *mut libc::c_double = REAL(sy);
    let mut tmp: libc::c_double = 0.;
    let mut phi: *mut libc::c_double = REAL(sPhi);
    let mut theta: *mut libc::c_double = REAL(sTheta);
    let mut w: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut resid: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut n: libc::c_int = LENGTH(sy);
    let mut arma: *mut libc::c_int = INTEGER(sarma);
    let mut p: libc::c_int = LENGTH(sPhi);
    let mut q: libc::c_int = LENGTH(sTheta);
    let mut ncond: libc::c_int = asInteger(sncond);
    let mut ns: libc::c_int = 0;
    let mut nu: libc::c_int = 0 as libc::c_int;
    let mut useResid: Rboolean = asLogical(giveResid) as Rboolean;
    w = R_alloc(
        n as size_t,
        ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
    ) as *mut libc::c_double;
    let mut l: libc::c_int = 0 as libc::c_int;
    while l < n {
        *w.offset(l as isize) = *y.offset(l as isize);
        l += 1
    }
    let mut i: libc::c_int = 0 as libc::c_int;
    while i < *arma.offset(5 as libc::c_int as isize) {
        let mut l_0: libc::c_int = n - 1 as libc::c_int;
        while l_0 > 0 as libc::c_int {
            *w.offset(l_0 as isize) -= *w.offset((l_0 - 1 as libc::c_int) as isize);
            l_0 -= 1
        }
        i += 1
    }
    ns = *arma.offset(4 as libc::c_int as isize);
    let mut i_0: libc::c_int = 0 as libc::c_int;
    while i_0 < *arma.offset(6 as libc::c_int as isize) {
        let mut l_1: libc::c_int = n - 1 as libc::c_int;
        while l_1 >= ns {
            *w.offset(l_1 as isize) -= *w.offset((l_1 - ns) as isize);
            l_1 -= 1
        }
        i_0 += 1
    }
    sResid = allocVector(14 as libc::c_int as SEXPTYPE, n as R_xlen_t);
    protect(sResid);
    resid = REAL(sResid);
    if useResid as u64 != 0 {
        let mut l_2: libc::c_int = 0 as libc::c_int;
        while l_2 < ncond {
            *resid.offset(l_2 as isize) = 0 as libc::c_int as libc::c_double;
            l_2 += 1
        }
    }
    let mut l_3: libc::c_int = ncond;
    while l_3 < n {
        tmp = *w.offset(l_3 as isize);
        let mut j: libc::c_int = 0 as libc::c_int;
        while j < p {
            tmp -= *phi.offset(j as isize) * *w.offset((l_3 - j - 1 as libc::c_int) as isize);
            j += 1
        }
        let mut j_0: libc::c_int = 0 as libc::c_int;
        while j_0 < (if l_3 - ncond < q { (l_3) - ncond } else { q }) {
            tmp -= *theta.offset(j_0 as isize)
                * *resid.offset((l_3 - j_0 - 1 as libc::c_int) as isize);
            j_0 += 1
        }
        *resid.offset(l_3 as isize) = tmp;
        if !(tmp.is_nan() as i32 != 0 as libc::c_int) {
            nu += 1;
            ssq += tmp * tmp
        }
        l_3 += 1
    }
    if useResid as u64 != 0 {
        res = allocVector(19 as libc::c_int as SEXPTYPE, 2 as libc::c_int as R_xlen_t);
        protect(res);
        SET_VECTOR_ELT(
            res,
            0 as libc::c_int as R_xlen_t,
            ScalarReal(ssq / nu as libc::c_double),
        );
        SET_VECTOR_ELT(res, 1 as libc::c_int as R_xlen_t, sResid);
        unprotect(2 as libc::c_int);
        return res;
    } else {
        unprotect(1 as libc::c_int);
        return ScalarReal(ssq / nu as libc::c_double);
    };
}
#[no_mangle]
pub unsafe extern "C" fn TSconv(mut a: SEXP, mut b: SEXP) -> SEXP {
    let mut na: libc::c_int = 0;
    let mut nb: libc::c_int = 0;
    let mut nab: libc::c_int = 0;
    let mut ab: SEXP = 0 as *mut SEXPREC;
    let mut ra: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut rb: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut rab: *mut libc::c_double = 0 as *mut libc::c_double;
    a = coerceVector(a, 14 as libc::c_int as SEXPTYPE);
    protect(a);
    b = coerceVector(b, 14 as libc::c_int as SEXPTYPE);
    protect(b);
    na = LENGTH(a);
    nb = LENGTH(b);
    nab = na + nb - 1 as libc::c_int;
    ab = allocVector(14 as libc::c_int as SEXPTYPE, nab as R_xlen_t);
    protect(ab);
    ra = REAL(a);
    rb = REAL(b);
    rab = REAL(ab);
    let mut i: libc::c_int = 0 as libc::c_int;
    while i < nab {
        *rab.offset(i as isize) = 0.0f64;
        i += 1
    }
    let mut i_0: libc::c_int = 0 as libc::c_int;
    while i_0 < na {
        let mut j: libc::c_int = 0 as libc::c_int;
        while j < nb {
            *rab.offset((i_0 + j) as isize) += *ra.offset(i_0 as isize) * *rb.offset(j as isize);
            j += 1
        }
        i_0 += 1
    }
    unprotect(3 as libc::c_int);
    return ab;
}
/* based on code from AS154 */
unsafe extern "C" fn inclu2(
    mut np: size_t,
    mut xnext: *mut libc::c_double,
    mut xrow: *mut libc::c_double,
    mut ynext: libc::c_double,
    mut d: *mut libc::c_double,
    mut rbar: *mut libc::c_double,
    mut thetab: *mut libc::c_double,
) {
    let mut cbar: libc::c_double = 0.;
    let mut sbar: libc::c_double = 0.;
    let mut di: libc::c_double = 0.;
    let mut xi: libc::c_double = 0.;
    let mut xk: libc::c_double = 0.;
    let mut rbthis: libc::c_double = 0.;
    let mut dpi: libc::c_double = 0.;
    let mut i: size_t = 0;
    let mut k: size_t = 0;
    let mut ithisr: size_t = 0;
    /*   This subroutine updates d, rbar, thetab by the inclusion
    of xnext and ynext. */
    i = 0 as libc::c_int as size_t;
    while i < np {
        *xrow.offset(i as isize) = *xnext.offset(i as isize);
        i = i.wrapping_add(1)
    }
    ithisr = 0 as libc::c_int as size_t;
    i = 0 as libc::c_int as size_t;
    while i < np {
        if *xrow.offset(i as isize) != 0.0f64 {
            xi = *xrow.offset(i as isize);
            di = *d.offset(i as isize);
            dpi = di + xi * xi;
            *d.offset(i as isize) = dpi;
            cbar = di / dpi;
            sbar = xi / dpi;
            k = i.wrapping_add(1 as libc::c_int as libc::c_ulong);
            while k < np {
                xk = *xrow.offset(k as isize);
                rbthis = *rbar.offset(ithisr as isize);
                *xrow.offset(k as isize) = xk - xi * rbthis;
                let fresh3 = ithisr;
                ithisr = ithisr.wrapping_add(1);
                *rbar.offset(fresh3 as isize) = cbar * rbthis + sbar * xk;
                k = k.wrapping_add(1)
            }
            xk = ynext;
            ynext = xk - xi * *thetab.offset(i as isize);
            *thetab.offset(i as isize) = cbar * *thetab.offset(i as isize) + sbar * xk;
            if di == 0.0f64 {
                return;
            }
        } else {
            ithisr = ithisr
                .wrapping_add(np)
                .wrapping_sub(i)
                .wrapping_sub(1 as libc::c_int as libc::c_ulong)
        }
        i = i.wrapping_add(1)
    }
}
/*
  Matwey V. Kornilov's implementation of algorithm by
  Dr. Raphael Rossignol
  See https://bugs.r-project.org/bugzilla3/show_bug.cgi?id=14682 for details.
*/
#[no_mangle]
pub unsafe extern "C" fn getQ0bis(mut sPhi: SEXP, mut sTheta: SEXP, mut sTol: SEXP) -> SEXP {
    let mut res: SEXP = 0 as *mut SEXPREC; // tol = REAL(sTol)[0];
    let mut p: libc::c_int = LENGTH(sPhi);
    let mut q: libc::c_int = LENGTH(sTheta);
    let mut phi: *mut libc::c_double = REAL(sPhi);
    let mut theta: *mut libc::c_double = REAL(sTheta);
    let mut i: libc::c_int = 0;
    let mut j: libc::c_int = 0;
    let mut r: libc::c_int = if p < q + 1 as libc::c_int {
        (q) + 1 as libc::c_int
    } else {
        p
    };
    /* Final result is block product
     *   Q0 = A1 SX A1^T + A1 SXZ A2^T + (A1 SXZ A2^T)^T + A2 A2^T ,
     * where A1 [i,j] = phi[i+j],
     *       A2 [i,j] = ttheta[i+j],  and SX, SXZ are defined below */
    res = allocMatrix(14 as libc::c_int as SEXPTYPE, r, r);
    protect(res);
    let mut P: *mut libc::c_double = REAL(res);
    /* Clean P */
    memset(
        P as *mut libc::c_void,
        0 as libc::c_int,
        ((r * r) as size_t).wrapping_mul(::std::mem::size_of::<libc::c_double>() as libc::c_ulong),
    );
    let mut ttheta: *mut libc::c_double = R_alloc(
        (q + 1 as libc::c_int) as size_t,
        ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
    ) as *mut libc::c_double;
    /* Init ttheta = c(1, theta) */
    *ttheta.offset(0 as libc::c_int as isize) = 1.0f64; // end if(p > 0)
    i = 1 as libc::c_int;
    while i < q + 1 as libc::c_int {
        *ttheta.offset(i as isize) = *theta.offset((i - 1 as libc::c_int) as isize);
        i += 1
    }
    if p > 0 as libc::c_int {
        let mut r2: libc::c_int = if p + q < p + 1 as libc::c_int {
            (p) + 1 as libc::c_int
        } else {
            (p) + q
        };
        let mut sgam: SEXP = protect(allocMatrix(14 as libc::c_int as SEXPTYPE, r2, r2));
        let mut sg: SEXP = protect(allocVector(14 as libc::c_int as SEXPTYPE, r2 as R_xlen_t));
        let mut gam: *mut libc::c_double = REAL(sgam);
        let mut g: *mut libc::c_double = REAL(sg);
        let mut tphi: *mut libc::c_double = R_alloc(
            (p + 1 as libc::c_int) as size_t,
            ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
        ) as *mut libc::c_double;
        /* Init tphi = c(1, -phi) */
        *tphi.offset(0 as libc::c_int as isize) = 1.0f64;
        i = 1 as libc::c_int;
        while i < p + 1 as libc::c_int {
            *tphi.offset(i as isize) = -*phi.offset((i - 1 as libc::c_int) as isize);
            i += 1
        }
        /* Compute the autocovariance function of U, the AR part of X */
        /* Gam := C1 + C2 ; initialize */
        memset(
            gam as *mut libc::c_void,
            0 as libc::c_int,
            ((r2 * r2) as size_t)
                .wrapping_mul(::std::mem::size_of::<libc::c_double>() as libc::c_ulong),
        );
        /* C1[E] */
        j = 0 as libc::c_int;
        while j < r2 {
            i = j;
            while i < r2 && i - j < p + 1 as libc::c_int {
                *gam.offset((j * r2 + i) as isize) += *tphi.offset((i - j) as isize);
                i += 1
            }
            j += 1
        }
        /* C2[E] */
        i = 0 as libc::c_int;
        while i < r2 {
            j = 1 as libc::c_int;
            while j < r2 && i + j < p + 1 as libc::c_int {
                *gam.offset((j * r2 + i) as isize) += *tphi.offset((i + j) as isize);
                j += 1
            }
            i += 1
        }
        /* Initialize g = (1 0 0 .... 0) */
        *g.offset(0 as libc::c_int as isize) = 1.0f64;
        i = 1 as libc::c_int;
        while i < r2 {
            *g.offset(i as isize) = 0.0f64;
            i += 1
        }
        /* rU = solve(Gam, g)  -> solve.default() -> .Internal(La_solve, .,)
         * --> fiddling with R-objects -> C and then F77_CALL(.) of dgesv, dlange, dgecon
         * FIXME: call these directly here, possibly even use 'info' instead of error(.)
         * e.g., in case of exact singularity.
         */
        let mut callS: SEXP = protect(lang4(
            install(b"solve.default\x00" as *const u8 as *const libc::c_char),
            sgam,
            sg,
            sTol,
        ));
        let mut su: SEXP = protect(eval(callS, R_BaseEnv));
        let mut u: *mut libc::c_double = REAL(su);
        /* SX = A SU A^T */
        /* A[i,j]  = ttheta[j-i] */
        /* SU[i,j] = u[abs(i-j)] */
        /* Q0 += ( A1 SX A1^T == A1 A SU A^T A1^T) */
        // (relying on good compiler optimization here:)
        i = 0 as libc::c_int;
        while i < r {
            j = i;
            while j < r {
                let mut k: libc::c_int = 0 as libc::c_int;
                while i + k < p {
                    let mut L: libc::c_int = k;
                    while L - k < q + 1 as libc::c_int {
                        let mut m: libc::c_int = 0 as libc::c_int;
                        while j + m < p {
                            let mut n: libc::c_int = m;
                            while n - m < q + 1 as libc::c_int {
                                *P.offset((r * i + j) as isize) += *phi.offset((i + k) as isize)
                                    * *phi.offset((j + m) as isize)
                                    * *ttheta.offset((L - k) as isize)
                                    * *ttheta.offset((n - m) as isize)
                                    * *u.offset((L - n).abs() as isize);
                                n += 1
                            }
                            m += 1
                        }
                        L += 1
                    }
                    k += 1
                }
                j += 1
            }
            i += 1
        }
        unprotect(4 as libc::c_int);
        /* Compute correlation matrix between X and Z */
        /* forwardsolve(C1, g) */
        /* C[i,j] = tphi[i-j] */
        /* g[i] = _ttheta(i) */
        let mut rrz: *mut libc::c_double = R_alloc(
            q as size_t,
            ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
        ) as *mut libc::c_double;
        if q > 0 as libc::c_int {
            i = 0 as libc::c_int;
            while i < q {
                *rrz.offset(i as isize) = *ttheta.offset(i as isize);
                j = if (0 as libc::c_int) < i - p {
                    (i) - p
                } else {
                    0 as libc::c_int
                };
                while j < i {
                    *rrz.offset(i as isize) -=
                        *rrz.offset(j as isize) * *tphi.offset((i - j) as isize);
                    j += 1
                }
                i += 1
            }
        }
        /* Q0 += A1 SXZ A2^T + (A1 SXZ A2^T)^T */
        /* SXZ[i,j] = rrz[j-i-1], j > 0 */
        i = 0 as libc::c_int;
        while i < r {
            j = i;
            while j < r {
                let mut k_0: libc::c_int = 0;
                let mut L_0: libc::c_int = 0;
                k_0 = 0 as libc::c_int;
                while i + k_0 < p {
                    L_0 = k_0 + 1 as libc::c_int;
                    while j + L_0 < q + 1 as libc::c_int {
                        *P.offset((r * i + j) as isize) += *phi.offset((i + k_0) as isize)
                            * *ttheta.offset((j + L_0) as isize)
                            * *rrz.offset((L_0 - k_0 - 1 as libc::c_int) as isize);
                        L_0 += 1
                    }
                    k_0 += 1
                }
                k_0 = 0 as libc::c_int;
                while j + k_0 < p {
                    L_0 = k_0 + 1 as libc::c_int;
                    while i + L_0 < q + 1 as libc::c_int {
                        *P.offset((r * i + j) as isize) += *phi.offset((j + k_0) as isize)
                            * *ttheta.offset((i + L_0) as isize)
                            * *rrz.offset((L_0 - k_0 - 1 as libc::c_int) as isize);
                        L_0 += 1
                    }
                    k_0 += 1
                }
                j += 1
            }
            i += 1
        }
    }
    /* Q0 += A2 A2^T */
    i = 0 as libc::c_int;
    while i < r {
        j = i;
        while j < r {
            let mut k_1: libc::c_int = 0 as libc::c_int;
            while j + k_1 < q + 1 as libc::c_int {
                *P.offset((r * i + j) as isize) +=
                    *ttheta.offset((i + k_1) as isize) * *ttheta.offset((j + k_1) as isize);
                k_1 += 1
            }
            j += 1
        }
        i += 1
    }
    /* Symmetrize result */
    i = 0 as libc::c_int;
    while i < r {
        j = i + 1 as libc::c_int;
        while j < r {
            *P.offset((r * j + i) as isize) = *P.offset((r * i + j) as isize);
            j += 1
        }
        i += 1
    }
    unprotect(1 as libc::c_int);
    return res;
}
#[no_mangle]
pub unsafe extern "C" fn getQ0(mut sPhi: SEXP, mut sTheta: SEXP) -> SEXP {
    let mut res: SEXP = 0 as *mut SEXPREC;
    let mut p: libc::c_int = LENGTH(sPhi);
    let mut q: libc::c_int = LENGTH(sTheta);
    let mut phi: *mut libc::c_double = REAL(sPhi);
    let mut theta: *mut libc::c_double = REAL(sTheta);
    /* thetab[np], xnext[np], xrow[np].  rbar[rbar] */
    /* NB: nrbar could overflow */
    let mut r: libc::c_int = if p < q + 1 as libc::c_int {
        (q) + 1 as libc::c_int
    } else {
        p
    };
    let mut np: size_t = (r * (r + 1 as libc::c_int) / 2 as libc::c_int) as size_t;
    let mut nrbar: size_t = np
        .wrapping_mul(np.wrapping_sub(1 as libc::c_int as libc::c_ulong))
        .wrapping_div(2 as libc::c_int as libc::c_ulong);
    let mut npr: size_t = 0;
    let mut npr1: size_t = 0;
    let mut indi: size_t = 0;
    let mut indj: size_t = 0;
    let mut indn: size_t = 0;
    let mut i: size_t = 0;
    let mut j: size_t = 0;
    let mut ithisr: size_t = 0;
    let mut ind: size_t = 0;
    let mut ind1: size_t = 0;
    let mut ind2: size_t = 0;
    let mut im: size_t = 0;
    let mut jm: size_t = 0;
    /* This is the limit using an int index.  We could use
    size_t and get more on a 64-bit system,
    but there seems no practical need. */
    if r > 350 as libc::c_int {
        error(dcgettext(
            b"stats\x00" as *const u8 as *const libc::c_char,
            b"maximum supported lag is 350\x00" as *const u8 as *const libc::c_char,
            5 as libc::c_int,
        )); // PR#16419
    }
    let mut xnext: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut xrow: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut rbar: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut thetab: *mut libc::c_double = 0 as *mut libc::c_double;
    let mut V: *mut libc::c_double = 0 as *mut libc::c_double;
    xnext = R_alloc(
        np,
        ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
    ) as *mut libc::c_double;
    xrow = R_alloc(
        np,
        ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
    ) as *mut libc::c_double;
    rbar = R_alloc(
        nrbar,
        ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
    ) as *mut libc::c_double;
    thetab = R_alloc(
        np,
        ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
    ) as *mut libc::c_double;
    V = R_alloc(
        np,
        ::std::mem::size_of::<libc::c_double>() as libc::c_ulong as libc::c_int,
    ) as *mut libc::c_double;
    ind = 0 as libc::c_int as size_t;
    j = 0 as libc::c_int as size_t;
    while j < r as libc::c_ulong {
        let mut vj: libc::c_double = 0.0f64;
        if j == 0 as libc::c_int as libc::c_ulong {
            vj = 1.0f64
        } else if j.wrapping_sub(1 as libc::c_int as libc::c_ulong) < q as libc::c_ulong {
            vj = *theta.offset(j.wrapping_sub(1 as libc::c_int as libc::c_ulong) as isize)
        }
        i = j;
        while i < r as libc::c_ulong {
            let mut vi: libc::c_double = 0.0f64;
            if i == 0 as libc::c_int as libc::c_ulong {
                vi = 1.0f64
            } else if i.wrapping_sub(1 as libc::c_int as libc::c_ulong) < q as libc::c_ulong {
                vi = *theta.offset(i.wrapping_sub(1 as libc::c_int as libc::c_ulong) as isize)
            }
            let fresh4 = ind;
            ind = ind.wrapping_add(1);
            *V.offset(fresh4 as isize) = vi * vj;
            i = i.wrapping_add(1)
        }
        j = j.wrapping_add(1)
    }
    res = allocMatrix(14 as libc::c_int as SEXPTYPE, r, r);
    protect(res);
    let mut P: *mut libc::c_double = REAL(res);
    if r == 1 as libc::c_int {
        if p == 0 as libc::c_int {
            *P.offset(0 as libc::c_int as isize) = 1.0f64
        } else {
            *P.offset(0 as libc::c_int as isize) = 1.0f64
                / (1.0f64
                    - *phi.offset(0 as libc::c_int as isize)
                        * *phi.offset(0 as libc::c_int as isize))
        }
        unprotect(1 as libc::c_int);
        return res;
    }
    if p > 0 as libc::c_int {
        /*      The set of equations s * vec(P0) = vec(v) is solved for
        vec(P0).  s is generated row by row in the array xnext.  The
        order of elements in P is changed, so as to bring more leading
        zeros into the rows of s. */
        i = 0 as libc::c_int as size_t;
        while i < nrbar {
            *rbar.offset(i as isize) = 0.0f64;
            i = i.wrapping_add(1)
        }
        i = 0 as libc::c_int as size_t;
        while i < np {
            *P.offset(i as isize) = 0.0f64;
            *thetab.offset(i as isize) = 0.0f64;
            *xnext.offset(i as isize) = 0.0f64;
            i = i.wrapping_add(1)
        }
        ind = 0 as libc::c_int as size_t;
        ind1 = -(1 as libc::c_int) as size_t;
        npr = np.wrapping_sub(r as libc::c_ulong);
        npr1 = npr.wrapping_add(1 as libc::c_int as libc::c_ulong);
        indj = npr;
        ind2 = npr.wrapping_sub(1 as libc::c_int as libc::c_ulong);
        j = 0 as libc::c_int as size_t;
        while j < r as libc::c_ulong {
            let mut phij: libc::c_double = if j < p as libc::c_ulong {
                *phi.offset(j as isize)
            } else {
                0.0f64
            };
            let fresh5 = indj;
            indj = indj.wrapping_add(1);
            *xnext.offset(fresh5 as isize) = 0.0f64;
            indi = npr1.wrapping_add(j);
            i = j;
            while i < r as libc::c_ulong {
                let fresh6 = ind;
                ind = ind.wrapping_add(1);
                let mut ynext: libc::c_double = *V.offset(fresh6 as isize);
                let mut phii: libc::c_double = if i < p as libc::c_ulong {
                    *phi.offset(i as isize)
                } else {
                    0.0f64
                };
                if j != (r - 1 as libc::c_int) as libc::c_ulong {
                    *xnext.offset(indj as isize) = -phii;
                    if i != (r - 1 as libc::c_int) as libc::c_ulong {
                        *xnext.offset(indi as isize) -= phij;
                        ind1 = ind1.wrapping_add(1);
                        *xnext.offset(ind1 as isize) = -1.0f64
                    }
                }
                *xnext.offset(npr as isize) = -phii * phij;
                ind2 = ind2.wrapping_add(1);
                if ind2 >= np {
                    ind2 = 0 as libc::c_int as size_t
                }
                *xnext.offset(ind2 as isize) += 1.0f64;
                inclu2(np, xnext, xrow, ynext, P, rbar, thetab);
                *xnext.offset(ind2 as isize) = 0.0f64;
                if i != (r - 1 as libc::c_int) as libc::c_ulong {
                    let fresh7 = indi;
                    indi = indi.wrapping_add(1);
                    *xnext.offset(fresh7 as isize) = 0.0f64;
                    *xnext.offset(ind1 as isize) = 0.0f64
                }
                i = i.wrapping_add(1)
            }
            j = j.wrapping_add(1)
        }
        ithisr = nrbar.wrapping_sub(1 as libc::c_int as libc::c_ulong);
        im = np.wrapping_sub(1 as libc::c_int as libc::c_ulong);
        i = 0 as libc::c_int as size_t;
        while i < np {
            let mut bi: libc::c_double = *thetab.offset(im as isize);
            jm = np.wrapping_sub(1 as libc::c_int as libc::c_ulong);
            j = 0 as libc::c_int as size_t;
            while j < i {
                let fresh8 = ithisr;
                ithisr = ithisr.wrapping_sub(1);
                let fresh9 = jm;
                jm = jm.wrapping_sub(1);
                bi -= *rbar.offset(fresh8 as isize) * *P.offset(fresh9 as isize);
                j = j.wrapping_add(1)
            }
            let fresh10 = im;
            im = im.wrapping_sub(1);
            *P.offset(fresh10 as isize) = bi;
            i = i.wrapping_add(1)
        }
        /*        now re-order p. */
        ind = npr;
        i = 0 as libc::c_int as size_t;
        while i < r as libc::c_ulong {
            let fresh11 = ind;
            ind = ind.wrapping_add(1);
            *xnext.offset(i as isize) = *P.offset(fresh11 as isize);
            i = i.wrapping_add(1)
        }
        ind = np.wrapping_sub(1 as libc::c_int as libc::c_ulong);
        ind1 = npr.wrapping_sub(1 as libc::c_int as libc::c_ulong);
        i = 0 as libc::c_int as size_t;
        while i < npr {
            let fresh12 = ind1;
            ind1 = ind1.wrapping_sub(1);
            let fresh13 = ind;
            ind = ind.wrapping_sub(1);
            *P.offset(fresh13 as isize) = *P.offset(fresh12 as isize);
            i = i.wrapping_add(1)
        }
        i = 0 as libc::c_int as size_t;
        while i < r as libc::c_ulong {
            *P.offset(i as isize) = *xnext.offset(i as isize);
            i = i.wrapping_add(1)
        }
    } else {
        /* P0 is obtained by backsubstitution for a moving average process. */
        indn = np;
        ind = np;
        i = 0 as libc::c_int as size_t;
        while i < r as libc::c_ulong {
            j = 0 as libc::c_int as size_t;
            while j <= i {
                ind = ind.wrapping_sub(1);
                *P.offset(ind as isize) = *V.offset(ind as isize);
                if j != 0 as libc::c_int as libc::c_ulong {
                    indn = indn.wrapping_sub(1);
                    *P.offset(ind as isize) += *P.offset(indn as isize)
                }
                j = j.wrapping_add(1)
            }
            i = i.wrapping_add(1)
        }
    }
    /* now unpack to a full matrix */
    i = (r - 1 as libc::c_int) as size_t;
    ind = np;
    while i > 0 as libc::c_int as libc::c_ulong {
        j = (r - 1 as libc::c_int) as size_t;
        while j >= i {
            ind = ind.wrapping_sub(1);
            *P.offset((r as libc::c_ulong).wrapping_mul(i).wrapping_add(j) as isize) =
                *P.offset(ind as isize);
            j = j.wrapping_sub(1)
        }
        i = i.wrapping_sub(1)
    }
    i = 0 as libc::c_int as size_t;
    while i < (r - 1 as libc::c_int) as libc::c_ulong {
        j = i.wrapping_add(1 as libc::c_int as libc::c_ulong);
        while j < r as libc::c_ulong {
            *P.offset(i.wrapping_add((r as libc::c_ulong).wrapping_mul(j)) as isize) =
                *P.offset(j.wrapping_add((r as libc::c_ulong).wrapping_mul(i)) as isize);
            j = j.wrapping_add(1)
        }
        i = i.wrapping_add(1)
    }
    unprotect(1 as libc::c_int);
    return res;
}