1C -------------------------------------------------------------------- 2C MISSING Simulates missing observations 3C Developed by A.Rossi, C.Planas and G.Fiorentini 4C 5C Copyright (C) 2010-2014 European Commission 6C 7C This file is part of Program DMM 8C 9C DMM is free software developed at the Joint Research Centre of the 10C European Commission: you can redistribute it and/or modify it under 11C the terms of the GNU General Public License as published by 12C the Free Software Foundation, either version 3 of the License, or 13C (at your option) any later version. 14C 15C DMM is distributed in the hope that it will be useful, 16C but WITHOUT ANY WARRANTY; without even the implied warranty of 17C MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 18C GNU General Public License for more details. 19C 20C You should have received a copy of the GNU General Public License 21C along with DMM. If not, see <http://www.gnu.org/licenses/>. 22C -------------------------------------------------------------------- 23 SUBROUTINE MISSING(yk,ny,nz,nx,nu,ns,nt,nmis,theta,S,STATE,pdll,ykmis) 24 25 USE dfwin 26 INTERFACE 27 SUBROUTINE DESIGN(ny,nz,nx,nu,ns,nt,theta,c,H,G,a,F,R) 28 INTEGER ny,nz,nx,nu,ns(6),nt 29 DOUBLE PRECISION theta(nt) 30 DOUBLE PRECISION c(ny,max(1,nz),ns(1)),H(ny,nx,ns(2)), 31 1 G(ny,nu,ns(3)),a(nx,ns(4)),F(nx,nx,ns(5)),R(nx,nu,ns(6)) 32 END SUBROUTINE 33 END INTERFACE 34 CHARACTER*1 fittizia 35 POINTER (pdll,fittizia) 36 POINTER (pdesign,DESIGN) 37 38C INPUT 39 INTEGER ny,nz,nx,nu,nt,ns(6),nmis 40 DOUBLE PRECISION yk(ny+nz),theta(nt),STATE(nx) 41 42C OUTPUT 43 DOUBLE PRECISION ykmis(nmis) 44C LOCALS 45 INTEGER S(6),I,J,K,IFAIL,NIYK(ny) 46 DOUBLE PRECISION U(nu) 47 DOUBLE PRECISION gennor 48 DOUBLE PRECISION,ALLOCATABLE::R(:,:,:),c(:,:,:),H(:,:,:), 49 1 G(:,:,:),a(:,:),F(:,:,:) 50 51 ALLOCATE(R(nx,nu,ns(6)),c(ny,max(nz,1),ns(1)),H(ny,nx,ns(2)), 52 1 G(ny,nu,ns(3)),a(nx,ns(4)),F(nx,nx,ns(5))) 53 pdesign = getprocaddress(pdll, "design_"C) 54 CALL DESIGN(ny,nz,nx,nu,ns,nt,theta,c,H,G,a,F,R) 55 56C NIYK = not(IYK) 57 K = 0 58 DO 10 J = 1,ny 59 IF(yk(J).EQ.-99999.D0) THEN 60 K = K+1 61 NIYK(K) = J 6210 ENDIF 63 64C SAMPLING U 65 IFAIL = -1 66 U(1:nu) = 0.D0 67 DO 20 I = 1,nu 68c CALL G05EAF(U(I),1,1.D0,1,1.D-14,WORKU,3,IFAIL) 69c20 CALL G05EZF(U(I),1,WORKU,3,IFAIL) 7020 U(I) = gennor(0.D0,1.D0) 71 72 73C DRAW yk ~ f(yk|x,S,zk,theta) 74 DO 30 I = 1,nmis 7530 ykmis(I) = SUM(c(NIYK(I),1:nz,S(1))*yk(ny+1:ny+nz)) 76 + + SUM(H(NIYK(I),1:nx,S(2))*STATE(1:nx)) 77 + + SUM(G(NIYK(I),1:nu,S(3))*U(1:nu)) 78 79 DEALLOCATE (R,c,H,G,a,F) 80 81 RETURN 82 END 83