ParaMonte Fortran 2.0.0
Parallel Monte Carlo and Machine Learning Library
See the latest version documentation.
pm_sampleCov::setCovMeanMerged Interface Reference

Return the merged covariance and mean of a sample resulting from the merger of two separate (potentially weighted) samples \(A\) and \(B\).
More...

Detailed Description

Return the merged covariance and mean of a sample resulting from the merger of two separate (potentially weighted) samples \(A\) and \(B\).

See the documentation of pm_sampleCov for more information and definition online updating of sample covariance.

Note
The input and output variances of this generic interface are all biased variances.
A biased covariance can be readily unbiased by multiplying it with the appropriate bias-correction factors detailed in the documentation of pm_sampleCov.
Parameters
[out]cov: The output object of the same type and kind and rank and shape as covA, containing the biased covariance of the sample resulting form the merger of the two samples.
(optional. If missing, the resulting merged covariance will be written to the argument covB.)
[out]mean: The output object of the same type and kind and rank and shape as meanA, containing the mean of the sample resulting form the merger of the two samples.
(optional. If missing, the resulting merged mean will be written to the argument meanB.)
[in,out]covB: The input or input/output object of the same type and kind and rank and shape as the input argument covA, containing the biased covariance of the second sample that must be merged with the first sample.
If the input argument cov is missing, then covB contains the updated biased covariance of the merged sample on return.
Otherwise, the contents of covB remain intact upon return.
[in]meanB: The input or input/output object of the same type and kind and rank as the input argument covA of size size(covA, 1), containing the mean of the second sample that must be merged with the mean of the first sample meanA.
If the input argument mean is missing, then meanB contains the updated mean of the merged sample on return.
Otherwise, the contents of meanB remain intact upon return.
[in]covA: The input contiguous matrix of shape (1 : ndim, 1 : ndim) of,
  1. type complex of kind any supported by the processor (e.g., CK, CK32, CK64, or CK128),
  2. type real of kind any supported by the processor (e.g., RK, RK32, RK64, or RK128),
containing the biased covariance of the first sample that must be merged with the second sample.
[in]meanA: The input object of the same type and kind as the input argument covA of size size(covA, 1), containing the mean of the first sample that must be merged with the mean of the second sample meanB.
[in]fracA: The input scalar of type real of the same kind as covA, containing the ratio of the sum of the weights of all points in sample \(A\) to sum of weights of all points in the merged sample.
If the sample is unweighted, then fracA is simply size(sampleA) / (size(sampleA) + size(sampleB)).
[in]subset: The input scalar constant that can be:
  1. The constant lowDia, implying that only the lower-diagonal subset of the input cov, covB, and covA matrices must be accessed and/or manipulated.
  2. The constant uppDia, implying that only the upper-diagonal subset of the input cov, covB, and covA matrices must be accessed and/or manipulated.
This input argument is merely serves to resolve the different procedures of this generic interface from each other at compile-time.


Possible calling interfaces

call setCovMeanMerged(covB(1:ndim, 1:ndim), meanB(1:ndim), covA(1:ndim, 1:ndim), meanA(1:ndim), fracA, subset)
call setCovMeanMerged(cov(1:ndim, 1:ndim), mean(1:ndim), covB(1:ndim, 1:ndim), meanB(1:ndim), covA(1:ndim, 1:ndim), meanA(1:ndim), fracA, subset)
Return the merged covariance and mean of a sample resulting from the merger of two separate (potentia...
This module contains classes and procedures for computing the properties related to the covariance ma...
Warning
The condition 0 < fracA .and. fracA < 1 must hold for the corresponding input arguments.
The condition all(shape(covB) == shape(covA)) must hold for the corresponding input arguments.
The condition all(size(meanB) == shape(covA)) must hold for the corresponding input arguments.
The condition all(size(meanA) == shape(covA)) must hold for the corresponding input arguments.
The condition all(size(mean) == shape(covA)) must hold for the corresponding input arguments.
The condition all(shape(cov) == shape(covA)) must hold for the corresponding input arguments.
These conditions are verified only if the library is built with the preprocessor macro CHECK_ENABLED=1.
The pure procedure(s) documented herein become impure when the ParaMonte library is compiled with preprocessor macro CHECK_ENABLED=1.
By default, these procedures are pure in release build and impure in debug and testing builds.
See also
getCor
setCor
getCov
setCov
getVar
setVar
getMean
setMean
getCovMerged
setCovMerged
setCovMeanMerged
getVarCorrection
getMeanMerged
setMeanMerged
getVarMerged
setVarMerged


Example usage

1program example
2
3 use pm_kind, only: SK, IK
4 use pm_kind, only: TKG => RKS ! All other real types are also supported.
5 use pm_sampleCov, only: getCov
6 use pm_sampleMean, only: getMean
7 use pm_sampleCov, only: uppDia, lowDia
10 use pm_arrayResize, only: setResized
11 use pm_distUnif, only: getUnifRand
12 use pm_arrayRange, only: getRange
13 use pm_io, only: display_type
14
15 implicit none
16
17 integer(IK) , parameter :: nsam = 2
18 integer(IK) , allocatable :: iweight(:)
19 integer(IK) :: isam, ndim, lb(nsam), ub(nsam)
20 integer(IK) :: dim, itry, ntry = 10
21 type(display_type) :: disp
22 disp = display_type(file = "main.out.F90")
23
24 call disp%skip()
25 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
26 call disp%show("!Compute the biased merged covariance of a multivariate sample.")
27 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
28 call disp%skip()
29
30 block
31 real(TKG), allocatable :: mean(:,:), cov(:,:,:), meanMerged(:), covMerged(:,:)
32 real(TKG), allocatable :: sample(:,:)
33 do itry = 1, ntry
34 call disp%skip()
35 call disp%show("dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)")
36 dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
37 call disp%show("do isam = 2, nsam")
38 call disp%show(" lb(isam) = ub(isam - 1) + 1")
39 call disp%show(" ub(isam) = ub(isam - 1) + getUnifRand(2, 7)")
40 call disp%show("end do")
41 do isam = 2, nsam
42 lb(isam) = ub(isam - 1) + 1
43 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
44 end do
45 call disp%show("lb")
46 call disp%show( lb )
47 call disp%show("ub")
48 call disp%show( ub )
49 call disp%show("ndim = getUnifRand(1, minval(ub - lb + 1, 1))")
50 ndim = getUnifRand(1, minval(ub - lb + 1, 1))
51 call disp%show("call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])")
52 call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
53 call disp%show("call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])")
54 call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
55 call disp%show("call setResized(covMerged, [ndim, ndim])")
56 call setResized(covMerged, [ndim, ndim])
57 call disp%show("call setResized(meanMerged, ndim)")
58 call setResized(meanMerged, ndim)
59 call disp%show("sample = getUnifRand(-1., +1., ndim, ub(nsam))")
60 sample = getUnifRand(-1., +1., ndim, ub(nsam))
61 call disp%show("sample")
62 call disp%show( sample )
63 call disp%show("cov(:,:,0) = getCov(sample, dim)")
64 cov(:,:,0) = getCov(sample, dim)
65 call disp%show("cov(:,:,0) ! reference")
66 call disp%show( cov(:,:,0) )
67 call disp%show("mean(:,0) = getMean(sample, dim)")
68 mean(:,0) = getMean(sample, dim)
69 call disp%show("mean(:,0) ! reference")
70 call disp%show( mean(:,0) )
71 call disp%show("do isam = 1, nsam")
72 call disp%show(" cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), dim)")
73 call disp%show(" mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), dim)")
74 call disp%show("end do")
75 do isam = 1, nsam
76 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), dim)
77 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), dim)
78 end do
79 call disp%show("call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), uppDia)")
80 call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), uppDia)
81 call disp%show("covMerged")
82 call disp%show( covMerged )
83 call disp%show("call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), lowDia)")
84 call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), lowDia)
85 call disp%show("covMerged")
86 call disp%show( covMerged )
87 call disp%show("call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), uppDia)")
88 call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), uppDia)
89 call disp%show("cov(:,:,2)")
90 call disp%show( cov(:,:,2) )
91 call disp%show("cov(:,:,0) ! reference")
92 call disp%show( cov(:,:,0) )
93 call disp%skip()
94 end do
95 end block
96
97 call disp%skip()
98 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
99 call disp%show("!Compute the biased merged covariance of a frequency weighted multivariate sample.")
100 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
101 call disp%skip()
102
103 block
104 real(TKG), allocatable :: mean(:,:), cov(:,:,:), meanMerged(:), covMerged(:,:)
105 real(TKG), allocatable :: sample(:,:)
106 do itry = 1, ntry
107 call disp%skip()
108 call disp%show("dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)")
109 dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
110 call disp%show("do isam = 2, nsam")
111 call disp%show(" lb(isam) = ub(isam - 1) + 1")
112 call disp%show(" ub(isam) = ub(isam - 1) + getUnifRand(2, 7)")
113 call disp%show("end do")
114 do isam = 2, nsam
115 lb(isam) = ub(isam - 1) + 1
116 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
117 end do
118 call disp%show("lb")
119 call disp%show( lb )
120 call disp%show("ub")
121 call disp%show( ub )
122 call disp%show("ndim = getUnifRand(1, minval(ub - lb + 1, 1))")
123 ndim = getUnifRand(1, minval(ub - lb + 1, 1))
124 call disp%show("call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])")
125 call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
126 call disp%show("call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])")
127 call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
128 call disp%show("call setResized(covMerged, [ndim, ndim])")
129 call setResized(covMerged, [ndim, ndim])
130 call disp%show("call setResized(meanMerged, ndim)")
131 call setResized(meanMerged, ndim)
132 call disp%show("sample = getUnifRand(-1., +1., ndim, ub(nsam))")
133 sample = getUnifRand(-1., +1., ndim, ub(nsam))
134 call disp%show("sample")
135 call disp%show( sample )
136 call disp%show("iweight = getUnifRand(1, 10, size(sample, dim, IK))")
137 iweight = getUnifRand(1, 10, size(sample, dim, IK))
138 call disp%show("iweight")
139 call disp%show( iweight )
140 call disp%show("cov(:,:,0) = getCov(sample, dim, iweight)")
141 cov(:,:,0) = getCov(sample, dim, iweight)
142 call disp%show("cov(:,:,0) ! reference")
143 call disp%show( cov(:,:,0) )
144 call disp%show("mean(:,0) = getMean(sample, dim, iweight)")
145 mean(:,0) = getMean(sample, dim, iweight)
146 call disp%show("mean(:,0) ! reference")
147 call disp%show( mean(:,0) )
148 call disp%show("do isam = 1, nsam")
149 call disp%show(" cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), dim, iweight(lb(isam):ub(isam)))")
150 call disp%show(" mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), dim, iweight(lb(isam):ub(isam)))")
151 call disp%show("end do")
152 do isam = 1, nsam
153 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), dim, iweight(lb(isam):ub(isam)))
154 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), dim, iweight(lb(isam):ub(isam)))
155 end do
156 call disp%show("call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), uppDia)")
157 call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), uppDia)
158 call disp%show("covMerged")
159 call disp%show( covMerged )
160 call disp%show("call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), lowDia)")
161 call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), lowDia)
162 call disp%show("covMerged")
163 call disp%show( covMerged )
164 call disp%show("call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), uppDia)")
165 call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), uppDia)
166 call disp%show("cov(:,:,2)")
167 call disp%show( cov(:,:,2) )
168 call disp%show("cov(:,:,0) ! reference")
169 call disp%show( cov(:,:,0) )
170 call disp%skip()
171 end do
172 end block
173
174 call disp%skip()
175 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
176 call disp%show("!Compute the biased merged covariance of a reliability weighted multivariate sample.")
177 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
178 call disp%skip()
179
180 block
181 real(TKG), allocatable :: mean(:,:), cov(:,:,:), meanMerged(:), covMerged(:,:)
182 real(TKG), allocatable :: sample(:,:)
183 real(TKG), allocatable :: rweight(:)
184 do itry = 1, ntry
185 call disp%skip()
186 call disp%show("dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)")
187 dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
188 call disp%show("do isam = 2, nsam")
189 call disp%show(" lb(isam) = ub(isam - 1) + 1")
190 call disp%show(" ub(isam) = ub(isam - 1) + getUnifRand(2, 7)")
191 call disp%show("end do")
192 do isam = 2, nsam
193 lb(isam) = ub(isam - 1) + 1
194 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
195 end do
196 call disp%show("lb")
197 call disp%show( lb )
198 call disp%show("ub")
199 call disp%show( ub )
200 call disp%show("ndim = getUnifRand(1, minval(ub - lb + 1, 1))")
201 ndim = getUnifRand(1, minval(ub - lb + 1, 1))
202 call disp%show("call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])")
203 call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
204 call disp%show("call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])")
205 call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
206 call disp%show("call setResized(covMerged, [ndim, ndim])")
207 call setResized(covMerged, [ndim, ndim])
208 call disp%show("call setResized(meanMerged, ndim)")
209 call setResized(meanMerged, ndim)
210 call disp%show("sample = getUnifRand(-1., +1., ndim, ub(nsam))")
211 sample = getUnifRand(-1., +1., ndim, ub(nsam))
212 call disp%show("sample")
213 call disp%show( sample )
214 call disp%show("rweight = getUnifRand(1., 2., size(sample, dim, IK))")
215 rweight = getUnifRand(1., 2., size(sample, dim, IK))
216 call disp%show("rweight")
217 call disp%show( rweight )
218 call disp%show("cov(:,:,0) = getCov(sample, 2_IK, rweight)")
219 cov(:,:,0) = getCov(sample, 2_IK, rweight)
220 call disp%show("cov(:,:,0) ! reference")
221 call disp%show( cov(:,:,0) )
222 call disp%show("mean(:,0) = getMean(sample, dim, rweight)")
223 mean(:,0) = getMean(sample, dim, rweight)
224 call disp%show("mean(:,0) ! reference")
225 call disp%show( mean(:,0) )
226 call disp%show("do isam = 1, nsam")
227 call disp%show(" cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), 2_IK, rweight(lb(isam):ub(isam)))")
228 call disp%show(" mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), 2_IK, rweight(lb(isam):ub(isam)))")
229 call disp%show("end do")
230 do isam = 1, nsam
231 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), 2_IK, rweight(lb(isam):ub(isam)))
232 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), 2_IK, rweight(lb(isam):ub(isam)))
233 end do
234 call disp%show("call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), uppDia)")
235 call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), uppDia)
236 call disp%show("covMerged")
237 call disp%show( covMerged )
238 call disp%show("call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), lowDia)")
239 call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), lowDia)
240 call disp%show("covMerged")
241 call disp%show( covMerged )
242 call disp%show("call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), uppDia)")
243 call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), uppDia)
244 call disp%show("cov(:,:,2)")
245 call disp%show( cov(:,:,2) )
246 call disp%show("cov(:,:,0) ! reference")
247 call disp%show( cov(:,:,0) )
248 call disp%skip()
249 end do
250 end block
251
252end program example
Generate minimally-spaced character, integer, real sequences or sequences at fixed intervals of size ...
Resize (shrink or expand) an input allocatable array of rank 1..3 to arbitrary new lower and upper bo...
Allocate or resize (shrink or expand) an input allocatable scalar string or array of rank 1....
Generate and return a scalar or a contiguous array of rank 1 of length s1 of randomly uniformly distr...
This is a generic method of the derived type display_type with pass attribute.
Definition: pm_io.F90:11726
This is a generic method of the derived type display_type with pass attribute.
Definition: pm_io.F90:11508
Generate and return the (optionally unbiased) covariance matrix of a pair of (potentially weighted) t...
Generate and return the (weighted) mean of an input sample of nsam observations with ndim = 1 or 2 at...
This module contains procedures and generic interfaces for generating ranges of discrete character,...
This module contains procedures and generic interfaces for resizing allocatable arrays of various typ...
This module contains procedures and generic interfaces for resizing allocatable arrays of various typ...
This module contains classes and procedures for computing various statistical quantities related to t...
This module contains classes and procedures for input/output (IO) or generic display operations on st...
Definition: pm_io.F90:252
type(display_type) disp
This is a scalar module variable an object of type display_type for general display.
Definition: pm_io.F90:11393
This module defines the relevant Fortran kind type-parameters frequently used in the ParaMonte librar...
Definition: pm_kind.F90:268
integer, parameter IK
The default integer kind in the ParaMonte library: int32 in Fortran, c_int32_t in C-Fortran Interoper...
Definition: pm_kind.F90:540
integer, parameter SK
The default character kind in the ParaMonte library: kind("a") in Fortran, c_char in C-Fortran Intero...
Definition: pm_kind.F90:539
integer, parameter RKS
The single-precision real kind in Fortran mode. On most platforms, this is an 32-bit real kind.
Definition: pm_kind.F90:567
This module contains classes and procedures for computing the first moment (i.e., the statistical mea...
Generate and return an object of type display_type.
Definition: pm_io.F90:10282

Example Unix compile command via Intel ifort compiler
1#!/usr/bin/env sh
2rm main.exe
3ifort -fpp -standard-semantics -O3 -Wl,-rpath,../../../lib -I../../../inc main.F90 ../../../lib/libparamonte* -o main.exe
4./main.exe

Example Windows Batch compile command via Intel ifort compiler
1del main.exe
2set PATH=..\..\..\lib;%PATH%
3ifort /fpp /standard-semantics /O3 /I:..\..\..\include main.F90 ..\..\..\lib\libparamonte*.lib /exe:main.exe
4main.exe

Example Unix / MinGW compile command via GNU gfortran compiler
1#!/usr/bin/env sh
2rm main.exe
3gfortran -cpp -ffree-line-length-none -O3 -Wl,-rpath,../../../lib -I../../../inc main.F90 ../../../lib/libparamonte* -o main.exe
4./main.exe

Example output
1
2!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3!Compute the biased merged covariance of a multivariate sample.
4!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
5
6
7dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
8do isam = 2, nsam
9 lb(isam) = ub(isam - 1) + 1
10 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
11end do
12lb
13+1, +6
14ub
15+5, +7
16ndim = getUnifRand(1, minval(ub - lb + 1, 1))
17call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
18call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
19call setResized(covMerged, [ndim, ndim])
20call setResized(meanMerged, ndim)
21sample = getUnifRand(-1., +1., ndim, ub(nsam))
22sample
23+0.401394963, +0.979276299, +0.534193277, -0.450837612, -0.838486075, +0.965392590E-2, +0.893495798
24+0.513590097, -0.936435103, -0.955090642, -0.972116470, +0.523586631, -0.497346282, +0.180513144
25cov(:,:,0) = getCov(sample, dim)
26cov(:,:,0) ! reference
27+0.396623105, -0.853262097E-1
28-0.853262097E-1, +0.413675517
29mean(:,0) = getMean(sample, dim)
30mean(:,0) ! reference
31+0.218384385, -0.306185544
32do isam = 1, nsam
33 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), dim)
34 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), dim)
35end do
36call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), uppDia)
37covMerged
38+0.396623135, -0.853261948E-1
39+0.307542974E-40, +0.413675517
40call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), lowDia)
41covMerged
42+0.396623135, -0.853261948E-1
43-0.853261948E-1, +0.413675517
44call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), uppDia)
45cov(:,:,2)
46+0.396623135, -0.853261948E-1
47+0.149780139, +0.413675517
48cov(:,:,0) ! reference
49+0.396623105, -0.853262097E-1
50-0.853262097E-1, +0.413675517
51
52
53dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
54do isam = 2, nsam
55 lb(isam) = ub(isam - 1) + 1
56 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
57end do
58lb
59+1, +5
60ub
61+4, +8
62ndim = getUnifRand(1, minval(ub - lb + 1, 1))
63call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
64call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
65call setResized(covMerged, [ndim, ndim])
66call setResized(meanMerged, ndim)
67sample = getUnifRand(-1., +1., ndim, ub(nsam))
68sample
69-0.203351617, +0.236961603, +0.381265283, +0.188914657, -0.636925578, +0.796635389, -0.868143082, -0.752728581
70+0.822335243, -0.457805395, +0.132377625, +0.770981908, +0.220284700, -0.338171721E-1, +0.530756593, -0.149759889
71-0.495133638, +0.865869164, +0.375614882, +0.918408394, -0.147484064, -0.844699740, -0.100102186, -0.708728790
72-0.293231487, -0.166109800E-1, +0.600808501, -0.688944459, +0.862112999, -0.589067459, -0.767288208, -0.496743679
73cov(:,:,0) = getCov(sample, dim)
74cov(:,:,0) ! reference
75+0.318405479, -0.497720167E-1, +0.811716318E-1, +0.342581235E-2
76-0.497720167E-1, +0.178811550, +0.432243943E-2, -0.611906350E-1
77+0.811716318E-1, +0.432243943E-2, +0.403089046, +0.624237731E-1
78+0.342581235E-2, -0.611906350E-1, +0.624237731E-1, +0.325805992
79mean(:,0) = getMean(sample, dim)
80mean(:,0) ! reference
81-0.107171491, +0.229419202, -0.170319974E-1, -0.173620597
82do isam = 1, nsam
83 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), dim)
84 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), dim)
85end do
86call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), uppDia)
87covMerged
88+0.318405479, -0.497720242E-1, +0.811716467E-1, +0.342581421E-2
89-0.853261948E-1, +0.178811550, +0.432246178E-2, -0.611906350E-1
90+0.00000000, +0.00000000, +0.403089046, +0.624237657E-1
91+0.00000000, +0.00000000, +0.00000000, +0.325805932
92call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), lowDia)
93covMerged
94+0.318405479, -0.497720242E-1, +0.811716467E-1, +0.342581421E-2
95-0.497720242E-1, +0.178811550, +0.432246178E-2, -0.611906350E-1
96+0.811716467E-1, +0.432246178E-2, +0.403089046, +0.624237657E-1
97+0.342581421E-2, -0.611906350E-1, +0.624237657E-1, +0.325805932
98call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), uppDia)
99cov(:,:,2)
100+0.318405479, -0.497720242E-1, +0.811716467E-1, +0.342581421E-2
101-0.769999474E-1, +0.178811550, +0.432246178E-2, -0.611906350E-1
102-0.154122844, +0.761472434E-1, +0.403089046, +0.624237657E-1
103-0.850856304E-1, +0.439180434E-2, +0.882764161E-1, +0.325805932
104cov(:,:,0) ! reference
105+0.318405479, -0.497720167E-1, +0.811716318E-1, +0.342581235E-2
106-0.497720167E-1, +0.178811550, +0.432243943E-2, -0.611906350E-1
107+0.811716318E-1, +0.432243943E-2, +0.403089046, +0.624237731E-1
108+0.342581235E-2, -0.611906350E-1, +0.624237731E-1, +0.325805992
109
110
111dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
112do isam = 2, nsam
113 lb(isam) = ub(isam - 1) + 1
114 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
115end do
116lb
117+1, +7
118ub
119+6, +10
120ndim = getUnifRand(1, minval(ub - lb + 1, 1))
121call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
122call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
123call setResized(covMerged, [ndim, ndim])
124call setResized(meanMerged, ndim)
125sample = getUnifRand(-1., +1., ndim, ub(nsam))
126sample
127+0.624962926, -0.127657890, -0.970286727, -0.968180299, +0.668535471, +0.501106858, -0.853803635, +0.305868983, -0.810999513, -0.363536358
128-0.563328266, +0.157229662, -0.184266567E-1, +0.781577706, -0.996863008, -0.805929184, +0.193685293E-1, -0.669980168, +0.466191769, -0.538042665
129+0.680117488, -0.558356285, +0.272832036, +0.855279326, +0.558489203, +0.770297170, -0.625860453, +0.632045031, -0.495243192, +0.996196985
130cov(:,:,0) = getCov(sample, dim)
131cov(:,:,0) ! reference
132+0.419856757, -0.301753849, +0.154538617
133-0.301753849, +0.308248430, -0.151841283
134+0.154538617, -0.151841283, +0.356217146
135mean(:,0) = getMean(sample, dim)
136mean(:,0) ! reference
137-0.199399024, -0.216820240, +0.308579773
138do isam = 1, nsam
139 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), dim)
140 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), dim)
141end do
142call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), uppDia)
143covMerged
144+0.419856817, -0.301753879, +0.154538602
145-0.497720242E-1, +0.308248371, -0.151841313
146+0.811716467E-1, +0.432246178E-2, +0.356217146
147call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), lowDia)
148covMerged
149+0.419856817, -0.301753879, +0.154538602
150-0.301753879, +0.308248371, -0.151841313
151+0.154538602, -0.151841313, +0.356217146
152call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), uppDia)
153cov(:,:,2)
154+0.419856817, -0.301753879, +0.154538602
155-0.178762853, +0.308248371, -0.151841313
156+0.246388942, -0.277714342, +0.356217146
157cov(:,:,0) ! reference
158+0.419856757, -0.301753849, +0.154538617
159-0.301753849, +0.308248430, -0.151841283
160+0.154538617, -0.151841283, +0.356217146
161
162
163dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
164do isam = 2, nsam
165 lb(isam) = ub(isam - 1) + 1
166 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
167end do
168lb
169+1, +5
170ub
171+4, +9
172ndim = getUnifRand(1, minval(ub - lb + 1, 1))
173call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
174call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
175call setResized(covMerged, [ndim, ndim])
176call setResized(meanMerged, ndim)
177sample = getUnifRand(-1., +1., ndim, ub(nsam))
178sample
179+0.181960821, -0.686401606, +0.763722539, +0.147343755, -0.626788020, +0.541644573, +0.115449786, +0.175318241, -0.708830237
180+0.932233572, +0.178673267, +0.900178432, -0.703949094, -0.655135632, -0.518699408, +0.518129945, +0.346210957, +0.430200100E-1
181+0.765156150, -0.833746195, -0.107750058, +0.464287281, +0.309900045E-1, -0.961233974, -0.632779241, -0.907546401, -0.949697256
182cov(:,:,0) = getCov(sample, dim)
183cov(:,:,0) ! reference
184+0.260105550, +0.957356989E-1, +0.627730861E-1
185+0.957356989E-1, +0.352771342, +0.331181437E-1
186+0.627730861E-1, +0.331181437E-1, +0.385390908
187mean(:,0) = getMean(sample, dim)
188mean(:,0) ! reference
189-0.107311280E-1, +0.115629114, -0.348035514
190do isam = 1, nsam
191 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), dim)
192 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), dim)
193end do
194call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), uppDia)
195covMerged
196+0.260105610, +0.957356989E-1, +0.627730638E-1
197-0.301753879, +0.352771223, +0.331181213E-1
198+0.154538602, -0.151841313, +0.385390788
199call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), lowDia)
200covMerged
201+0.260105610, +0.957356989E-1, +0.627730638E-1
202+0.957356989E-1, +0.352771223, +0.331181213E-1
203+0.627730638E-1, +0.331181213E-1, +0.385390788
204call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), uppDia)
205cov(:,:,2)
206+0.260105610, +0.957356989E-1, +0.627730638E-1
207+0.385766290E-1, +0.352771223, +0.331181213E-1
208-0.886560902E-1, -0.773828030E-1, +0.385390788
209cov(:,:,0) ! reference
210+0.260105550, +0.957356989E-1, +0.627730861E-1
211+0.957356989E-1, +0.352771342, +0.331181437E-1
212+0.627730861E-1, +0.331181437E-1, +0.385390908
213
214
215dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
216do isam = 2, nsam
217 lb(isam) = ub(isam - 1) + 1
218 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
219end do
220lb
221+1, +8
222ub
223+7, +9
224ndim = getUnifRand(1, minval(ub - lb + 1, 1))
225call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
226call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
227call setResized(covMerged, [ndim, ndim])
228call setResized(meanMerged, ndim)
229sample = getUnifRand(-1., +1., ndim, ub(nsam))
230sample
231+0.496601105, +0.771960616, -0.545367956, +0.342359900, +0.116581798, -0.579060316, -0.746165395, -0.928222656, -0.196714878
232cov(:,:,0) = getCov(sample, dim)
233cov(:,:,0) ! reference
234+0.320497394
235mean(:,0) = getMean(sample, dim)
236mean(:,0) ! reference
237-0.140891984
238do isam = 1, nsam
239 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), dim)
240 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), dim)
241end do
242call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), uppDia)
243covMerged
244+0.320497364
245call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), lowDia)
246covMerged
247+0.320497364
248call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), uppDia)
249cov(:,:,2)
250+0.320497364
251cov(:,:,0) ! reference
252+0.320497394
253
254
255dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
256do isam = 2, nsam
257 lb(isam) = ub(isam - 1) + 1
258 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
259end do
260lb
261+1, +6
262ub
263+5, +7
264ndim = getUnifRand(1, minval(ub - lb + 1, 1))
265call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
266call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
267call setResized(covMerged, [ndim, ndim])
268call setResized(meanMerged, ndim)
269sample = getUnifRand(-1., +1., ndim, ub(nsam))
270sample
271-0.264775157, +0.259539127, +0.402713537, -0.878988385, -0.187849283, -0.149602890, +0.157950640
272-0.169064045, -0.792035818, +0.915624857, +0.210077763, -0.282754064, -0.894224644, -0.325002909
273cov(:,:,0) = getCov(sample, dim)
274cov(:,:,0) ! reference
275+0.156066045, +0.464977697E-2
276+0.464977697E-2, +0.324015439
277mean(:,0) = getMean(sample, dim)
278mean(:,0) ! reference
279-0.944303498E-1, -0.191054136
280do isam = 1, nsam
281 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), dim)
282 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), dim)
283end do
284call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), uppDia)
285covMerged
286+0.156066030, +0.464977697E-2
287+0.307542974E-40, +0.324015379
288call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), lowDia)
289covMerged
290+0.156066030, +0.464977697E-2
291+0.464977697E-2, +0.324015379
292call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), uppDia)
293cov(:,:,2)
294+0.156066030, +0.464977697E-2
295+0.437665395E-1, +0.324015379
296cov(:,:,0) ! reference
297+0.156066045, +0.464977697E-2
298+0.464977697E-2, +0.324015439
299
300
301dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
302do isam = 2, nsam
303 lb(isam) = ub(isam - 1) + 1
304 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
305end do
306lb
307+1, +8
308ub
309+7, +11
310ndim = getUnifRand(1, minval(ub - lb + 1, 1))
311call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
312call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
313call setResized(covMerged, [ndim, ndim])
314call setResized(meanMerged, ndim)
315sample = getUnifRand(-1., +1., ndim, ub(nsam))
316sample
317+0.680088997, +0.856828094, +0.618274808, +0.935130596, +0.841377854, +0.744893074, +0.379459620, -0.134975195, +0.127558351, +0.221639037, -0.991853237
318-0.470487833, +0.567016482, -0.510230660, -0.450091243, +0.720406771, +0.488170505, -0.115825653, -0.432596207, -0.758372426, -0.263102651, -0.746771097E-1
319cov(:,:,0) = getCov(sample, dim)
320cov(:,:,0) ! reference
321+0.296680480, +0.762073770E-1
322+0.762073770E-1, +0.223634332
323mean(:,0) = getMean(sample, dim)
324mean(:,0) ! reference
325+0.388947457, -0.118162736
326do isam = 1, nsam
327 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), dim)
328 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), dim)
329end do
330call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), uppDia)
331covMerged
332+0.296680480, +0.762073845E-1
333+0.464977697E-2, +0.223634273
334call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), lowDia)
335covMerged
336+0.296680480, +0.762073845E-1
337+0.762073845E-1, +0.223634273
338call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), uppDia)
339cov(:,:,2)
340+0.296680480, +0.762073845E-1
341-0.799481496E-1, +0.223634273
342cov(:,:,0) ! reference
343+0.296680480, +0.762073770E-1
344+0.762073770E-1, +0.223634332
345
346
347dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
348do isam = 2, nsam
349 lb(isam) = ub(isam - 1) + 1
350 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
351end do
352lb
353+1, +8
354ub
355+7, +13
356ndim = getUnifRand(1, minval(ub - lb + 1, 1))
357call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
358call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
359call setResized(covMerged, [ndim, ndim])
360call setResized(meanMerged, ndim)
361sample = getUnifRand(-1., +1., ndim, ub(nsam))
362sample
363+0.352897048, -0.625674367, +0.769135237, -0.107919097, +0.280082226, -0.278752685, -0.608349562, +0.600519776, +0.405497074, +0.741304159, -0.499057293, -0.325697422, -0.781980634
364-0.694901466, +0.871395946, +0.504343867, +0.553126574, -0.502991438, -0.261566043, -0.342728138, +0.821214318, +0.220376253E-1, +0.654630184, +0.924949408, -0.322317839, +0.324513316
365-0.978813648, +0.257770300, -0.579929829, -0.538145423, +0.966723800, -0.719702721, +0.267160177, +0.563987613, -0.191644073, +0.923226237, -0.569140911E-1, +0.117556810, -0.152100086
366cov(:,:,0) = getCov(sample, dim)
367cov(:,:,0) ! reference
368+0.283554226, +0.543091865E-2, +0.359909013E-1
369+0.543091865E-2, +0.300668955, +0.670931935E-1
370+0.359909013E-1, +0.670931935E-1, +0.340043515
371mean(:,0) = getMean(sample, dim)
372mean(:,0) ! reference
373-0.599965686E-2, +0.196285114, -0.929422583E-2
374do isam = 1, nsam
375 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), dim)
376 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), dim)
377end do
378call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), uppDia)
379covMerged
380+0.283554256, +0.543093588E-2, +0.359908566E-1
381+0.762073845E-1, +0.300668985, +0.670932829E-1
382+0.00000000, +0.331181213E-1, +0.340043724
383call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), lowDia)
384covMerged
385+0.283554256, +0.543093588E-2, +0.359908566E-1
386+0.543093588E-2, +0.300668985, +0.670932829E-1
387+0.359908566E-1, +0.670932829E-1, +0.340043724
388call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), uppDia)
389cov(:,:,2)
390+0.283554256, +0.543093588E-2, +0.359908566E-1
391+0.533605441E-1, +0.300668985, +0.670932829E-1
392+0.171034560, +0.727906153E-1, +0.340043724
393cov(:,:,0) ! reference
394+0.283554226, +0.543091865E-2, +0.359909013E-1
395+0.543091865E-2, +0.300668955, +0.670931935E-1
396+0.359909013E-1, +0.670931935E-1, +0.340043515
397
398
399dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
400do isam = 2, nsam
401 lb(isam) = ub(isam - 1) + 1
402 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
403end do
404lb
405+1, +5
406ub
407+4, +10
408ndim = getUnifRand(1, minval(ub - lb + 1, 1))
409call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
410call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
411call setResized(covMerged, [ndim, ndim])
412call setResized(meanMerged, ndim)
413sample = getUnifRand(-1., +1., ndim, ub(nsam))
414sample
415-0.220445395E-1, +0.101730824E-1, +0.319839120, -0.790473342, -0.979961276, +0.523224354, +0.361259222, -0.985275745, +0.806435108, +0.418053746
416+0.806536674E-1, -0.528206944, +0.816679001E-2, -0.585064888, -0.757138729, +0.746034622, -0.862583160, +0.560268998, -0.294501543, +0.642369270
417-0.357433915, +0.894447803, +0.903111100, -0.908998728, +0.777324438E-1, +0.873299837, +0.632227063, +0.111449003, -0.929121614, +0.984597325
418cov(:,:,0) = getCov(sample, dim)
419cov(:,:,0) ! reference
420+0.387672901, +0.724123493E-1, +0.126333743
421+0.724123493E-1, +0.321701169, +0.127962053
422+0.126333743, +0.127962053, +0.506277025
423mean(:,0) = getMean(sample, dim)
424mean(:,0) ! reference
425-0.338770263E-1, -0.990001932E-1, +0.228131041
426do isam = 1, nsam
427 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), dim)
428 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), dim)
429end do
430call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), uppDia)
431covMerged
432+0.387672842, +0.724123493E-1, +0.126333728
433+0.543093588E-2, +0.321701199, +0.127962068
434+0.359908566E-1, +0.670932829E-1, +0.506277204
435call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), lowDia)
436covMerged
437+0.387672842, +0.724123493E-1, +0.126333728
438+0.724123493E-1, +0.321701199, +0.127962068
439+0.126333728, +0.127962068, +0.506277204
440call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), uppDia)
441cov(:,:,2)
442+0.387672842, +0.724123493E-1, +0.126333728
443+0.498162918E-1, +0.321701199, +0.127962068
444+0.199596100E-1, +0.167634130, +0.506277204
445cov(:,:,0) ! reference
446+0.387672901, +0.724123493E-1, +0.126333743
447+0.724123493E-1, +0.321701169, +0.127962053
448+0.126333743, +0.127962053, +0.506277025
449
450
451dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
452do isam = 2, nsam
453 lb(isam) = ub(isam - 1) + 1
454 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
455end do
456lb
457+1, +3
458ub
459+2, +9
460ndim = getUnifRand(1, minval(ub - lb + 1, 1))
461call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
462call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
463call setResized(covMerged, [ndim, ndim])
464call setResized(meanMerged, ndim)
465sample = getUnifRand(-1., +1., ndim, ub(nsam))
466sample
467+0.403508544, +0.520493984, -0.125038624, +0.593883276, +0.493101239, -0.395837665, -0.303887606, +0.722774386, +0.105382919
468cov(:,:,0) = getCov(sample, dim)
469cov(:,:,0) ! reference
470+0.152988747
471mean(:,0) = getMean(sample, dim)
472mean(:,0) ! reference
473+0.223820046
474do isam = 1, nsam
475 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), dim)
476 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), dim)
477end do
478call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), uppDia)
479covMerged
480+0.152988762
481call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), lowDia)
482covMerged
483+0.152988762
484call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(ub(1), TKG) / real(ub(2), TKG), uppDia)
485cov(:,:,2)
486+0.152988762
487cov(:,:,0) ! reference
488+0.152988747
489
490
491!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
492!Compute the biased merged covariance of a frequency weighted multivariate sample.
493!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
494
495
496dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
497do isam = 2, nsam
498 lb(isam) = ub(isam - 1) + 1
499 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
500end do
501lb
502+1, +3
503ub
504+2, +9
505ndim = getUnifRand(1, minval(ub - lb + 1, 1))
506call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
507call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
508call setResized(covMerged, [ndim, ndim])
509call setResized(meanMerged, ndim)
510sample = getUnifRand(-1., +1., ndim, ub(nsam))
511sample
512+0.988867044, +0.862279415, +0.291793227, +0.295010805, +0.802412033, +0.140203595, +0.631700754, +0.820457101, +0.116583347
513iweight = getUnifRand(1, 10, size(sample, dim, IK))
514iweight
515+4, +1, +7, +10, +9, +1, +3, +8, +3
516cov(:,:,0) = getCov(sample, dim, iweight)
517cov(:,:,0) ! reference
518+0.844530016E-1
519mean(:,0) = getMean(sample, dim, iweight)
520mean(:,0) ! reference
521+0.564800620
522do isam = 1, nsam
523 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), dim, iweight(lb(isam):ub(isam)))
524 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), dim, iweight(lb(isam):ub(isam)))
525end do
526call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), uppDia)
527covMerged
528+0.844530538E-1
529call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), lowDia)
530covMerged
531+0.844530538E-1
532call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), uppDia)
533cov(:,:,2)
534+0.844530538E-1
535cov(:,:,0) ! reference
536+0.844530016E-1
537
538
539dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
540do isam = 2, nsam
541 lb(isam) = ub(isam - 1) + 1
542 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
543end do
544lb
545+1, +7
546ub
547+6, +11
548ndim = getUnifRand(1, minval(ub - lb + 1, 1))
549call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
550call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
551call setResized(covMerged, [ndim, ndim])
552call setResized(meanMerged, ndim)
553sample = getUnifRand(-1., +1., ndim, ub(nsam))
554sample
555-0.285781384, -0.375538111, +0.799287915, -0.833766460, -0.325526595, +0.790790319E-1, +0.122757912, -0.928408980, +0.960289836, +0.735855579, +0.305666804
556+0.518711209, -0.551659346, -0.326926589, +0.836314917, -0.141353846, +0.624673367, -0.825314760, -0.127631426E-1, +0.699443460, -0.769469857, +0.466434956E-1
557iweight = getUnifRand(1, 10, size(sample, dim, IK))
558iweight
559+10, +3, +6, +4, +3, +6, +9, +5, +7, +6, +9
560cov(:,:,0) = getCov(sample, dim, iweight)
561cov(:,:,0) ! reference
562+0.329561979, -0.641065612E-1
563-0.641065612E-1, +0.331481755
564mean(:,0) = getMean(sample, dim, iweight)
565mean(:,0) ! reference
566+0.107721582, +0.212829579E-1
567do isam = 1, nsam
568 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), dim, iweight(lb(isam):ub(isam)))
569 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), dim, iweight(lb(isam):ub(isam)))
570end do
571call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), uppDia)
572covMerged
573+0.329561949, -0.641065687E-1
574+0.307542974E-40, +0.331481785
575call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), lowDia)
576covMerged
577+0.329561949, -0.641065687E-1
578-0.641065687E-1, +0.331481785
579call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), uppDia)
580cov(:,:,2)
581+0.329561949, -0.641065687E-1
582+0.703652129E-1, +0.331481785
583cov(:,:,0) ! reference
584+0.329561979, -0.641065612E-1
585-0.641065612E-1, +0.331481755
586
587
588dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
589do isam = 2, nsam
590 lb(isam) = ub(isam - 1) + 1
591 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
592end do
593lb
594+1, +8
595ub
596+7, +10
597ndim = getUnifRand(1, minval(ub - lb + 1, 1))
598call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
599call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
600call setResized(covMerged, [ndim, ndim])
601call setResized(meanMerged, ndim)
602sample = getUnifRand(-1., +1., ndim, ub(nsam))
603sample
604-0.781103492, -0.193438172, -0.248659492, -0.353055000, -0.694214106, -0.232177138, +0.348670483, -0.623793602, -0.994449615, -0.746173382
605iweight = getUnifRand(1, 10, size(sample, dim, IK))
606iweight
607+9, +5, +9, +6, +3, +1, +2, +6, +10, +2
608cov(:,:,0) = getCov(sample, dim, iweight)
609cov(:,:,0) ! reference
610+0.113772787
611mean(:,0) = getMean(sample, dim, iweight)
612mean(:,0) ! reference
613-0.550008893
614do isam = 1, nsam
615 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), dim, iweight(lb(isam):ub(isam)))
616 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), dim, iweight(lb(isam):ub(isam)))
617end do
618call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), uppDia)
619covMerged
620+0.113772780
621call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), lowDia)
622covMerged
623+0.113772780
624call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), uppDia)
625cov(:,:,2)
626+0.113772780
627cov(:,:,0) ! reference
628+0.113772787
629
630
631dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
632do isam = 2, nsam
633 lb(isam) = ub(isam - 1) + 1
634 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
635end do
636lb
637+1, +5
638ub
639+4, +6
640ndim = getUnifRand(1, minval(ub - lb + 1, 1))
641call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
642call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
643call setResized(covMerged, [ndim, ndim])
644call setResized(meanMerged, ndim)
645sample = getUnifRand(-1., +1., ndim, ub(nsam))
646sample
647+0.785812497, -0.824784756, +0.666745663, +0.365077257, -0.185629249, +0.420706153
648iweight = getUnifRand(1, 10, size(sample, dim, IK))
649iweight
650+10, +4, +6, +6, +1, +1
651cov(:,:,0) = getCov(sample, dim, iweight)
652cov(:,:,0) ! reference
653+0.295174271
654mean(:,0) = getMean(sample, dim, iweight)
655mean(:,0) ! reference
656+0.392321467
657do isam = 1, nsam
658 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), dim, iweight(lb(isam):ub(isam)))
659 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), dim, iweight(lb(isam):ub(isam)))
660end do
661call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), uppDia)
662covMerged
663+0.295174271
664call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), lowDia)
665covMerged
666+0.295174271
667call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), uppDia)
668cov(:,:,2)
669+0.295174271
670cov(:,:,0) ! reference
671+0.295174271
672
673
674dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
675do isam = 2, nsam
676 lb(isam) = ub(isam - 1) + 1
677 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
678end do
679lb
680+1, +3
681ub
682+2, +6
683ndim = getUnifRand(1, minval(ub - lb + 1, 1))
684call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
685call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
686call setResized(covMerged, [ndim, ndim])
687call setResized(meanMerged, ndim)
688sample = getUnifRand(-1., +1., ndim, ub(nsam))
689sample
690+0.570910573, -0.127544522, +0.428713202, -0.257845879, -0.195810795, -0.793071866
691iweight = getUnifRand(1, 10, size(sample, dim, IK))
692iweight
693+4, +4, +7, +10, +3, +5
694cov(:,:,0) = getCov(sample, dim, iweight)
695cov(:,:,0) ! reference
696+0.194295719
697mean(:,0) = getMean(sample, dim, iweight)
698mean(:,0) ! reference
699-0.714179948E-1
700do isam = 1, nsam
701 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), dim, iweight(lb(isam):ub(isam)))
702 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), dim, iweight(lb(isam):ub(isam)))
703end do
704call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), uppDia)
705covMerged
706+0.194295719
707call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), lowDia)
708covMerged
709+0.194295719
710call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), uppDia)
711cov(:,:,2)
712+0.194295719
713cov(:,:,0) ! reference
714+0.194295719
715
716
717dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
718do isam = 2, nsam
719 lb(isam) = ub(isam - 1) + 1
720 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
721end do
722lb
723+1, +3
724ub
725+2, +9
726ndim = getUnifRand(1, minval(ub - lb + 1, 1))
727call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
728call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
729call setResized(covMerged, [ndim, ndim])
730call setResized(meanMerged, ndim)
731sample = getUnifRand(-1., +1., ndim, ub(nsam))
732sample
733+0.340297818, -0.314716816, +0.447024822, -0.253781915, +0.531564593, -0.360658050, -0.883807182, -0.783993602, +0.303907156
734iweight = getUnifRand(1, 10, size(sample, dim, IK))
735iweight
736+1, +6, +5, +3, +6, +2, +4, +7, +10
737cov(:,:,0) = getCov(sample, dim, iweight)
738cov(:,:,0) ! reference
739+0.270808876
740mean(:,0) = getMean(sample, dim, iweight)
741mean(:,0) ! reference
742-0.815969482E-1
743do isam = 1, nsam
744 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), dim, iweight(lb(isam):ub(isam)))
745 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), dim, iweight(lb(isam):ub(isam)))
746end do
747call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), uppDia)
748covMerged
749+0.270808846
750call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), lowDia)
751covMerged
752+0.270808846
753call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), uppDia)
754cov(:,:,2)
755+0.270808846
756cov(:,:,0) ! reference
757+0.270808876
758
759
760dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
761do isam = 2, nsam
762 lb(isam) = ub(isam - 1) + 1
763 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
764end do
765lb
766+1, +7
767ub
768+6, +10
769ndim = getUnifRand(1, minval(ub - lb + 1, 1))
770call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
771call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
772call setResized(covMerged, [ndim, ndim])
773call setResized(meanMerged, ndim)
774sample = getUnifRand(-1., +1., ndim, ub(nsam))
775sample
776+0.526983738, +0.734173536, -0.750278234E-1, +0.138120294, -0.859922290, -0.789736509, -0.469662428, +0.126173973, -0.336406827, -0.912330151E-1
777+0.595440626, +0.790780544, -0.110277057, +0.260060906, -0.881801963, -0.463621616E-1, +0.955803394, -0.987929702, -0.481496334, -0.146150231
778-0.400498033, +0.673115015, -0.780857801E-1, +0.560301185, +0.240870476, +0.764454842, +0.842976093, -0.204507709, +0.836968541, -0.985610485
779iweight = getUnifRand(1, 10, size(sample, dim, IK))
780iweight
781+6, +10, +8, +7, +8, +6, +2, +7, +10, +2
782cov(:,:,0) = getCov(sample, dim, iweight)
783cov(:,:,0) ! reference
784+0.277378827, +0.209649324, -0.492883772E-1
785+0.209649324, +0.391249955, +0.596357696E-1
786-0.492883772E-1, +0.596357696E-1, +0.239725173
787mean(:,0) = getMean(sample, dim, iweight)
788mean(:,0) ! reference
789-0.659118518E-1, -0.761377886E-1, +0.315032393
790do isam = 1, nsam
791 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), dim, iweight(lb(isam):ub(isam)))
792 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), dim, iweight(lb(isam):ub(isam)))
793end do
794call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), uppDia)
795covMerged
796+0.277378827, +0.209649414, -0.492884852E-1
797+0.307542974E-40, +0.391250044, +0.596356615E-1
798+0.00000000, +0.140129846E-44, +0.239725292
799call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), lowDia)
800covMerged
801+0.277378827, +0.209649414, -0.492884852E-1
802+0.209649414, +0.391250044, +0.596356615E-1
803-0.492884852E-1, +0.596356615E-1, +0.239725292
804call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), uppDia)
805cov(:,:,2)
806+0.277378827, +0.209649414, -0.492884852E-1
807-0.885012969E-1, +0.391250044, +0.596356615E-1
808-0.117471650, +0.118432455, +0.239725292
809cov(:,:,0) ! reference
810+0.277378827, +0.209649324, -0.492883772E-1
811+0.209649324, +0.391249955, +0.596357696E-1
812-0.492883772E-1, +0.596357696E-1, +0.239725173
813
814
815dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
816do isam = 2, nsam
817 lb(isam) = ub(isam - 1) + 1
818 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
819end do
820lb
821+1, +4
822ub
823+3, +7
824ndim = getUnifRand(1, minval(ub - lb + 1, 1))
825call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
826call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
827call setResized(covMerged, [ndim, ndim])
828call setResized(meanMerged, ndim)
829sample = getUnifRand(-1., +1., ndim, ub(nsam))
830sample
831-0.842973709, -0.172415376, +0.549086452, +0.757463813, -0.449250698, -0.330590129, -0.171786427
832-0.987430811E-1, -0.832129240, -0.561766267, +0.797159672, -0.800980210, -0.812227726, +0.451050162
833+0.260770321E-1, +0.187863469, +0.721985102, +0.766978621, -0.580928326, -0.888453722, +0.660766840
834iweight = getUnifRand(1, 10, size(sample, dim, IK))
835iweight
836+2, +10, +2, +5, +6, +1, +1
837cov(:,:,0) = getCov(sample, dim, iweight)
838cov(:,:,0) ! reference
839+0.238151833, +0.209738970, +0.196257070
840+0.209738970, +0.412076622, +0.218169644
841+0.196257070, +0.218169644, +0.264301926
842mean(:,0) = getMean(sample, dim, iweight)
843mean(:,0) ! reference
844-0.637959167E-1, -0.400873035, +0.129496112
845do isam = 1, nsam
846 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), dim, iweight(lb(isam):ub(isam)))
847 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), dim, iweight(lb(isam):ub(isam)))
848end do
849call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), uppDia)
850covMerged
851+0.238151759, +0.209738970, +0.196257055
852+0.209649414, +0.412076652, +0.218169659
853-0.492884852E-1, +0.596356615E-1, +0.264301926
854call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), lowDia)
855covMerged
856+0.238151759, +0.209738970, +0.196257055
857+0.209738970, +0.412076652, +0.218169659
858+0.196257055, +0.218169659, +0.264301926
859call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), uppDia)
860cov(:,:,2)
861+0.238151759, +0.209738970, +0.196257055
862+0.417133272, +0.412076652, +0.218169659
863+0.357336581, +0.529201031, +0.264301926
864cov(:,:,0) ! reference
865+0.238151833, +0.209738970, +0.196257070
866+0.209738970, +0.412076622, +0.218169644
867+0.196257070, +0.218169644, +0.264301926
868
869
870dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
871do isam = 2, nsam
872 lb(isam) = ub(isam - 1) + 1
873 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
874end do
875lb
876+1, +4
877ub
878+3, +10
879ndim = getUnifRand(1, minval(ub - lb + 1, 1))
880call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
881call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
882call setResized(covMerged, [ndim, ndim])
883call setResized(meanMerged, ndim)
884sample = getUnifRand(-1., +1., ndim, ub(nsam))
885sample
886-0.899042368, +0.415559888, +0.658540368, -0.353461504E-1, -0.461444497, +0.355355263, +0.992943525, -0.881807804, -0.572038889E-1, +0.793961287
887+0.375169516E-1, +0.421079159, +0.912170768, +0.649054885, +0.140879393, -0.322211623, +0.166822314, -0.508674026, -0.735155463, +0.518001437
888+0.701062679E-1, -0.610464692, +0.537359357, +0.243846416, -0.413061261, +0.626007318, -0.272456288, +0.724389315, -0.332440138E-1, -0.903913736
889iweight = getUnifRand(1, 10, size(sample, dim, IK))
890iweight
891+1, +3, +6, +5, +8, +7, +6, +2, +7, +9
892cov(:,:,0) = getCov(sample, dim, iweight)
893cov(:,:,0) ! reference
894+0.298760593, +0.118387401, -0.699651837E-1
895+0.118387401, +0.261046827, -0.596563146E-1
896-0.699651837E-1, -0.596563146E-1, +0.290581763
897mean(:,0) = getMean(sample, dim, iweight)
898mean(:,0) ! reference
899+0.256617576, +0.155372545, -0.887820199E-1
900do isam = 1, nsam
901 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), dim, iweight(lb(isam):ub(isam)))
902 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), dim, iweight(lb(isam):ub(isam)))
903end do
904call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), uppDia)
905covMerged
906+0.298760474, +0.118387401, -0.699651986E-1
907+0.209738970, +0.261046827, -0.596563220E-1
908+0.196257055, +0.218169659, +0.290581733
909call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), lowDia)
910covMerged
911+0.298760474, +0.118387401, -0.699651986E-1
912+0.118387401, +0.261046827, -0.596563220E-1
913-0.699651986E-1, -0.596563220E-1, +0.290581733
914call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), uppDia)
915cov(:,:,2)
916+0.298760474, +0.118387401, -0.699651986E-1
917+0.931684598E-1, +0.261046827, -0.596563220E-1
918-0.112460896, -0.134293810, +0.290581733
919cov(:,:,0) ! reference
920+0.298760593, +0.118387401, -0.699651837E-1
921+0.118387401, +0.261046827, -0.596563146E-1
922-0.699651837E-1, -0.596563146E-1, +0.290581763
923
924
925dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
926do isam = 2, nsam
927 lb(isam) = ub(isam - 1) + 1
928 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
929end do
930lb
931+1, +7
932ub
933+6, +8
934ndim = getUnifRand(1, minval(ub - lb + 1, 1))
935call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
936call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
937call setResized(covMerged, [ndim, ndim])
938call setResized(meanMerged, ndim)
939sample = getUnifRand(-1., +1., ndim, ub(nsam))
940sample
941+0.888599515, -0.585449457, +0.648303151, -0.978484273, +0.484613538, -0.912868023, -0.652709126, +0.463062048
942+0.143523693, +0.640499353, +0.379489541, -0.746734619, +0.141961813, +0.387100220, +0.328666449, -0.565282702
943iweight = getUnifRand(1, 10, size(sample, dim, IK))
944iweight
945+1, +3, +3, +5, +6, +1, +1, +5
946cov(:,:,0) = getCov(sample, dim, iweight)
947cov(:,:,0) ! reference
948+0.464220256, +0.767516866E-1
949+0.767516866E-1, +0.252795905
950mean(:,0) = getMean(sample, dim, iweight)
951mean(:,0) ! reference
952-0.631384831E-2, -0.715623572E-1
953do isam = 1, nsam
954 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), dim, iweight(lb(isam):ub(isam)))
955 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), dim, iweight(lb(isam):ub(isam)))
956end do
957call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), uppDia)
958covMerged
959+0.464220226, +0.767517015E-1
960+0.118387401, +0.252795935
961call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), lowDia)
962covMerged
963+0.464220226, +0.767517015E-1
964+0.767517015E-1, +0.252795935
965call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(iweight(:ub(1))), TKG) / real(sum(iweight), TKG), uppDia)
966cov(:,:,2)
967+0.464220226, +0.767517015E-1
968-0.138533592, +0.252795935
969cov(:,:,0) ! reference
970+0.464220256, +0.767516866E-1
971+0.767516866E-1, +0.252795905
972
973
974!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
975!Compute the biased merged covariance of a reliability weighted multivariate sample.
976!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
977
978
979dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
980do isam = 2, nsam
981 lb(isam) = ub(isam - 1) + 1
982 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
983end do
984lb
985+1, +3
986ub
987+2, +8
988ndim = getUnifRand(1, minval(ub - lb + 1, 1))
989call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
990call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
991call setResized(covMerged, [ndim, ndim])
992call setResized(meanMerged, ndim)
993sample = getUnifRand(-1., +1., ndim, ub(nsam))
994sample
995-0.571417451, +0.931302309, +0.977839231, -0.977189541, +0.124298334E-1, -0.929996610, -0.237989902, -0.665398240
996-0.721635222, -0.580464959, +0.696216345, +0.908992887, +0.697550774, +0.947196960, -0.743318677, +0.271038175
997rweight = getUnifRand(1., 2., size(sample, dim, IK))
998rweight
999+1.74241102, +1.07246280, +1.34209347, +1.74612665, +1.54444861, +1.87357593, +1.24614286, +1.18447518
1000cov(:,:,0) = getCov(sample, 2_IK, rweight)
1001cov(:,:,0) ! reference
1002+0.493098676, -0.114245765
1003-0.114245765, +0.498654991
1004mean(:,0) = getMean(sample, dim, rweight)
1005mean(:,0) ! reference
1006-0.272193044, +0.245787114
1007do isam = 1, nsam
1008 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), 2_IK, rweight(lb(isam):ub(isam)))
1009 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), 2_IK, rweight(lb(isam):ub(isam)))
1010end do
1011call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), uppDia)
1012covMerged
1013+0.493098617, -0.114245795
1014+0.307542974E-40, +0.498655081
1015call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), lowDia)
1016covMerged
1017+0.493098617, -0.114245795
1018-0.114245795, +0.498655081
1019call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), uppDia)
1020cov(:,:,2)
1021+0.493098617, -0.114245795
1022-0.625649914E-1, +0.498655081
1023cov(:,:,0) ! reference
1024+0.493098676, -0.114245765
1025-0.114245765, +0.498654991
1026
1027
1028dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
1029do isam = 2, nsam
1030 lb(isam) = ub(isam - 1) + 1
1031 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
1032end do
1033lb
1034+1, +4
1035ub
1036+3, +7
1037ndim = getUnifRand(1, minval(ub - lb + 1, 1))
1038call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
1039call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
1040call setResized(covMerged, [ndim, ndim])
1041call setResized(meanMerged, ndim)
1042sample = getUnifRand(-1., +1., ndim, ub(nsam))
1043sample
1044-0.100600958, +0.465979815, +0.937115550, +0.947117805E-2, -0.884447813, +0.688210726E-1, -0.206218958E-1
1045-0.982659101, -0.230036139, -0.986904025, -0.227110505, +0.249925733, +0.933954835, -0.398386002
1046rweight = getUnifRand(1., 2., size(sample, dim, IK))
1047rweight
1048+1.35797143, +1.48018575, +1.08264220, +1.86524141, +1.07807255, +1.06541657, +1.98629522
1049cov(:,:,0) = getCov(sample, 2_IK, rweight)
1050cov(:,:,0) ! reference
1051+0.210854188, -0.101127788
1052-0.101127788, +0.314686835
1053mean(:,0) = getMean(sample, dim, rweight)
1054mean(:,0) ! reference
1055+0.669851825E-1, -0.271669000
1056do isam = 1, nsam
1057 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), 2_IK, rweight(lb(isam):ub(isam)))
1058 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), 2_IK, rweight(lb(isam):ub(isam)))
1059end do
1060call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), uppDia)
1061covMerged
1062+0.210854143, -0.101127759
1063-0.114245795, +0.314686865
1064call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), lowDia)
1065covMerged
1066+0.210854143, -0.101127759
1067-0.101127759, +0.314686865
1068call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), uppDia)
1069cov(:,:,2)
1070+0.210854143, -0.101127759
1071-0.250287931E-1, +0.314686865
1072cov(:,:,0) ! reference
1073+0.210854188, -0.101127788
1074-0.101127788, +0.314686835
1075
1076
1077dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
1078do isam = 2, nsam
1079 lb(isam) = ub(isam - 1) + 1
1080 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
1081end do
1082lb
1083+1, +4
1084ub
1085+3, +7
1086ndim = getUnifRand(1, minval(ub - lb + 1, 1))
1087call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
1088call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
1089call setResized(covMerged, [ndim, ndim])
1090call setResized(meanMerged, ndim)
1091sample = getUnifRand(-1., +1., ndim, ub(nsam))
1092sample
1093+0.828327656, +0.342382431, -0.813050389, +0.176590443, +0.405714512, -0.922640324, +0.829164505
1094rweight = getUnifRand(1., 2., size(sample, dim, IK))
1095rweight
1096+1.14408064, +1.77327693, +1.16854942, +1.93906152, +1.79191124, +1.28518343, +1.46284866
1097cov(:,:,0) = getCov(sample, 2_IK, rweight)
1098cov(:,:,0) ! reference
1099+0.373551428
1100mean(:,0) = getMean(sample, dim, rweight)
1101mean(:,0) ! reference
1102+0.161035627
1103do isam = 1, nsam
1104 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), 2_IK, rweight(lb(isam):ub(isam)))
1105 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), 2_IK, rweight(lb(isam):ub(isam)))
1106end do
1107call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), uppDia)
1108covMerged
1109+0.373551577
1110call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), lowDia)
1111covMerged
1112+0.373551577
1113call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), uppDia)
1114cov(:,:,2)
1115+0.373551577
1116cov(:,:,0) ! reference
1117+0.373551428
1118
1119
1120dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
1121do isam = 2, nsam
1122 lb(isam) = ub(isam - 1) + 1
1123 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
1124end do
1125lb
1126+1, +4
1127ub
1128+3, +6
1129ndim = getUnifRand(1, minval(ub - lb + 1, 1))
1130call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
1131call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
1132call setResized(covMerged, [ndim, ndim])
1133call setResized(meanMerged, ndim)
1134sample = getUnifRand(-1., +1., ndim, ub(nsam))
1135sample
1136+0.888656020, +0.824515104, +0.657945395, -0.835144401, -0.136049867, +0.254159451
1137+0.439628601, +0.224233150, -0.341976762, -0.504526258, +0.470494628, +0.399492979E-1
1138+0.936748981E-1, +0.357176065E-1, -0.518700480, +0.749618888, +0.589431524, +0.201505542
1139rweight = getUnifRand(1., 2., size(sample, dim, IK))
1140rweight
1141+1.11516237, +1.84511876, +1.37388158, +1.87129498, +1.12207687, +1.71832609
1142cov(:,:,0) = getCov(sample, 2_IK, rweight)
1143cov(:,:,0) ! reference
1144+0.404780746, +0.130590633, -0.214089751
1145+0.130590633, +0.132172212, -0.106344931E-1
1146-0.214089751, -0.106344931E-1, +0.166609764
1147mean(:,0) = getMean(sample, dim, rweight)
1148mean(:,0) ! reference
1149+0.236299470, +0.957544707E-2, +0.206517741
1150do isam = 1, nsam
1151 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), 2_IK, rweight(lb(isam):ub(isam)))
1152 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), 2_IK, rweight(lb(isam):ub(isam)))
1153end do
1154call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), uppDia)
1155covMerged
1156+0.404780746, +0.130590633, -0.214089751
1157+0.307542974E-40, +0.132172182, -0.106344968E-1
1158+0.00000000, +0.140203595, +0.166609764
1159call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), lowDia)
1160covMerged
1161+0.404780746, +0.130590633, -0.214089751
1162+0.130590633, +0.132172182, -0.106344968E-1
1163-0.214089751, -0.106344968E-1, +0.166609764
1164call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), uppDia)
1165cov(:,:,2)
1166+0.404780746, +0.130590633, -0.214089751
1167+0.135784522, +0.132172182, -0.106344968E-1
1168-0.110218026, -0.434923284E-1, +0.166609764
1169cov(:,:,0) ! reference
1170+0.404780746, +0.130590633, -0.214089751
1171+0.130590633, +0.132172212, -0.106344931E-1
1172-0.214089751, -0.106344931E-1, +0.166609764
1173
1174
1175dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
1176do isam = 2, nsam
1177 lb(isam) = ub(isam - 1) + 1
1178 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
1179end do
1180lb
1181+1, +7
1182ub
1183+6, +13
1184ndim = getUnifRand(1, minval(ub - lb + 1, 1))
1185call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
1186call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
1187call setResized(covMerged, [ndim, ndim])
1188call setResized(meanMerged, ndim)
1189sample = getUnifRand(-1., +1., ndim, ub(nsam))
1190sample
1191+0.832262516, -0.127016664, +0.675074220, -0.692545176E-1, -0.312504053, -0.147016168, -0.656745791, -0.496482849, -0.112376213, -0.515277267, +0.206839681, +0.527799249, +0.666803122E-1
1192-0.189897180, +0.834278584, +0.197322369E-1, +0.380127549, +0.828397870, +0.599107385, -0.388283610, -0.889161468, +0.149561167E-1, +0.280150890, +0.230778813, +0.975948691, -0.463267922
1193-0.672663331, -0.887493372, +0.931710243, -0.195828080, -0.598317385, -0.292548776, +0.232463956, -0.939492345, -0.344697952, -0.279079795, -0.455149651, +0.930317283, +0.960207701
1194rweight = getUnifRand(1., 2., size(sample, dim, IK))
1195rweight
1196+1.17317104, +1.08978057, +1.03113854, +1.65376806, +1.15174007, +1.87458849, +1.07473564, +1.15723777, +1.65951502, +1.98242521, +1.76555133, +1.12642741, +1.52484775
1197cov(:,:,0) = getCov(sample, 2_IK, rweight)
1198cov(:,:,0) ! reference
1199+0.171400830, +0.239375029E-1, +0.849001184E-1
1200+0.239375029E-1, +0.252791524, -0.173558444E-1
1201+0.849001184E-1, -0.173558444E-1, +0.368183792
1202mean(:,0) = getMean(sample, dim, rweight)
1203mean(:,0) ! reference
1204-0.352020413E-1, +0.183240026, -0.142958015
1205do isam = 1, nsam
1206 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), 2_IK, rweight(lb(isam):ub(isam)))
1207 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), 2_IK, rweight(lb(isam):ub(isam)))
1208end do
1209call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), uppDia)
1210covMerged
1211+0.171400696, +0.239375569E-1, +0.849001408E-1
1212+0.130590633, +0.252791554, -0.173558611E-1
1213-0.214089751, -0.106344968E-1, +0.368183732
1214call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), lowDia)
1215covMerged
1216+0.171400696, +0.239375569E-1, +0.849001408E-1
1217+0.239375588E-1, +0.252791554, -0.173558611E-1
1218+0.849001408E-1, -0.173558593E-1, +0.368183732
1219call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), uppDia)
1220cov(:,:,2)
1221+0.171400696, +0.239375569E-1, +0.849001408E-1
1222+0.107305847, +0.252791554, -0.173558611E-1
1223+0.114116281, +0.839109495E-1, +0.368183732
1224cov(:,:,0) ! reference
1225+0.171400830, +0.239375029E-1, +0.849001184E-1
1226+0.239375029E-1, +0.252791524, -0.173558444E-1
1227+0.849001184E-1, -0.173558444E-1, +0.368183792
1228
1229
1230dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
1231do isam = 2, nsam
1232 lb(isam) = ub(isam - 1) + 1
1233 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
1234end do
1235lb
1236+1, +5
1237ub
1238+4, +11
1239ndim = getUnifRand(1, minval(ub - lb + 1, 1))
1240call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
1241call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
1242call setResized(covMerged, [ndim, ndim])
1243call setResized(meanMerged, ndim)
1244sample = getUnifRand(-1., +1., ndim, ub(nsam))
1245sample
1246+0.404199958, +0.385584474, +0.962077379, +0.989631295, +0.906902075, +0.369826794, +0.214111447, +0.719902635, +0.955797791, +0.859641910, +0.374192119
1247+0.916485667, +0.685104728, -0.771127343, -0.129209995, +0.691563368, -0.682598472, -0.353188515E-1, -0.270760059, +0.834522605, +0.212278366, -0.663943291E-1
1248-0.448723078, +0.341016531, -0.316381097, +0.161957860, -0.454055309, -0.789866090, -0.953052044E-1, -0.617180824, -0.371762156, +0.395909309, -0.237938404
1249rweight = getUnifRand(1., 2., size(sample, dim, IK))
1250rweight
1251+1.25719142, +1.88753462, +1.78973460, +1.63784313, +1.07069731, +1.45528889, +1.52669919, +1.99930871, +1.52823269, +1.76109016, +1.27821636
1252cov(:,:,0) = getCov(sample, 2_IK, rweight)
1253cov(:,:,0) ! reference
1254+0.798297375E-1, -0.112736505E-1, +0.807016809E-2
1255-0.112736505E-1, +0.312828720, +0.651230067E-1
1256+0.807016809E-2, +0.651230067E-1, +0.144341275
1257mean(:,0) = getMean(sample, dim, rweight)
1258mean(:,0) ! reference
1259+0.657695055, +0.913084820E-1, -0.198439956
1260do isam = 1, nsam
1261 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), 2_IK, rweight(lb(isam):ub(isam)))
1262 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), 2_IK, rweight(lb(isam):ub(isam)))
1263end do
1264call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), uppDia)
1265covMerged
1266+0.798297450E-1, -0.112736579E-1, +0.807017367E-2
1267+0.239375588E-1, +0.312828660, +0.651229993E-1
1268+0.849001408E-1, -0.173558593E-1, +0.144341305
1269call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), lowDia)
1270covMerged
1271+0.798297450E-1, -0.112736579E-1, +0.807017367E-2
1272-0.112736579E-1, +0.312828660, +0.651229993E-1
1273+0.807017367E-2, +0.651229993E-1, +0.144341305
1274call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), uppDia)
1275cov(:,:,2)
1276+0.798297450E-1, -0.112736579E-1, +0.807017367E-2
1277+0.901869610E-1, +0.312828660, +0.651229993E-1
1278+0.167856105E-1, +0.656217933E-1, +0.144341305
1279cov(:,:,0) ! reference
1280+0.798297375E-1, -0.112736505E-1, +0.807016809E-2
1281-0.112736505E-1, +0.312828720, +0.651230067E-1
1282+0.807016809E-2, +0.651230067E-1, +0.144341275
1283
1284
1285dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
1286do isam = 2, nsam
1287 lb(isam) = ub(isam - 1) + 1
1288 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
1289end do
1290lb
1291+1, +8
1292ub
1293+7, +13
1294ndim = getUnifRand(1, minval(ub - lb + 1, 1))
1295call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
1296call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
1297call setResized(covMerged, [ndim, ndim])
1298call setResized(meanMerged, ndim)
1299sample = getUnifRand(-1., +1., ndim, ub(nsam))
1300sample
1301-0.863177776E-1, +0.153170824, -0.867692471, +0.265655756, +0.258558869, +0.170848370E-1, +0.649551392, -0.763564944, -0.666157007, -0.504313231, +0.824491262, +0.474032640, +0.132736564
1302+0.615460396, +0.999525785E-1, +0.774540186, +0.651717782, +0.994927049, -0.955973148, +0.117524505, +0.375386596, -0.455659151, -0.404311299, +0.775979519, +0.570208192, +0.464315414
1303+0.677419186, -0.516756177, -0.764681458, +0.816069722, -0.223775864, -0.983417034E-2, -0.323994160, +0.495145798, +0.383114815, +0.322707653, +0.425106883, +0.324522257E-1, -0.161771059
1304+0.717651963, +0.416641474, -0.685375571, +0.877829790, -0.771937609, -0.419826031, +0.795935512, -0.852182627, +0.297023416, -0.382172465, -0.691636205, -0.291558981, +0.744650722
1305-0.110471249E-1, -0.543419719, -0.533681154, +0.893100619, -0.272954702, -0.570667624, -0.508905411, -0.869824409, -0.593479633, -0.367645383, -0.568282604E-1, +0.868855476, -0.526344299
1306-0.669735193, -0.908916593, +0.354645371, +0.552516103, +0.761848688, -0.593866110, +0.405617833, -0.139012694, +0.181772470, +0.539301634, -0.897852182E-1, -0.131992459, -0.417092800
1307rweight = getUnifRand(1., 2., size(sample, dim, IK))
1308rweight
1309+1.58328199, +1.16560173, +1.15800667, +1.47889853, +1.47997761, +1.28861511, +1.32251430, +1.47728014, +1.29625607, +1.83890510, +1.22764826, +1.77725387, +1.29353333
1310cov(:,:,0) = getCov(sample, 2_IK, rweight)
1311cov(:,:,0) ! reference
1312+0.262752771, +0.880049840E-1, -0.105997641E-1, +0.883777365E-1, +0.142275095, -0.193362106E-1
1313+0.880049840E-1, +0.301021129, -0.737055438E-2, -0.137093151E-1, +0.120418891, +0.496521890E-1
1314-0.105997641E-1, -0.737055438E-2, +0.198423773, +0.504211709E-1, +0.874293372E-1, +0.475427136E-2
1315+0.883777365E-1, -0.137093151E-1, +0.504211709E-1, +0.410423517, +0.821109265E-1, -0.693959594E-1
1316+0.142275095, +0.120418891, +0.874293372E-1, +0.821109265E-1, +0.299637049, +0.404402576E-1
1317-0.193362106E-1, +0.496521890E-1, +0.475427136E-2, -0.693959594E-1, +0.404402576E-1, +0.254095823
1318mean(:,0) = getMean(sample, dim, rweight)
1319mean(:,0) ! reference
1320-0.108147785E-1, +0.279249042, +0.120246388, -0.263803434E-1, -0.201071024, +0.784909632E-2
1321do isam = 1, nsam
1322 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), 2_IK, rweight(lb(isam):ub(isam)))
1323 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), 2_IK, rweight(lb(isam):ub(isam)))
1324end do
1325call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), uppDia)
1326covMerged
1327+0.262752801, +0.880050063E-1, -0.105997575E-1, +0.883777142E-1, +0.142275095, -0.193362311E-1
1328-0.112736579E-1, +0.301021159, -0.737057161E-2, -0.137093123E-1, +0.120418914, +0.496522300E-1
1329+0.807017367E-2, +0.651229993E-1, +0.198423848, +0.504211709E-1, +0.874293670E-1, +0.475430535E-2
1330+0.00000000, +0.00000000, +0.00000000, +0.410423577, +0.821109340E-1, -0.693959594E-1
1331+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.299637049, +0.404402614E-1
1332+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.254095703
1333call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), lowDia)
1334covMerged
1335+0.262752801, +0.880050063E-1, -0.105997575E-1, +0.883777142E-1, +0.142275095, -0.193362311E-1
1336+0.880050063E-1, +0.301021159, -0.737057161E-2, -0.137093123E-1, +0.120418914, +0.496522300E-1
1337-0.105997585E-1, -0.737057161E-2, +0.198423848, +0.504211709E-1, +0.874293670E-1, +0.475430535E-2
1338+0.883777142E-1, -0.137093123E-1, +0.504211709E-1, +0.410423577, +0.821109340E-1, -0.693959594E-1
1339+0.142275095, +0.120418914, +0.874293670E-1, +0.821109340E-1, +0.299637049, +0.404402614E-1
1340-0.193362311E-1, +0.496522300E-1, +0.475430535E-2, -0.693959594E-1, +0.404402614E-1, +0.254095703
1341call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), uppDia)
1342cov(:,:,2)
1343+0.262752801, +0.880050063E-1, -0.105997575E-1, +0.883777142E-1, +0.142275095, -0.193362311E-1
1344+0.209802270, +0.301021159, -0.737057161E-2, -0.137093123E-1, +0.120418914, +0.496522300E-1
1345-0.605140515E-1, -0.326006003E-1, +0.198423848, +0.504211709E-1, +0.874293670E-1, +0.475430535E-2
1346+0.482197665E-2, -0.531651005E-1, -0.840576440E-1, +0.410423577, +0.821109340E-1, -0.693959594E-1
1347+0.249660119, +0.116730630, -0.608433075E-1, -0.240780152E-1, +0.299637049, +0.404402614E-1
1348-0.853117481E-1, -0.121517599, +0.345167629E-1, -0.404614992E-1, -0.235077329E-1, +0.254095703
1349cov(:,:,0) ! reference
1350+0.262752771, +0.880049840E-1, -0.105997641E-1, +0.883777365E-1, +0.142275095, -0.193362106E-1
1351+0.880049840E-1, +0.301021129, -0.737055438E-2, -0.137093151E-1, +0.120418891, +0.496521890E-1
1352-0.105997641E-1, -0.737055438E-2, +0.198423773, +0.504211709E-1, +0.874293372E-1, +0.475427136E-2
1353+0.883777365E-1, -0.137093151E-1, +0.504211709E-1, +0.410423517, +0.821109265E-1, -0.693959594E-1
1354+0.142275095, +0.120418891, +0.874293372E-1, +0.821109265E-1, +0.299637049, +0.404402576E-1
1355-0.193362106E-1, +0.496521890E-1, +0.475427136E-2, -0.693959594E-1, +0.404402576E-1, +0.254095823
1356
1357
1358dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
1359do isam = 2, nsam
1360 lb(isam) = ub(isam - 1) + 1
1361 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
1362end do
1363lb
1364+1, +3
1365ub
1366+2, +9
1367ndim = getUnifRand(1, minval(ub - lb + 1, 1))
1368call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
1369call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
1370call setResized(covMerged, [ndim, ndim])
1371call setResized(meanMerged, ndim)
1372sample = getUnifRand(-1., +1., ndim, ub(nsam))
1373sample
1374+0.457989335, +0.810496807E-1, -0.654122233, +0.952951312, -0.175797582, -0.288653135, +0.610257149, -0.605512977, +0.164504528
1375-0.379562378E-3, +0.612921000, -0.233135819, +0.323578358, +0.360003829, +0.951012969, -0.248130918, +0.459004998, +0.377686262
1376rweight = getUnifRand(1., 2., size(sample, dim, IK))
1377rweight
1378+1.99604058, +1.03367066, +1.79838407, +1.52714586, +1.87341881, +1.16860247, +1.90648484, +1.74027216, +1.75829124
1379cov(:,:,0) = getCov(sample, 2_IK, rweight)
1380cov(:,:,0) ! reference
1381+0.274544597, -0.367375016E-1
1382-0.367375016E-1, +0.125712201
1383mean(:,0) = getMean(sample, dim, rweight)
1384mean(:,0) ! reference
1385+0.681750923E-1, +0.235321149
1386do isam = 1, nsam
1387 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), 2_IK, rweight(lb(isam):ub(isam)))
1388 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), 2_IK, rweight(lb(isam):ub(isam)))
1389end do
1390call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), uppDia)
1391covMerged
1392+0.274544597, -0.367375202E-1
1393+0.880050063E-1, +0.125712141
1394call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), lowDia)
1395covMerged
1396+0.274544597, -0.367375202E-1
1397-0.367375202E-1, +0.125712141
1398call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), uppDia)
1399cov(:,:,2)
1400+0.274544597, -0.367375202E-1
1401-0.305830371E-1, +0.125712141
1402cov(:,:,0) ! reference
1403+0.274544597, -0.367375016E-1
1404-0.367375016E-1, +0.125712201
1405
1406
1407dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
1408do isam = 2, nsam
1409 lb(isam) = ub(isam - 1) + 1
1410 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
1411end do
1412lb
1413+1, +8
1414ub
1415+7, +9
1416ndim = getUnifRand(1, minval(ub - lb + 1, 1))
1417call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
1418call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
1419call setResized(covMerged, [ndim, ndim])
1420call setResized(meanMerged, ndim)
1421sample = getUnifRand(-1., +1., ndim, ub(nsam))
1422sample
1423-0.165395021, +0.555638313, -0.414738774, +0.955900669, +0.472225070, +0.573829889, -0.314169765, +0.488227963, +0.341623902
1424rweight = getUnifRand(1., 2., size(sample, dim, IK))
1425rweight
1426+1.10429144, +1.66944027, +1.59734917, +1.46626925, +1.16047907, +1.60012817, +1.47653294, +1.06352961, +1.09223926
1427cov(:,:,0) = getCov(sample, 2_IK, rweight)
1428cov(:,:,0) ! reference
1429+0.207610339
1430mean(:,0) = getMean(sample, dim, rweight)
1431mean(:,0) ! reference
1432+0.276264995
1433do isam = 1, nsam
1434 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), 2_IK, rweight(lb(isam):ub(isam)))
1435 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), 2_IK, rweight(lb(isam):ub(isam)))
1436end do
1437call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), uppDia)
1438covMerged
1439+0.207610354
1440call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), lowDia)
1441covMerged
1442+0.207610354
1443call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), uppDia)
1444cov(:,:,2)
1445+0.207610354
1446cov(:,:,0) ! reference
1447+0.207610339
1448
1449
1450dim = 2; lb(1) = 1; ub(1) = getUnifRand(2, 7)
1451do isam = 2, nsam
1452 lb(isam) = ub(isam - 1) + 1
1453 ub(isam) = ub(isam - 1) + getUnifRand(2, 7)
1454end do
1455lb
1456+1, +8
1457ub
1458+7, +12
1459ndim = getUnifRand(1, minval(ub - lb + 1, 1))
1460call setRebound(cov, [1_IK, 1_IK, 0_IK], [ndim, ndim, nsam])
1461call setRebound(mean, [1_IK, 0_IK], [ndim, nsam])
1462call setResized(covMerged, [ndim, ndim])
1463call setResized(meanMerged, ndim)
1464sample = getUnifRand(-1., +1., ndim, ub(nsam))
1465sample
1466+0.661407113, +0.921583056, +0.328138232, -0.931953192, -0.989997864, -0.397196651, +0.530970097E-1, +0.857881069, +0.103220701, -0.173329234, +0.921910048, -0.611423254E-1
1467+0.162328124, -0.912275314, +0.334692001E-1, -0.149823785, +0.845754743, +0.403283119, +0.381153345, +0.251645207, +0.405446649, -0.707323909, +0.597941875E-2, +0.683746815
1468-0.376863718, -0.151410222, -0.118592024, -0.125443101, -0.994588852, -0.908079505, -0.937530160, -0.672951102, -0.983699679, -0.565156698, +0.949061751, +0.775063634
1469+0.313614607, +0.163513303, +0.915796041, +0.971271992E-1, -0.107908368, -0.381834149, +0.244665146, -0.621290803, -0.405495882, -0.538779974, -0.101119280E-1, +0.921919703
1470rweight = getUnifRand(1., 2., size(sample, dim, IK))
1471rweight
1472+1.83676660, +1.25321245, +1.51819849, +1.94060791, +1.59117222, +1.13043523, +1.51070869, +1.50350833, +1.81314301, +1.87763572, +1.20856524, +1.23270035
1473cov(:,:,0) = getCov(sample, 2_IK, rweight)
1474cov(:,:,0) ! reference
1475+0.411606401, -0.870113969E-1, +0.112263426, +0.210920442E-1
1476-0.870113969E-1, +0.238117605, -0.557373241E-1, +0.350888185E-1
1477+0.112263426, -0.557373241E-1, +0.342847049, +0.151569068
1478+0.210920442E-1, +0.350888185E-1, +0.151569068, +0.229862705
1479mean(:,0) = getMean(sample, dim, rweight)
1480mean(:,0) ! reference
1481+0.708988383E-1, +0.104679950, -0.384704679, +0.309142582E-1
1482do isam = 1, nsam
1483 cov(:,:,isam) = getCov(sample(:,lb(isam):ub(isam)), 2_IK, rweight(lb(isam):ub(isam)))
1484 mean(:,isam) = getMean(sample(:,lb(isam):ub(isam)), 2_IK, rweight(lb(isam):ub(isam)))
1485end do
1486call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), uppDia)
1487covMerged
1488+0.411606431, -0.870114118E-1, +0.112263419, +0.210920833E-1
1489+0.307542974E-40, +0.238117635, -0.557373203E-1, +0.350888222E-1
1490+0.00000000, +0.731962323E-1, +0.342847049, +0.151569054
1491+0.00000000, +0.353250235, +0.704036851E-3, +0.229862720
1492call setCovMeanMerged(covMerged, meanMerged, cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), lowDia)
1493covMerged
1494+0.411606431, -0.870114118E-1, +0.112263419, +0.210920833E-1
1495-0.870114118E-1, +0.238117635, -0.557373203E-1, +0.350888222E-1
1496+0.112263419, -0.557373203E-1, +0.342847049, +0.151569054
1497+0.210920833E-1, +0.350888222E-1, +0.151569054, +0.229862720
1498call setCovMeanMerged(cov(:,:,2), mean(:,2), cov(:,:,1), mean(:,1), real(sum(rweight(:ub(1))), TKG) / real(sum(rweight), TKG), uppDia)
1499cov(:,:,2)
1500+0.411606431, -0.870114118E-1, +0.112263419, +0.210920833E-1
1501+0.528442077E-1, +0.238117635, -0.557373203E-1, +0.350888222E-1
1502+0.830376446E-1, +0.758205354E-1, +0.342847049, +0.151569054
1503-0.440159403E-1, +0.142612323, +0.318926066, +0.229862720
1504cov(:,:,0) ! reference
1505+0.411606401, -0.870113969E-1, +0.112263426, +0.210920442E-1
1506-0.870113969E-1, +0.238117605, -0.557373241E-1, +0.350888185E-1
1507+0.112263426, -0.557373241E-1, +0.342847049, +0.151569068
1508+0.210920442E-1, +0.350888185E-1, +0.151569068, +0.229862705
1509
1510
Test:
test_pm_sampleCov


Final Remarks


If you believe this algorithm or its documentation can be improved, we appreciate your contribution and help to edit this page's documentation and source file on GitHub.
For details on the naming abbreviations, see this page.
For details on the naming conventions, see this page.
This software is distributed under the MIT license with additional terms outlined below.

  1. If you use any parts or concepts from this library to any extent, please acknowledge the usage by citing the relevant publications of the ParaMonte library.
  2. If you regenerate any parts/ideas from this library in a programming environment other than those currently supported by this ParaMonte library (i.e., other than C, C++, Fortran, MATLAB, Python, R), please also ask the end users to cite this original ParaMonte library.

This software is available to the public under a highly permissive license.
Help us justify its continued development and maintenance by acknowledging its benefit to society, distributing it, and contributing to it.

Author:
Amir Shahmoradi, April 21, 2017, 1:54 AM, Institute for Computational Engineering and Sciences (ICES), The University of Texas Austin

Definition at line 7883 of file pm_sampleCov.F90.


The documentation for this interface was generated from the following file: