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

Return an affine-transformation of the input sample of shape (1:ndim) or (1:ndim, 1:nsam) or (1:nsam, 1:ndim) based on the specified values for the translation tlate and tranformations tform along the specified axis dim.
More...

Detailed Description

Return an affine-transformation of the input sample of shape (1:ndim) or (1:ndim, 1:nsam) or (1:nsam, 1:ndim) based on the specified values for the translation tlate and tranformations tform along the specified axis dim.

Here, ndim stands for the number of dimensions (data attributes) of the input sample and nsam represents the number of data points in the sample.
If the input tlate is the negative of the mean of the sample, then the returned sample will have a mean of zero.

Parameters
[out]affinity: The output object of the same type and kind and rank and shape as sample, containing the affine-transformed sample.
[in]sample: The input contiguous array of shape (nsam), (ndim, nsam), or (nsam, ndim) of,
  1. type real of kind any supported by the processor (e.g., RK, RK32, RK64, or RK128),
containing the sample to be affine-transformed.
[in]dim: The input scalar of type integer of default kind IK, whose value represents the dimension of the input sample containing different nsam observations:
  1. If dim = 1, the input sample is assumed to have the shape (1:nsam, 1:ndim).
  2. If dim = 2, the input sample is assumed to have the shape (1:ndim, 1:nsam).
(optional. It must be present if and only if the input arguments the condition rank(sample) > 1 holds.)
[in]tform: The input contiguous matrix of the same type and kind as sample of shape (1:ndim, 1:ndim), containing the transformation matrix of the affine-transformation.
[in]class: The input scalar constant that can be any of the following:
  1. The scalar constant genrecmat or equivalently an object of type genrecmat_type signifying the general rectangular matrix as the class of the input transformation matrix.
    All elements of the input matrix will be used in constructing the output sample.
  2. The scalar constant upperDiag or equivalently an object of type upperDiag_type signifying the upper-diagonal matrix as the class of the input transformation matrix.
    Only the upper-diagonal elements of tform will be used in constructing the output sample.
    The lower triangle elements of tform are assumed to be zero.
  3. The scalar constant lowerDiag or equivalently an object of type lowerDiag_type signifying the lower-diagonal matrix as the class of the input transformation matrix.
    Only the lower-diagonal elements of tform will be used in constructing the output sample.
    The upper triangle elements of tform are assumed to be zero.
  4. The scalar constant upperUnit or equivalently an object of type upperUnit_type signifying the unit-triangular matrix as the class of the input transformation matrix.
    Only the upper elements of tform will be used in constructing the output sample.
    The lower triangle elements of tform are assumed to be zero.
    The diagonal elements of tform are assumed to be one.
    The diagonal elements of tform are assumed to be zero.
  5. The scalar constant lowerUnit or equivalently an object of type lowerUnit_type signifying the unit-triangular matrix as the class of the input transformation matrix.
    Only the lower elements of tform will be used in constructing the output sample.
    The upper triangle elements of tform are assumed to be zero.
    The diagonal elements of tform are assumed to be zero.
[in]tlate: The input contiguous vector of the same type and kind as sample of shape (1:ndim), containing the amount by which the input sample must be translated (i.e., shifted).
  1. If the input argument dim = 1 then, size(tlate) == size(sample, 2) == ndim must hold.
  2. If the input argument dim = 2 then, size(tlate) == size(sample, 1) == ndim must hold.
(optional, default = [(0, i = 1, ndim)].)


Possible calling interfaces

call setAffinity(affinity(1:ndim), sample(1:ndim), tform(1:ndim, 1:ndim), class)
call setAffinity(affinity(1:ndim), sample(1:ndim), tform(1:ndim, 1:ndim), class, tlate(1:ndim))
call setAffinity(affinity(@shape(sample)), sample(:,:), dim, tform(1:size(sample, 3 - dim), 1:size(sample, 3 - dim)), class)
call setAffinity(affinity(@shape(sample)), sample(:,:), dim, tform(1:size(sample, 3 - dim), 1:size(sample, 3 - dim)), class, tlate(1:size(sample, 3 - dim)))
Return an affine-transformation of the input sample of shape (1:ndim) or (1:ndim, 1:nsam) or (1:nsam,...
This module contains classes and procedures for affine transformation of multivariate samples.
Warning
The condition 1 <= dim .and. dim <= rank(sample) must hold for the corresponding input arguments.
The condition size(tlate) == size(sample, 3 - dim) must hold for the corresponding input arguments.
The condition all(shape(tform) == size(sample, 3 - dim)) 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
getShifted
setShifted
setAffinity
setAffinity


Example usage

1program example
2
3 use pm_kind, only: SK, IK, LK
4 use pm_matrixChol, only: getMatChol, uppDia, lowDia
5 use pm_sampleAffinity, only: upperDiag, upperUnit
6 use pm_sampleAffinity, only: lowerDiag, lowerUnit
7 use pm_matrixInit, only: setMatInit, lowDia
9 use pm_sampleAffinity, only: genrecmat
10 use pm_arrayResize, only: setResized
11 use pm_distUnif, only: getUnifRand
12 use pm_arrayFill, only: getFilled
13 use pm_matrixInv, only: getMatInv
14 use pm_distCov, only: getCovRand
15 use pm_io, only: display_type
16
17 implicit none
18
19 type(display_type) :: disp
20 integer(IK) :: itry, ntry = 5
21 integer(IK) :: dim, isam, ndim, nsam
22 disp = display_type(file = "main.out.F90")
23
24 call disp%skip()
25 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
26 call disp%show("!Transform 2D sample along the second dimension.")
27 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
28 call disp%skip()
29
30 block
31 use pm_kind, only: RKG => RKS
32 real(RKG), allocatable :: tlate(:), tform(:,:), sample(:,:), affinInv(:,:), affinity(:,:)
33 do itry = 1, ntry
34 call disp%skip()
35 call disp%show("dim = 2; ndim = getUnifRand(1, 3); nsam = getUnifRand(4, 5)")
36 dim = 2; ndim = getUnifRand(1, 3); nsam = getUnifRand(4, 5)
37 call disp%show("[dim, ndim, nsam]")
38 call disp%show( [dim, ndim, nsam] )
39 call disp%show("tlate = getUnifRand(-5, +5, ndim)")
40 tlate = getUnifRand(-5, +5, ndim)
41 call disp%show("tlate")
42 call disp%show( tlate )
43 call disp%show("sample = getUnifRand(-9, +9, ndim, nsam)")
44 sample = getUnifRand(-9, +9, ndim, nsam)
45 call disp%show("sample")
46 call disp%show( sample )
47 call disp%show("call setResized(affinity, shape(sample, IK))")
48 call setResized(affinity, shape(sample, IK))
49 call disp%show("call setResized(affinInv, shape(sample, IK))")
50 call setResized(affinInv, shape(sample, IK))
51 call disp%show("tform = getCovRand(1., ndim)")
52 tform = getCovRand(1., ndim)
53 call disp%show("tform")
54 call disp%show( tform )
55 call disp%show("call setAffinity(affinity, sample, dim, tform, genrecmat) ! whole sample transformation.")
56 call setAffinity(affinity, sample, dim, tform, genrecmat) ! whole sample transformation.
57 call disp%show("transpose(affinity)")
58 call disp%show( transpose(affinity) )
59 call disp%show("do isam = 1, nsam")
60 call disp%show(" call setAffinity(affinity(:,isam), sample(:,isam), tform, genrecmat) ! point-wise transformation")
61 call disp%show("end do")
62 do isam = 1, nsam
63 call setAffinity(affinity(:,isam), sample(:,isam), tform, genrecmat) ! point-wise transformation
64 end do
65 call disp%show("transpose(affinity)")
66 call disp%show( transpose(affinity) )
67 call disp%show("call setAffinity(affinInv, affinity, dim, getMatInv(tform), genrecmat) ! inverse transformation.")
68 call setAffinity(affinInv, affinity, dim, getMatInv(tform), genrecmat) ! inverse transformation.
69 call disp%show("transpose(affinInv)")
70 call disp%show( transpose(affinInv) )
71 call disp%show("transpose(sample) ! for comparison with affinInv.")
72 call disp%show( transpose(sample) )
73 call disp%skip()
74 call disp%show("call setMatInit(tform(2:,1:ndim-1), lowDia, vlow = 0._RKG, vdia = 0._RKG) ! make tform an upper-diagonal matrix.")
75 call setMatInit(tform(2:,1:ndim-1), lowDia, vlow = 0._RKG, vdia = 0._RKG) ! make tform an upper-diagonal matrix.
76 call disp%show("tform")
77 call disp%show( tform )
78 call disp%show("call setAffinity(affinity, sample, dim, tform, class = upperDiag) ! whole sample transformation.")
79 call setAffinity(affinity, sample, dim, tform, class = upperDiag) ! whole sample transformation.
80 call disp%show("transpose(affinity)")
81 call disp%show( transpose(affinity) )
82 call disp%show("do isam = 1, nsam")
83 call disp%show(" call setAffinity(affinity(:,isam), sample(:,isam), tform, upperDiag) ! point-wise transformation")
84 call disp%show("end do")
85 do isam = 1, nsam
86 call setAffinity(affinity(:,isam), sample(:,isam), tform, upperDiag) ! point-wise transformation
87 end do
88 call disp%show("transpose(affinity)")
89 call disp%show( transpose(affinity) )
90 call disp%show("call setAffinity(affinInv, affinity, dim, getMatInv(tform, upperDiag), class = upperDiag) ! inverse transformation.")
91 call setAffinity(affinInv, affinity, dim, getMatInv(tform, upperDiag), class = upperDiag) ! inverse transformation.
92 call disp%show("transpose(affinInv)")
93 call disp%show( transpose(affinInv) )
94 call disp%show("transpose(sample) ! for comparison with affinInv.")
95 call disp%show( transpose(sample) )
96 call disp%skip()
97 call disp%show("call setAffinity(affinity, sample, dim, tform, class = upperUnit) ! whole sample transformation.")
98 call setAffinity(affinity, sample, dim, tform, class = upperUnit) ! whole sample transformation.
99 call disp%show("transpose(affinity)")
100 call disp%show( transpose(affinity) )
101 call disp%show("do isam = 1, nsam")
102 call disp%show(" call setAffinity(affinity(:,isam), sample(:,isam), tform, upperUnit) ! point-wise transformation")
103 call disp%show("end do")
104 do isam = 1, nsam
105 call setAffinity(affinity(:,isam), sample(:,isam), tform, upperUnit) ! point-wise transformation
106 end do
107 call disp%show("transpose(affinity)")
108 call disp%show( transpose(affinity) )
109 call disp%show("call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = upperUnit) ! inverse transformation.")
110 call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = upperUnit) ! inverse transformation.
111 call disp%show("transpose(affinInv)")
112 call disp%show( transpose(affinInv) )
113 call disp%show("transpose(sample) ! for comparison with affinInv.")
114 call disp%show( transpose(sample) )
115 call disp%skip()
116 call disp%show("call setAffinity(affinity, sample, dim, tform, class = lowerDiag) ! whole sample transformation.")
117 call setAffinity(affinity, sample, dim, tform, class = lowerDiag) ! whole sample transformation.
118 call disp%show("transpose(affinity)")
119 call disp%show( transpose(affinity) )
120 call disp%show("do isam = 1, nsam")
121 call disp%show(" call setAffinity(affinity(:,isam), sample(:,isam), tform, lowerDiag) ! point-wise transformation")
122 call disp%show("end do")
123 do isam = 1, nsam
124 call setAffinity(affinity(:,isam), sample(:,isam), tform, lowerDiag) ! point-wise transformation
125 end do
126 call disp%show("transpose(affinity)")
127 call disp%show( transpose(affinity) )
128 call disp%show("call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerDiag) ! inverse transformation.")
129 call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerDiag) ! inverse transformation.
130 call disp%show("transpose(affinInv)")
131 call disp%show( transpose(affinInv) )
132 call disp%show("transpose(sample) ! for comparison with affinInv.")
133 call disp%show( transpose(sample) )
134 call disp%skip()
135 call disp%show("call setAffinity(affinity, sample, dim, tform, class = lowerUnit) ! whole sample transformation.")
136 call setAffinity(affinity, sample, dim, tform, class = lowerUnit) ! whole sample transformation.
137 call disp%show("transpose(affinity)")
138 call disp%show( transpose(affinity) )
139 call disp%show("do isam = 1, nsam")
140 call disp%show(" call setAffinity(affinity(:,isam), sample(:,isam), tform, lowerUnit) ! point-wise transformation")
141 call disp%show("end do")
142 do isam = 1, nsam
143 call setAffinity(affinity(:,isam), sample(:,isam), tform, lowerUnit) ! point-wise transformation
144 end do
145 call disp%show("transpose(affinity)")
146 call disp%show( transpose(affinity) )
147 call disp%show("call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerUnit) ! inverse transformation.")
148 call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerUnit) ! inverse transformation.
149 call disp%show("transpose(affinInv)")
150 call disp%show( transpose(affinInv) )
151 call disp%show("transpose(sample) ! for comparison with affinInv.")
152 call disp%show( transpose(sample) )
153 end do
154 end block
155
156 call disp%skip()
157 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
158 call disp%show("!Transform 2D sample along the first dimension.")
159 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
160 call disp%skip()
161
162 block
163 use pm_kind, only: RKG => RKS
164 real(RKG), allocatable :: tlate(:), tform(:,:), sample(:,:), affinInv(:,:), affinity(:,:)
165 do itry = 1, ntry
166 call disp%skip()
167 call disp%show("dim = 1; ndim = getUnifRand(1, 3); nsam = getUnifRand(4, 5)")
168 dim = 1; ndim = getUnifRand(1, 3); nsam = getUnifRand(4, 5)
169 call disp%show("[dim, ndim, nsam]")
170 call disp%show( [dim, ndim, nsam] )
171 call disp%show("tlate = getUnifRand(-5, +5, ndim)")
172 tlate = getUnifRand(-5, +5, ndim)
173 call disp%show("tlate")
174 call disp%show( tlate )
175 call disp%show("sample = getUnifRand(-9, +9, nsam, ndim)")
176 sample = getUnifRand(-9, +9, nsam, ndim)
177 call disp%show("sample")
178 call disp%show( sample )
179 call disp%show("call setResized(affinity, shape(sample, IK))")
180 call setResized(affinity, shape(sample, IK))
181 call disp%show("call setResized(affinInv, shape(sample, IK))")
182 call setResized(affinInv, shape(sample, IK))
183 call disp%show("tform = getCovRand(1., ndim)")
184 tform = getCovRand(1., ndim)
185 call disp%show("tform")
186 call disp%show( tform )
187 call disp%show("call setAffinity(affinity, sample, dim, tform) ! whole sample transformation.")
188 call setAffinity(affinity, sample, dim, tform, genrecmat) ! whole sample transformation.
189 call disp%show("affinity")
190 call disp%show( affinity )
191 call disp%show("do isam = 1, nsam")
192 call disp%show(" call setAffinity(affinity(isam,:), sample(isam,:), tform, genrecmat) ! point-wise transformation")
193 call disp%show("end do")
194 do isam = 1, nsam
195 call setAffinity(affinity(isam,:), sample(isam,:), tform, genrecmat) ! point-wise transformation
196 end do
197 call disp%show("affinity")
198 call disp%show( affinity )
199 call disp%show("call setAffinity(affinInv, affinity, dim, getMatInv(tform), genrecmat) ! inverse transformation.")
200 call setAffinity(affinInv, affinity, dim, getMatInv(tform), genrecmat) ! inverse transformation.
201 call disp%show("affinInv")
202 call disp%show( affinInv )
203 call disp%show("sample ! for comparison with affinInv.")
204 call disp%show( sample )
205 call disp%skip()
206 call disp%show("call setMatInit(tform(2:,1:ndim-1), lowDia, vlow = 0._RKG, vdia = 0._RKG) ! make tform an upper-diagonal matrix.")
207 call setMatInit(tform(2:,1:ndim-1), lowDia, vlow = 0._RKG, vdia = 0._RKG) ! make tform an upper-diagonal matrix.
208 call disp%show("tform")
209 call disp%show( tform )
210 call disp%show("call setAffinity(affinity, sample, dim, tform, class = upperDiag) ! whole sample transformation.")
211 call setAffinity(affinity, sample, dim, tform, class = upperDiag) ! whole sample transformation.
212 call disp%show("affinity")
213 call disp%show( affinity )
214 call disp%show("do isam = 1, nsam")
215 call disp%show(" call setAffinity(affinity(isam,:), sample(isam,:), tform, upperDiag) ! point-wise transformation")
216 call disp%show("end do")
217 do isam = 1, nsam
218 call setAffinity(affinity(isam,:), sample(isam,:), tform, upperDiag) ! point-wise transformation
219 end do
220 call disp%show("affinity")
221 call disp%show( affinity )
222 call disp%show("call setAffinity(affinInv, affinity, dim, getMatInv(tform, upperDiag), class = upperDiag) ! inverse transformation.")
223 call setAffinity(affinInv, affinity, dim, getMatInv(tform, upperDiag), class = upperDiag) ! inverse transformation.
224 call disp%show("affinInv")
225 call disp%show( affinInv )
226 call disp%show("sample ! for comparison with affinInv.")
227 call disp%show( sample )
228 call disp%skip()
229 call disp%show("call setAffinity(affinity, sample, dim, tform, class = upperUnit) ! whole sample transformation.")
230 call setAffinity(affinity, sample, dim, tform, class = upperUnit) ! whole sample transformation.
231 call disp%show("affinity")
232 call disp%show( affinity )
233 call disp%show("do isam = 1, nsam")
234 call disp%show(" call setAffinity(affinity(isam,:), sample(isam,:), tform, upperUnit) ! point-wise transformation")
235 call disp%show("end do")
236 do isam = 1, nsam
237 call setAffinity(affinity(isam,:), sample(isam,:), tform, upperUnit) ! point-wise transformation
238 end do
239 call disp%show("affinity")
240 call disp%show( affinity )
241 call disp%show("call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = upperUnit) ! inverse transformation.")
242 call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = upperUnit) ! inverse transformation.
243 call disp%show("affinInv")
244 call disp%show( affinInv )
245 call disp%show("sample ! for comparison with affinInv.")
246 call disp%show( sample )
247 call disp%skip()
248 call disp%show("call setAffinity(affinity, sample, dim, tform, class = lowerDiag) ! whole sample transformation.")
249 call setAffinity(affinity, sample, dim, tform, class = lowerDiag) ! whole sample transformation.
250 call disp%show("affinity")
251 call disp%show( affinity )
252 call disp%show("do isam = 1, nsam")
253 call disp%show(" call setAffinity(affinity(isam,:), sample(isam,:), tform, lowerDiag) ! point-wise transformation")
254 call disp%show("end do")
255 do isam = 1, nsam
256 call setAffinity(affinity(isam,:), sample(isam,:), tform, lowerDiag) ! point-wise transformation
257 end do
258 call disp%show("affinity")
259 call disp%show( affinity )
260 call disp%show("call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerDiag) ! inverse transformation.")
261 call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerDiag) ! inverse transformation.
262 call disp%show("affinInv")
263 call disp%show( affinInv )
264 call disp%show("sample ! for comparison with affinInv.")
265 call disp%show( sample )
266 call disp%skip()
267 call disp%show("call setAffinity(affinity, sample, dim, tform, class = lowerUnit) ! whole sample transformation.")
268 call setAffinity(affinity, sample, dim, tform, class = lowerUnit) ! whole sample transformation.
269 call disp%show("affinity")
270 call disp%show( affinity )
271 call disp%show("do isam = 1, nsam")
272 call disp%show(" call setAffinity(affinity(isam,:), sample(isam,:), tform, lowerUnit) ! point-wise transformation")
273 call disp%show("end do")
274 do isam = 1, nsam
275 call setAffinity(affinity(isam,:), sample(isam,:), tform, lowerUnit) ! point-wise transformation
276 end do
277 call disp%show("affinity")
278 call disp%show( affinity )
279 call disp%show("call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerUnit) ! inverse transformation.")
280 call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerUnit) ! inverse transformation.
281 call disp%show("affinInv")
282 call disp%show( affinInv )
283 call disp%show("sample ! for comparison with affinInv.")
284 call disp%show( sample )
285 end do
286 end block
287
288 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
289 ! Output an example rand array for visualization.
290 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
291
292 block
293 use pm_mathConst, only: PI
294 use pm_container, only: css_type
295 use pm_matrixInit, only: uppLowDia
296 use pm_matrixInit, only: getMatInit
297 use pm_arraySpace, only: getLinSpace
298 use pm_sampleAffinity, only: genrecmat
299 use pm_io, only: getErrTableWrite, trans
300 use pm_distUnifEll, only: getUnifEllRand, uppDia
301 use pm_sampleShift, only: getShifted
302 use pm_matrixInv, only: getMatInv
303 use pm_kind, only: RKG => RKS
304 type(css_type), allocatable :: shapes(:), tforms(:)
305 integer(IK) :: ishape, itform
306 integer(IK), parameter :: nsam = 1000, ndim = 2, dim = 2
307 real(RKG), allocatable :: angle, tlate(:), cov(:,:), tform(:,:), sample(:,:), affinity(:,:), affinInv(:,:)
308 tforms = [css_type(SK_"warp"), css_type(SK_"rotation")]
309 shapes = [css_type(SK_"circle"), css_type(SK_"square")]
310 tlate = getUnifRand(0, 5, ndim)
311 do ishape = 1, size(shapes)
312 if (shapes(ishape)%val == "circle") then
313 sample = getUnifEllRand(getMatInit([ndim, ndim], uppLowDia, 0., 0., getLinSpace(1., real(ndim), ndim)), uppDia, nsam)
314 elseif (shapes(ishape)%val == "square") then
315 sample = getUnifRand(0., 1., ndim, nsam)
316 else
317 error stop "Unrecognized shape."
318 end if
319 call setResized(affinity, shape(sample, IK))
320 call setResized(affinInv, shape(sample, IK))
321 do itform = 1, size(tforms)
322 if (tforms(itform)%val == "rotation") then
323 angle = PI / 4
324 tform = reshape([cos(angle), -sin(angle), sin(angle), cos(angle)], [ndim, ndim])
325 call setAffinity(affinity, sample, dim, tform, genrecmat, tlate)
326 elseif (tforms(itform)%val == "warp") then
327 cov = reshape([1., .5, .5, 1.], shape = [ndim, ndim]) ! getCovRand(1., ndim)
328 tform = getMatChol(cov, lowDia)
329 call setAffinity(affinity, sample, dim, tform, lowerDiag, tlate)
330 else
331 error stop "Unrecognized tform."
332 end if
333 call setAffinity(affinInv, getShifted(affinity, dim, -tlate), dim, getMatInv(tform), genrecmat)
334 if (0 /= getErrTableWrite(SK_"setAffinity."//tforms(itform)%val//SK_"."//shapes(ishape)%val//SK_".sample.txt", sample, trans)) error stop "Failed table-write."
335 if (0 /= getErrTableWrite(SK_"setAffinity."//tforms(itform)%val//SK_"."//shapes(ishape)%val//SK_".affinity.txt", affinity, trans)) error stop "Failed table-write."
336 if (0 /= getErrTableWrite(SK_"setAffinity."//tforms(itform)%val//SK_"."//shapes(ishape)%val//SK_".affinInv.txt", affinInv, trans)) error stop "Failed table-write."
337 end do
338 end do
339 end block
340
341end program example
Generate and return an array of the specified rank and shape of arbitrary intrinsic type and kind wit...
Allocate or resize (shrink or expand) an input allocatable scalar string or array of rank 1....
Generate count evenly spaced points over the interval [x1, x2] if x1 < x2, or [x2,...
Generate and return a random positive-definite (correlation or covariance) matrix using the Gram meth...
Definition: pm_distCov.F90:394
Generate and return a (collection) of random vector(s) of size ndim from the ndim-dimensional MultiVa...
Generate and return a scalar or a contiguous array of rank 1 of length s1 of randomly uniformly distr...
Generate and return the iostat code resulting from writing the input table of rank 1 or 2 to the spec...
Definition: pm_io.F90:5940
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 upper or the lower Cholesky factorization of the input symmetric positive-def...
Generate and return a matrix of shape (shape(1), shape(2)) with the upper/lower triangle and the diag...
Set the upper/lower triangle and the diagonal elements of the input matrix of arbitrary shape (:,...
Generate and return the full inverse of an input matrix of general or triangular form directly or thr...
Generate a sample of shape (nsam), or (ndim, nsam) or (nsam, ndim) that is shifted by the specified i...
This module contains procedures and generic interfaces for convenient allocation and filling of array...
This module contains procedures and generic interfaces for resizing allocatable arrays of various typ...
This module contains procedures and generic interfaces for generating arrays with linear or logarithm...
This module contains the derived types for generating allocatable containers of scalar,...
This module contains classes and procedures for generating random matrices distributed on the space o...
Definition: pm_distCov.F90:72
This module contains classes and procedures for computing various statistical quantities related to t...
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 LK
The default logical kind in the ParaMonte library: kind(.true.) in Fortran, kind(....
Definition: pm_kind.F90:541
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 relevant mathematical constants.
real(RKB), parameter PI
The scalar real constant of kind with highest available precision RKB representing the irrational num...
This module contains procedures and generic interfaces for computing the Cholesky factorization of po...
This module contains procedures and generic interfaces relevant to generating and initializing matric...
This module contains abstract and concrete derived types and procedures related to the inversion of s...
This module contains classes and procedures for shifting univariate or multivariate samples by arbitr...
This is the css_type type for generating instances of container of scalar of string objects.
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!Transform 2D sample along the second dimension.
4!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
5
6
7dim = 2; ndim = getUnifRand(1, 3); nsam = getUnifRand(4, 5)
8[dim, ndim, nsam]
9+2, +3, +5
10tlate = getUnifRand(-5, +5, ndim)
11tlate
12+3.00000000, +2.00000000, +2.00000000
13sample = getUnifRand(-9, +9, ndim, nsam)
14sample
15+1.00000000, +1.00000000, -3.00000000, +1.00000000, -9.00000000
16+3.00000000, +9.00000000, -3.00000000, +0.00000000, +1.00000000
17-7.00000000, -4.00000000, -3.00000000, +0.00000000, +9.00000000
18call setResized(affinity, shape(sample, IK))
19call setResized(affinInv, shape(sample, IK))
20tform = getCovRand(1., ndim)
21tform
22+1.00000000, -0.941470921, -0.876164198
23-0.941470921, +1.00000000, +0.983183086
24-0.876164198, +0.983183086, +1.00000000
25call setAffinity(affinity, sample, dim, tform, genrecmat) ! whole sample transformation.
26transpose(affinity)
27+4.30873632, -4.82375240, -4.92661476
28-3.96858120, +4.12579632, +3.97248363
29+2.45290542, -3.12513638, -3.32105660
30+1.00000000, -0.941470921, -0.876164198
31-17.8269501, +18.3218861, +17.8686600
32do isam = 1, nsam
33 call setAffinity(affinity(:,isam), sample(:,isam), tform, genrecmat) ! point-wise transformation
34end do
35transpose(affinity)
36+4.30873632, -4.82375240, -4.92661476
37-3.96858120, +4.12579632, +3.97248363
38+2.45290542, -3.12513638, -3.32105660
39+1.00000000, -0.941470921, -0.876164198
40-17.8269501, +18.3218861, +17.8686600
41call setAffinity(affinInv, affinity, dim, getMatInv(tform), genrecmat) ! inverse transformation.
42transpose(affinInv)
43+1.00000000, +3.00000000, -6.99996948
44+0.999984741, +9.00000000, -4.00003052
45-2.99997711, -2.99996948, -2.99996948
46+1.00000000, -0.762939453E-5, +0.762939453E-5
47-9.00006104, +1.00000000, +8.99987793
48transpose(sample) ! for comparison with affinInv.
49+1.00000000, +3.00000000, -7.00000000
50+1.00000000, +9.00000000, -4.00000000
51-3.00000000, -3.00000000, -3.00000000
52+1.00000000, +0.00000000, +0.00000000
53-9.00000000, +1.00000000, +9.00000000
54
55call setMatInit(tform(2:,1:ndim-1), lowDia, vlow = 0._RKG, vdia = 0._RKG) ! make tform an upper-diagonal matrix.
56tform
57+1.00000000, -0.941470921, -0.876164198
58+0.00000000, +1.00000000, +0.983183086
59+0.00000000, +0.00000000, +1.00000000
60call setAffinity(affinity, sample, dim, tform, class = upperDiag) ! whole sample transformation.
61transpose(affinity)
62+4.30873632, -3.88228178, -7.00000000
63-3.96858120, +5.06726742, -4.00000000
64+2.45290565, -5.94954920, -3.00000000
65+1.00000000, +0.00000000, +0.00000000
66-17.8269501, +9.84864807, +9.00000000
67do isam = 1, nsam
68 call setAffinity(affinity(:,isam), sample(:,isam), tform, upperDiag) ! point-wise transformation
69end do
70transpose(affinity)
71+4.30873632, -3.88228178, -7.00000000
72-3.96858120, +5.06726742, -4.00000000
73+2.45290565, -5.94954920, -3.00000000
74+1.00000000, +0.00000000, +0.00000000
75-17.8269501, +9.84864807, +9.00000000
76call setAffinity(affinInv, affinity, dim, getMatInv(tform, upperDiag), class = upperDiag) ! inverse transformation.
77transpose(affinInv)
78+0.999999523, +3.00000000, -7.00000000
79+1.00000000, +9.00000000, -4.00000000
80-2.99999952, -3.00000000, -3.00000000
81+1.00000000, +0.00000000, +0.00000000
82-9.00000095, +1.00000000, +9.00000000
83transpose(sample) ! for comparison with affinInv.
84+1.00000000, +3.00000000, -7.00000000
85+1.00000000, +9.00000000, -4.00000000
86-3.00000000, -3.00000000, -3.00000000
87+1.00000000, +0.00000000, +0.00000000
88-9.00000000, +1.00000000, +9.00000000
89
90call setAffinity(affinity, sample, dim, tform, class = upperUnit) ! whole sample transformation.
91transpose(affinity)
92+4.30873632, -3.88228178, -7.00000000
93-3.96858120, +5.06726742, -4.00000000
94+2.45290565, -5.94954920, -3.00000000
95+1.00000000, +0.00000000, +0.00000000
96-17.8269501, +9.84864807, +9.00000000
97do isam = 1, nsam
98 call setAffinity(affinity(:,isam), sample(:,isam), tform, upperUnit) ! point-wise transformation
99end do
100transpose(affinity)
101+4.30873632, -3.88228178, -7.00000000
102-3.96858120, +5.06726742, -4.00000000
103+2.45290565, -5.94954920, -3.00000000
104+1.00000000, +0.00000000, +0.00000000
105-17.8269501, +9.84864807, +9.00000000
106call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = upperUnit) ! inverse transformation.
107transpose(affinInv)
108+0.999999523, +3.00000000, -7.00000000
109+1.00000000, +9.00000000, -4.00000000
110-2.99999952, -3.00000000, -3.00000000
111+1.00000000, +0.00000000, +0.00000000
112-9.00000095, +1.00000000, +9.00000000
113transpose(sample) ! for comparison with affinInv.
114+1.00000000, +3.00000000, -7.00000000
115+1.00000000, +9.00000000, -4.00000000
116-3.00000000, -3.00000000, -3.00000000
117+1.00000000, +0.00000000, +0.00000000
118-9.00000000, +1.00000000, +9.00000000
119
120call setAffinity(affinity, sample, dim, tform, class = lowerDiag) ! whole sample transformation.
121transpose(affinity)
122+1.00000000, +3.00000000, -7.00000000
123+1.00000000, +9.00000000, -4.00000000
124-3.00000000, -3.00000000, -3.00000000
125+1.00000000, +0.00000000, +0.00000000
126-9.00000000, +1.00000000, +9.00000000
127do isam = 1, nsam
128 call setAffinity(affinity(:,isam), sample(:,isam), tform, lowerDiag) ! point-wise transformation
129end do
130transpose(affinity)
131+1.00000000, +3.00000000, -7.00000000
132+1.00000000, +9.00000000, -4.00000000
133-3.00000000, -3.00000000, -3.00000000
134+1.00000000, +0.00000000, +0.00000000
135-9.00000000, +1.00000000, +9.00000000
136call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerDiag) ! inverse transformation.
137transpose(affinInv)
138+1.00000000, +3.00000000, -7.00000000
139+1.00000000, +9.00000000, -4.00000000
140-3.00000000, -3.00000000, -3.00000000
141+1.00000000, +0.00000000, +0.00000000
142-9.00000000, +1.00000000, +9.00000000
143transpose(sample) ! for comparison with affinInv.
144+1.00000000, +3.00000000, -7.00000000
145+1.00000000, +9.00000000, -4.00000000
146-3.00000000, -3.00000000, -3.00000000
147+1.00000000, +0.00000000, +0.00000000
148-9.00000000, +1.00000000, +9.00000000
149
150call setAffinity(affinity, sample, dim, tform, class = lowerUnit) ! whole sample transformation.
151transpose(affinity)
152+1.00000000, +3.00000000, -7.00000000
153+1.00000000, +9.00000000, -4.00000000
154-3.00000000, -3.00000000, -3.00000000
155+1.00000000, +0.00000000, +0.00000000
156-9.00000000, +1.00000000, +9.00000000
157do isam = 1, nsam
158 call setAffinity(affinity(:,isam), sample(:,isam), tform, lowerUnit) ! point-wise transformation
159end do
160transpose(affinity)
161+1.00000000, +3.00000000, -7.00000000
162+1.00000000, +9.00000000, -4.00000000
163-3.00000000, -3.00000000, -3.00000000
164+1.00000000, +0.00000000, +0.00000000
165-9.00000000, +1.00000000, +9.00000000
166call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerUnit) ! inverse transformation.
167transpose(affinInv)
168+1.00000000, +3.00000000, -7.00000000
169+1.00000000, +9.00000000, -4.00000000
170-3.00000000, -3.00000000, -3.00000000
171+1.00000000, +0.00000000, +0.00000000
172-9.00000000, +1.00000000, +9.00000000
173transpose(sample) ! for comparison with affinInv.
174+1.00000000, +3.00000000, -7.00000000
175+1.00000000, +9.00000000, -4.00000000
176-3.00000000, -3.00000000, -3.00000000
177+1.00000000, +0.00000000, +0.00000000
178-9.00000000, +1.00000000, +9.00000000
179
180dim = 2; ndim = getUnifRand(1, 3); nsam = getUnifRand(4, 5)
181[dim, ndim, nsam]
182+2, +2, +5
183tlate = getUnifRand(-5, +5, ndim)
184tlate
185+2.00000000, -3.00000000
186sample = getUnifRand(-9, +9, ndim, nsam)
187sample
188+7.00000000, -4.00000000, +0.00000000, +1.00000000, -8.00000000
189-9.00000000, +9.00000000, -2.00000000, +4.00000000, -6.00000000
190call setResized(affinity, shape(sample, IK))
191call setResized(affinInv, shape(sample, IK))
192tform = getCovRand(1., ndim)
193tform
194+1.00000000, +0.319688737
195+0.319688737, +1.00000000
196call setAffinity(affinity, sample, dim, tform, genrecmat) ! whole sample transformation.
197transpose(affinity)
198+4.12280130, -6.76217890
199-1.12280130, +7.72124481
200-0.639377475, -2.00000000
201+2.27875495, +4.31968880
202-9.91813278, -8.55751038
203do isam = 1, nsam
204 call setAffinity(affinity(:,isam), sample(:,isam), tform, genrecmat) ! point-wise transformation
205end do
206transpose(affinity)
207+4.12280130, -6.76217890
208-1.12280130, +7.72124481
209-0.639377475, -2.00000000
210+2.27875495, +4.31968880
211-9.91813278, -8.55751038
212call setAffinity(affinInv, affinity, dim, getMatInv(tform), genrecmat) ! inverse transformation.
213transpose(affinInv)
214+6.99999952, -9.00000000
215-3.99999952, +9.00000000
216+0.00000000, -1.99999988
217+0.999999881, +3.99999976
218-8.00000000, -6.00000000
219transpose(sample) ! for comparison with affinInv.
220+7.00000000, -9.00000000
221-4.00000000, +9.00000000
222+0.00000000, -2.00000000
223+1.00000000, +4.00000000
224-8.00000000, -6.00000000
225
226call setMatInit(tform(2:,1:ndim-1), lowDia, vlow = 0._RKG, vdia = 0._RKG) ! make tform an upper-diagonal matrix.
227tform
228+1.00000000, +0.319688737
229+0.00000000, +1.00000000
230call setAffinity(affinity, sample, dim, tform, class = upperDiag) ! whole sample transformation.
231transpose(affinity)
232+4.12280130, -9.00000000
233-1.12280130, +9.00000000
234-0.639377475, -2.00000000
235+2.27875495, +4.00000000
236-9.91813278, -6.00000000
237do isam = 1, nsam
238 call setAffinity(affinity(:,isam), sample(:,isam), tform, upperDiag) ! point-wise transformation
239end do
240transpose(affinity)
241+4.12280130, -9.00000000
242-1.12280130, +9.00000000
243-0.639377475, -2.00000000
244+2.27875495, +4.00000000
245-9.91813278, -6.00000000
246call setAffinity(affinInv, affinity, dim, getMatInv(tform, upperDiag), class = upperDiag) ! inverse transformation.
247transpose(affinInv)
248+7.00000000, -9.00000000
249-4.00000000, +9.00000000
250+0.00000000, -2.00000000
251+1.00000000, +4.00000000
252-8.00000000, -6.00000000
253transpose(sample) ! for comparison with affinInv.
254+7.00000000, -9.00000000
255-4.00000000, +9.00000000
256+0.00000000, -2.00000000
257+1.00000000, +4.00000000
258-8.00000000, -6.00000000
259
260call setAffinity(affinity, sample, dim, tform, class = upperUnit) ! whole sample transformation.
261transpose(affinity)
262+4.12280130, -9.00000000
263-1.12280130, +9.00000000
264-0.639377475, -2.00000000
265+2.27875495, +4.00000000
266-9.91813278, -6.00000000
267do isam = 1, nsam
268 call setAffinity(affinity(:,isam), sample(:,isam), tform, upperUnit) ! point-wise transformation
269end do
270transpose(affinity)
271+4.12280130, -9.00000000
272-1.12280130, +9.00000000
273-0.639377475, -2.00000000
274+2.27875495, +4.00000000
275-9.91813278, -6.00000000
276call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = upperUnit) ! inverse transformation.
277transpose(affinInv)
278+7.00000000, -9.00000000
279-4.00000000, +9.00000000
280+0.00000000, -2.00000000
281+1.00000000, +4.00000000
282-8.00000000, -6.00000000
283transpose(sample) ! for comparison with affinInv.
284+7.00000000, -9.00000000
285-4.00000000, +9.00000000
286+0.00000000, -2.00000000
287+1.00000000, +4.00000000
288-8.00000000, -6.00000000
289
290call setAffinity(affinity, sample, dim, tform, class = lowerDiag) ! whole sample transformation.
291transpose(affinity)
292+7.00000000, -9.00000000
293-4.00000000, +9.00000000
294+0.00000000, -2.00000000
295+1.00000000, +4.00000000
296-8.00000000, -6.00000000
297do isam = 1, nsam
298 call setAffinity(affinity(:,isam), sample(:,isam), tform, lowerDiag) ! point-wise transformation
299end do
300transpose(affinity)
301+7.00000000, -9.00000000
302-4.00000000, +9.00000000
303+0.00000000, -2.00000000
304+1.00000000, +4.00000000
305-8.00000000, -6.00000000
306call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerDiag) ! inverse transformation.
307transpose(affinInv)
308+7.00000000, -9.00000000
309-4.00000000, +9.00000000
310+0.00000000, -2.00000000
311+1.00000000, +4.00000000
312-8.00000000, -6.00000000
313transpose(sample) ! for comparison with affinInv.
314+7.00000000, -9.00000000
315-4.00000000, +9.00000000
316+0.00000000, -2.00000000
317+1.00000000, +4.00000000
318-8.00000000, -6.00000000
319
320call setAffinity(affinity, sample, dim, tform, class = lowerUnit) ! whole sample transformation.
321transpose(affinity)
322+7.00000000, -9.00000000
323-4.00000000, +9.00000000
324+0.00000000, -2.00000000
325+1.00000000, +4.00000000
326-8.00000000, -6.00000000
327do isam = 1, nsam
328 call setAffinity(affinity(:,isam), sample(:,isam), tform, lowerUnit) ! point-wise transformation
329end do
330transpose(affinity)
331+7.00000000, -9.00000000
332-4.00000000, +9.00000000
333+0.00000000, -2.00000000
334+1.00000000, +4.00000000
335-8.00000000, -6.00000000
336call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerUnit) ! inverse transformation.
337transpose(affinInv)
338+7.00000000, -9.00000000
339-4.00000000, +9.00000000
340+0.00000000, -2.00000000
341+1.00000000, +4.00000000
342-8.00000000, -6.00000000
343transpose(sample) ! for comparison with affinInv.
344+7.00000000, -9.00000000
345-4.00000000, +9.00000000
346+0.00000000, -2.00000000
347+1.00000000, +4.00000000
348-8.00000000, -6.00000000
349
350dim = 2; ndim = getUnifRand(1, 3); nsam = getUnifRand(4, 5)
351[dim, ndim, nsam]
352+2, +3, +4
353tlate = getUnifRand(-5, +5, ndim)
354tlate
355-5.00000000, -2.00000000, -1.00000000
356sample = getUnifRand(-9, +9, ndim, nsam)
357sample
358+7.00000000, -5.00000000, +9.00000000, +7.00000000
359-2.00000000, -3.00000000, +0.00000000, +1.00000000
360+9.00000000, +1.00000000, +2.00000000, +1.00000000
361call setResized(affinity, shape(sample, IK))
362call setResized(affinInv, shape(sample, IK))
363tform = getCovRand(1., ndim)
364tform
365+1.00000000, +0.737021267, -0.490478843
366+0.737021267, +1.00000000, +0.224671006
367-0.490478843, +0.224671006, +1.00000000
368call setAffinity(affinity, sample, dim, tform, genrecmat) ! whole sample transformation.
369transpose(affinity)
370+1.11164808, +5.18118763, +5.11730576
371-7.70154285, -6.46043539, +2.77838135
372+8.01904202, +7.08253384, -2.41430950
373+7.24654245, +6.38381958, -2.20868111
374do isam = 1, nsam
375 call setAffinity(affinity(:,isam), sample(:,isam), tform, genrecmat) ! point-wise transformation
376end do
377transpose(affinity)
378+1.11164808, +5.18118763, +5.11730576
379-7.70154285, -6.46043539, +2.77838135
380+8.01904202, +7.08253384, -2.41430950
381+7.24654245, +6.38381958, -2.20868111
382call setAffinity(affinInv, affinity, dim, getMatInv(tform), genrecmat) ! inverse transformation.
383transpose(affinInv)
384+7.00006104, -2.00000000, +9.00000000
385-4.99987793, -3.00024414, +0.999969482
386+8.99954224, +0.335693359E-3, +1.99993896
387+6.99981689, +1.00012207, +1.00006104
388transpose(sample) ! for comparison with affinInv.
389+7.00000000, -2.00000000, +9.00000000
390-5.00000000, -3.00000000, +1.00000000
391+9.00000000, +0.00000000, +2.00000000
392+7.00000000, +1.00000000, +1.00000000
393
394call setMatInit(tform(2:,1:ndim-1), lowDia, vlow = 0._RKG, vdia = 0._RKG) ! make tform an upper-diagonal matrix.
395tform
396+1.00000000, +0.737021267, -0.490478843
397+0.00000000, +1.00000000, +0.224671006
398+0.00000000, +0.00000000, +1.00000000
399call setAffinity(affinity, sample, dim, tform, class = upperDiag) ! whole sample transformation.
400transpose(affinity)
401+1.11164808, +0.220389366E-1, +9.00000000
402-7.70154285, -2.77532911, +1.00000000
403+8.01904202, +0.449342012, +2.00000000
404+7.24654245, +1.22467101, +1.00000000
405do isam = 1, nsam
406 call setAffinity(affinity(:,isam), sample(:,isam), tform, upperDiag) ! point-wise transformation
407end do
408transpose(affinity)
409+1.11164808, +0.220389366E-1, +9.00000000
410-7.70154285, -2.77532911, +1.00000000
411+8.01904202, +0.449342012, +2.00000000
412+7.24654245, +1.22467101, +1.00000000
413call setAffinity(affinInv, affinity, dim, getMatInv(tform, upperDiag), class = upperDiag) ! inverse transformation.
414transpose(affinInv)
415+7.00000048, -2.00000000, +9.00000000
416-5.00000000, -3.00000000, +1.00000000
417+9.00000000, +0.00000000, +2.00000000
418+7.00000000, +1.00000000, +1.00000000
419transpose(sample) ! for comparison with affinInv.
420+7.00000000, -2.00000000, +9.00000000
421-5.00000000, -3.00000000, +1.00000000
422+9.00000000, +0.00000000, +2.00000000
423+7.00000000, +1.00000000, +1.00000000
424
425call setAffinity(affinity, sample, dim, tform, class = upperUnit) ! whole sample transformation.
426transpose(affinity)
427+1.11164808, +0.220389366E-1, +9.00000000
428-7.70154285, -2.77532911, +1.00000000
429+8.01904202, +0.449342012, +2.00000000
430+7.24654245, +1.22467101, +1.00000000
431do isam = 1, nsam
432 call setAffinity(affinity(:,isam), sample(:,isam), tform, upperUnit) ! point-wise transformation
433end do
434transpose(affinity)
435+1.11164808, +0.220389366E-1, +9.00000000
436-7.70154285, -2.77532911, +1.00000000
437+8.01904202, +0.449342012, +2.00000000
438+7.24654245, +1.22467101, +1.00000000
439call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = upperUnit) ! inverse transformation.
440transpose(affinInv)
441+7.00000048, -2.00000000, +9.00000000
442-5.00000000, -3.00000000, +1.00000000
443+9.00000000, +0.00000000, +2.00000000
444+7.00000000, +1.00000000, +1.00000000
445transpose(sample) ! for comparison with affinInv.
446+7.00000000, -2.00000000, +9.00000000
447-5.00000000, -3.00000000, +1.00000000
448+9.00000000, +0.00000000, +2.00000000
449+7.00000000, +1.00000000, +1.00000000
450
451call setAffinity(affinity, sample, dim, tform, class = lowerDiag) ! whole sample transformation.
452transpose(affinity)
453+7.00000000, -2.00000000, +9.00000000
454-5.00000000, -3.00000000, +1.00000000
455+9.00000000, +0.00000000, +2.00000000
456+7.00000000, +1.00000000, +1.00000000
457do isam = 1, nsam
458 call setAffinity(affinity(:,isam), sample(:,isam), tform, lowerDiag) ! point-wise transformation
459end do
460transpose(affinity)
461+7.00000000, -2.00000000, +9.00000000
462-5.00000000, -3.00000000, +1.00000000
463+9.00000000, +0.00000000, +2.00000000
464+7.00000000, +1.00000000, +1.00000000
465call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerDiag) ! inverse transformation.
466transpose(affinInv)
467+7.00000000, -2.00000000, +9.00000000
468-5.00000000, -3.00000000, +1.00000000
469+9.00000000, +0.00000000, +2.00000000
470+7.00000000, +1.00000000, +1.00000000
471transpose(sample) ! for comparison with affinInv.
472+7.00000000, -2.00000000, +9.00000000
473-5.00000000, -3.00000000, +1.00000000
474+9.00000000, +0.00000000, +2.00000000
475+7.00000000, +1.00000000, +1.00000000
476
477call setAffinity(affinity, sample, dim, tform, class = lowerUnit) ! whole sample transformation.
478transpose(affinity)
479+7.00000000, -2.00000000, +9.00000000
480-5.00000000, -3.00000000, +1.00000000
481+9.00000000, +0.00000000, +2.00000000
482+7.00000000, +1.00000000, +1.00000000
483do isam = 1, nsam
484 call setAffinity(affinity(:,isam), sample(:,isam), tform, lowerUnit) ! point-wise transformation
485end do
486transpose(affinity)
487+7.00000000, -2.00000000, +9.00000000
488-5.00000000, -3.00000000, +1.00000000
489+9.00000000, +0.00000000, +2.00000000
490+7.00000000, +1.00000000, +1.00000000
491call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerUnit) ! inverse transformation.
492transpose(affinInv)
493+7.00000000, -2.00000000, +9.00000000
494-5.00000000, -3.00000000, +1.00000000
495+9.00000000, +0.00000000, +2.00000000
496+7.00000000, +1.00000000, +1.00000000
497transpose(sample) ! for comparison with affinInv.
498+7.00000000, -2.00000000, +9.00000000
499-5.00000000, -3.00000000, +1.00000000
500+9.00000000, +0.00000000, +2.00000000
501+7.00000000, +1.00000000, +1.00000000
502
503dim = 2; ndim = getUnifRand(1, 3); nsam = getUnifRand(4, 5)
504[dim, ndim, nsam]
505+2, +3, +4
506tlate = getUnifRand(-5, +5, ndim)
507tlate
508-2.00000000, -4.00000000, +3.00000000
509sample = getUnifRand(-9, +9, ndim, nsam)
510sample
511-3.00000000, -8.00000000, -3.00000000, +5.00000000
512+4.00000000, +5.00000000, +4.00000000, -7.00000000
513-1.00000000, -9.00000000, +5.00000000, +3.00000000
514call setResized(affinity, shape(sample, IK))
515call setResized(affinInv, shape(sample, IK))
516tform = getCovRand(1., ndim)
517tform
518+1.00000000, +0.447494864, -0.525756478
519+0.447494864, +1.00000000, -0.206653520
520-0.525756478, -0.206653520, +1.00000000
521call setAffinity(affinity, sample, dim, tform, genrecmat) ! whole sample transformation.
522transpose(affinity)
523-0.684264064, +2.86416912, -0.249344647
524-1.03071737, +3.27992272, -5.82721567
525-3.83880281, +1.62424791, +5.75065517
526+0.290266633, -5.38248634, +1.81779242
527do isam = 1, nsam
528 call setAffinity(affinity(:,isam), sample(:,isam), tform, genrecmat) ! point-wise transformation
529end do
530transpose(affinity)
531-0.684264064, +2.86416912, -0.249344647
532-1.03071737, +3.27992272, -5.82721567
533-3.83880281, +1.62424791, +5.75065517
534+0.290266633, -5.38248634, +1.81779242
535call setAffinity(affinInv, affinity, dim, getMatInv(tform), genrecmat) ! inverse transformation.
536transpose(affinInv)
537-3.00000000, +3.99999976, -1.00000000
538-7.99999952, +4.99999905, -9.00000000
539-3.00000000, +3.99999976, +4.99999952
540+5.00000000, -6.99999952, +3.00000048
541transpose(sample) ! for comparison with affinInv.
542-3.00000000, +4.00000000, -1.00000000
543-8.00000000, +5.00000000, -9.00000000
544-3.00000000, +4.00000000, +5.00000000
545+5.00000000, -7.00000000, +3.00000000
546
547call setMatInit(tform(2:,1:ndim-1), lowDia, vlow = 0._RKG, vdia = 0._RKG) ! make tform an upper-diagonal matrix.
548tform
549+1.00000000, +0.447494864, -0.525756478
550+0.00000000, +1.00000000, -0.206653520
551+0.00000000, +0.00000000, +1.00000000
552call setAffinity(affinity, sample, dim, tform, class = upperDiag) ! whole sample transformation.
553transpose(affinity)
554-0.684264183, +4.20665359, -1.00000000
555-1.03071737, +6.85988140, -9.00000000
556-3.83880281, +2.96673250, +5.00000000
557+0.290266514, -7.61996078, +3.00000000
558do isam = 1, nsam
559 call setAffinity(affinity(:,isam), sample(:,isam), tform, upperDiag) ! point-wise transformation
560end do
561transpose(affinity)
562-0.684264183, +4.20665359, -1.00000000
563-1.03071737, +6.85988140, -9.00000000
564-3.83880281, +2.96673250, +5.00000000
565+0.290266514, -7.61996078, +3.00000000
566call setAffinity(affinInv, affinity, dim, getMatInv(tform, upperDiag), class = upperDiag) ! inverse transformation.
567transpose(affinInv)
568-3.00000000, +4.00000000, -1.00000000
569-7.99999952, +5.00000000, -9.00000000
570-3.00000000, +4.00000000, +5.00000000
571+5.00000000, -7.00000000, +3.00000000
572transpose(sample) ! for comparison with affinInv.
573-3.00000000, +4.00000000, -1.00000000
574-8.00000000, +5.00000000, -9.00000000
575-3.00000000, +4.00000000, +5.00000000
576+5.00000000, -7.00000000, +3.00000000
577
578call setAffinity(affinity, sample, dim, tform, class = upperUnit) ! whole sample transformation.
579transpose(affinity)
580-0.684264183, +4.20665359, -1.00000000
581-1.03071737, +6.85988140, -9.00000000
582-3.83880281, +2.96673250, +5.00000000
583+0.290266514, -7.61996078, +3.00000000
584do isam = 1, nsam
585 call setAffinity(affinity(:,isam), sample(:,isam), tform, upperUnit) ! point-wise transformation
586end do
587transpose(affinity)
588-0.684264183, +4.20665359, -1.00000000
589-1.03071737, +6.85988140, -9.00000000
590-3.83880281, +2.96673250, +5.00000000
591+0.290266514, -7.61996078, +3.00000000
592call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = upperUnit) ! inverse transformation.
593transpose(affinInv)
594-3.00000000, +4.00000000, -1.00000000
595-7.99999952, +5.00000000, -9.00000000
596-3.00000000, +4.00000000, +5.00000000
597+5.00000000, -7.00000000, +3.00000000
598transpose(sample) ! for comparison with affinInv.
599-3.00000000, +4.00000000, -1.00000000
600-8.00000000, +5.00000000, -9.00000000
601-3.00000000, +4.00000000, +5.00000000
602+5.00000000, -7.00000000, +3.00000000
603
604call setAffinity(affinity, sample, dim, tform, class = lowerDiag) ! whole sample transformation.
605transpose(affinity)
606-3.00000000, +4.00000000, -1.00000000
607-8.00000000, +5.00000000, -9.00000000
608-3.00000000, +4.00000000, +5.00000000
609+5.00000000, -7.00000000, +3.00000000
610do isam = 1, nsam
611 call setAffinity(affinity(:,isam), sample(:,isam), tform, lowerDiag) ! point-wise transformation
612end do
613transpose(affinity)
614-3.00000000, +4.00000000, -1.00000000
615-8.00000000, +5.00000000, -9.00000000
616-3.00000000, +4.00000000, +5.00000000
617+5.00000000, -7.00000000, +3.00000000
618call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerDiag) ! inverse transformation.
619transpose(affinInv)
620-3.00000000, +4.00000000, -1.00000000
621-8.00000000, +5.00000000, -9.00000000
622-3.00000000, +4.00000000, +5.00000000
623+5.00000000, -7.00000000, +3.00000000
624transpose(sample) ! for comparison with affinInv.
625-3.00000000, +4.00000000, -1.00000000
626-8.00000000, +5.00000000, -9.00000000
627-3.00000000, +4.00000000, +5.00000000
628+5.00000000, -7.00000000, +3.00000000
629
630call setAffinity(affinity, sample, dim, tform, class = lowerUnit) ! whole sample transformation.
631transpose(affinity)
632-3.00000000, +4.00000000, -1.00000000
633-8.00000000, +5.00000000, -9.00000000
634-3.00000000, +4.00000000, +5.00000000
635+5.00000000, -7.00000000, +3.00000000
636do isam = 1, nsam
637 call setAffinity(affinity(:,isam), sample(:,isam), tform, lowerUnit) ! point-wise transformation
638end do
639transpose(affinity)
640-3.00000000, +4.00000000, -1.00000000
641-8.00000000, +5.00000000, -9.00000000
642-3.00000000, +4.00000000, +5.00000000
643+5.00000000, -7.00000000, +3.00000000
644call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerUnit) ! inverse transformation.
645transpose(affinInv)
646-3.00000000, +4.00000000, -1.00000000
647-8.00000000, +5.00000000, -9.00000000
648-3.00000000, +4.00000000, +5.00000000
649+5.00000000, -7.00000000, +3.00000000
650transpose(sample) ! for comparison with affinInv.
651-3.00000000, +4.00000000, -1.00000000
652-8.00000000, +5.00000000, -9.00000000
653-3.00000000, +4.00000000, +5.00000000
654+5.00000000, -7.00000000, +3.00000000
655
656dim = 2; ndim = getUnifRand(1, 3); nsam = getUnifRand(4, 5)
657[dim, ndim, nsam]
658+2, +2, +4
659tlate = getUnifRand(-5, +5, ndim)
660tlate
661-2.00000000, -2.00000000
662sample = getUnifRand(-9, +9, ndim, nsam)
663sample
664+5.00000000, -4.00000000, +9.00000000, +6.00000000
665+9.00000000, -2.00000000, +4.00000000, +6.00000000
666call setResized(affinity, shape(sample, IK))
667call setResized(affinInv, shape(sample, IK))
668tform = getCovRand(1., ndim)
669tform
670+1.00000000, -0.981980860
671-0.981980860, +1.00000000
672call setAffinity(affinity, sample, dim, tform, genrecmat) ! whole sample transformation.
673transpose(affinity)
674-3.83782768, +4.09009552
675-2.03603840, +1.92792344
676+5.07207680, -4.83782768
677+0.108114719, +0.108114719
678do isam = 1, nsam
679 call setAffinity(affinity(:,isam), sample(:,isam), tform, genrecmat) ! point-wise transformation
680end do
681transpose(affinity)
682-3.83782768, +4.09009552
683-2.03603840, +1.92792344
684+5.07207680, -4.83782768
685+0.108114719, +0.108114719
686call setAffinity(affinInv, affinity, dim, getMatInv(tform), genrecmat) ! inverse transformation.
687transpose(affinInv)
688+5.00000000, +8.99999237
689-3.99999619, -1.99999619
690+8.99998474, +4.00000000
691+5.99998951, +5.99998951
692transpose(sample) ! for comparison with affinInv.
693+5.00000000, +9.00000000
694-4.00000000, -2.00000000
695+9.00000000, +4.00000000
696+6.00000000, +6.00000000
697
698call setMatInit(tform(2:,1:ndim-1), lowDia, vlow = 0._RKG, vdia = 0._RKG) ! make tform an upper-diagonal matrix.
699tform
700+1.00000000, -0.981980860
701+0.00000000, +1.00000000
702call setAffinity(affinity, sample, dim, tform, class = upperDiag) ! whole sample transformation.
703transpose(affinity)
704-3.83782768, +9.00000000
705-2.03603840, -2.00000000
706+5.07207680, +4.00000000
707+0.108114719, +6.00000000
708do isam = 1, nsam
709 call setAffinity(affinity(:,isam), sample(:,isam), tform, upperDiag) ! point-wise transformation
710end do
711transpose(affinity)
712-3.83782768, +9.00000000
713-2.03603840, -2.00000000
714+5.07207680, +4.00000000
715+0.108114719, +6.00000000
716call setAffinity(affinInv, affinity, dim, getMatInv(tform, upperDiag), class = upperDiag) ! inverse transformation.
717transpose(affinInv)
718+5.00000000, +9.00000000
719-4.00000000, -2.00000000
720+9.00000000, +4.00000000
721+6.00000000, +6.00000000
722transpose(sample) ! for comparison with affinInv.
723+5.00000000, +9.00000000
724-4.00000000, -2.00000000
725+9.00000000, +4.00000000
726+6.00000000, +6.00000000
727
728call setAffinity(affinity, sample, dim, tform, class = upperUnit) ! whole sample transformation.
729transpose(affinity)
730-3.83782768, +9.00000000
731-2.03603840, -2.00000000
732+5.07207680, +4.00000000
733+0.108114719, +6.00000000
734do isam = 1, nsam
735 call setAffinity(affinity(:,isam), sample(:,isam), tform, upperUnit) ! point-wise transformation
736end do
737transpose(affinity)
738-3.83782768, +9.00000000
739-2.03603840, -2.00000000
740+5.07207680, +4.00000000
741+0.108114719, +6.00000000
742call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = upperUnit) ! inverse transformation.
743transpose(affinInv)
744+5.00000000, +9.00000000
745-4.00000000, -2.00000000
746+9.00000000, +4.00000000
747+6.00000000, +6.00000000
748transpose(sample) ! for comparison with affinInv.
749+5.00000000, +9.00000000
750-4.00000000, -2.00000000
751+9.00000000, +4.00000000
752+6.00000000, +6.00000000
753
754call setAffinity(affinity, sample, dim, tform, class = lowerDiag) ! whole sample transformation.
755transpose(affinity)
756+5.00000000, +9.00000000
757-4.00000000, -2.00000000
758+9.00000000, +4.00000000
759+6.00000000, +6.00000000
760do isam = 1, nsam
761 call setAffinity(affinity(:,isam), sample(:,isam), tform, lowerDiag) ! point-wise transformation
762end do
763transpose(affinity)
764+5.00000000, +9.00000000
765-4.00000000, -2.00000000
766+9.00000000, +4.00000000
767+6.00000000, +6.00000000
768call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerDiag) ! inverse transformation.
769transpose(affinInv)
770+5.00000000, +9.00000000
771-4.00000000, -2.00000000
772+9.00000000, +4.00000000
773+6.00000000, +6.00000000
774transpose(sample) ! for comparison with affinInv.
775+5.00000000, +9.00000000
776-4.00000000, -2.00000000
777+9.00000000, +4.00000000
778+6.00000000, +6.00000000
779
780call setAffinity(affinity, sample, dim, tform, class = lowerUnit) ! whole sample transformation.
781transpose(affinity)
782+5.00000000, +9.00000000
783-4.00000000, -2.00000000
784+9.00000000, +4.00000000
785+6.00000000, +6.00000000
786do isam = 1, nsam
787 call setAffinity(affinity(:,isam), sample(:,isam), tform, lowerUnit) ! point-wise transformation
788end do
789transpose(affinity)
790+5.00000000, +9.00000000
791-4.00000000, -2.00000000
792+9.00000000, +4.00000000
793+6.00000000, +6.00000000
794call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerUnit) ! inverse transformation.
795transpose(affinInv)
796+5.00000000, +9.00000000
797-4.00000000, -2.00000000
798+9.00000000, +4.00000000
799+6.00000000, +6.00000000
800transpose(sample) ! for comparison with affinInv.
801+5.00000000, +9.00000000
802-4.00000000, -2.00000000
803+9.00000000, +4.00000000
804+6.00000000, +6.00000000
805
806!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
807!Transform 2D sample along the first dimension.
808!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
809
810
811dim = 1; ndim = getUnifRand(1, 3); nsam = getUnifRand(4, 5)
812[dim, ndim, nsam]
813+1, +2, +5
814tlate = getUnifRand(-5, +5, ndim)
815tlate
816+5.00000000, -5.00000000
817sample = getUnifRand(-9, +9, nsam, ndim)
818sample
819-4.00000000, +7.00000000
820+2.00000000, +8.00000000
821-4.00000000, +4.00000000
822-1.00000000, +5.00000000
823-8.00000000, -6.00000000
824call setResized(affinity, shape(sample, IK))
825call setResized(affinInv, shape(sample, IK))
826tform = getCovRand(1., ndim)
827tform
828+1.00000000, +0.165688619E-1
829+0.165688619E-1, +1.00000000
830call setAffinity(affinity, sample, dim, tform) ! whole sample transformation.
831affinity
832-3.88401794, +6.93372440
833+2.13255095, +8.03313732
834-3.93372464, +3.93372464
835-0.917155683, +4.98343134
836-8.09941292, -6.13255072
837do isam = 1, nsam
838 call setAffinity(affinity(isam,:), sample(isam,:), tform, genrecmat) ! point-wise transformation
839end do
840affinity
841-3.88401794, +6.93372440
842+2.13255095, +8.03313732
843-3.93372464, +3.93372464
844-0.917155683, +4.98343134
845-8.09941292, -6.13255072
846call setAffinity(affinInv, affinity, dim, getMatInv(tform), genrecmat) ! inverse transformation.
847affinInv
848-4.00000000, +7.00000048
849+2.00000024, +8.00000000
850-4.00000048, +4.00000048
851-1.00000000, +5.00000048
852-8.00000095, -6.00000000
853sample ! for comparison with affinInv.
854-4.00000000, +7.00000000
855+2.00000000, +8.00000000
856-4.00000000, +4.00000000
857-1.00000000, +5.00000000
858-8.00000000, -6.00000000
859
860call setMatInit(tform(2:,1:ndim-1), lowDia, vlow = 0._RKG, vdia = 0._RKG) ! make tform an upper-diagonal matrix.
861tform
862+1.00000000, +0.165688619E-1
863+0.00000000, +1.00000000
864call setAffinity(affinity, sample, dim, tform, class = upperDiag) ! whole sample transformation.
865affinity
866-3.88401794, +7.00000000
867+2.13255095, +8.00000000
868-3.93372464, +4.00000000
869-0.917155683, +5.00000000
870-8.09941292, -6.00000000
871do isam = 1, nsam
872 call setAffinity(affinity(isam,:), sample(isam,:), tform, upperDiag) ! point-wise transformation
873end do
874affinity
875-3.88401794, +7.00000000
876+2.13255095, +8.00000000
877-3.93372464, +4.00000000
878-0.917155683, +5.00000000
879-8.09941292, -6.00000000
880call setAffinity(affinInv, affinity, dim, getMatInv(tform, upperDiag), class = upperDiag) ! inverse transformation.
881affinInv
882-4.00000000, +7.00000000
883+2.00000000, +8.00000000
884-4.00000000, +4.00000000
885-1.00000000, +5.00000000
886-7.99999952, -6.00000000
887sample ! for comparison with affinInv.
888-4.00000000, +7.00000000
889+2.00000000, +8.00000000
890-4.00000000, +4.00000000
891-1.00000000, +5.00000000
892-8.00000000, -6.00000000
893
894call setAffinity(affinity, sample, dim, tform, class = upperUnit) ! whole sample transformation.
895affinity
896-3.88401794, +7.00000000
897+2.13255095, +8.00000000
898-3.93372464, +4.00000000
899-0.917155683, +5.00000000
900-8.09941292, -6.00000000
901do isam = 1, nsam
902 call setAffinity(affinity(isam,:), sample(isam,:), tform, upperUnit) ! point-wise transformation
903end do
904affinity
905-3.88401794, +7.00000000
906+2.13255095, +8.00000000
907-3.93372464, +4.00000000
908-0.917155683, +5.00000000
909-8.09941292, -6.00000000
910call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = upperUnit) ! inverse transformation.
911affinInv
912-4.00000000, +7.00000000
913+2.00000000, +8.00000000
914-4.00000000, +4.00000000
915-1.00000000, +5.00000000
916-7.99999952, -6.00000000
917sample ! for comparison with affinInv.
918-4.00000000, +7.00000000
919+2.00000000, +8.00000000
920-4.00000000, +4.00000000
921-1.00000000, +5.00000000
922-8.00000000, -6.00000000
923
924call setAffinity(affinity, sample, dim, tform, class = lowerDiag) ! whole sample transformation.
925affinity
926-4.00000000, +7.00000000
927+2.00000000, +8.00000000
928-4.00000000, +4.00000000
929-1.00000000, +5.00000000
930-8.00000000, -6.00000000
931do isam = 1, nsam
932 call setAffinity(affinity(isam,:), sample(isam,:), tform, lowerDiag) ! point-wise transformation
933end do
934affinity
935-4.00000000, +7.00000000
936+2.00000000, +8.00000000
937-4.00000000, +4.00000000
938-1.00000000, +5.00000000
939-8.00000000, -6.00000000
940call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerDiag) ! inverse transformation.
941affinInv
942-4.00000000, +7.00000000
943+2.00000000, +8.00000000
944-4.00000000, +4.00000000
945-1.00000000, +5.00000000
946-8.00000000, -6.00000000
947sample ! for comparison with affinInv.
948-4.00000000, +7.00000000
949+2.00000000, +8.00000000
950-4.00000000, +4.00000000
951-1.00000000, +5.00000000
952-8.00000000, -6.00000000
953
954call setAffinity(affinity, sample, dim, tform, class = lowerUnit) ! whole sample transformation.
955affinity
956-4.00000000, +7.00000000
957+2.00000000, +8.00000000
958-4.00000000, +4.00000000
959-1.00000000, +5.00000000
960-8.00000000, -6.00000000
961do isam = 1, nsam
962 call setAffinity(affinity(isam,:), sample(isam,:), tform, lowerUnit) ! point-wise transformation
963end do
964affinity
965-4.00000000, +7.00000000
966+2.00000000, +8.00000000
967-4.00000000, +4.00000000
968-1.00000000, +5.00000000
969-8.00000000, -6.00000000
970call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerUnit) ! inverse transformation.
971affinInv
972-4.00000000, +7.00000000
973+2.00000000, +8.00000000
974-4.00000000, +4.00000000
975-1.00000000, +5.00000000
976-8.00000000, -6.00000000
977sample ! for comparison with affinInv.
978-4.00000000, +7.00000000
979+2.00000000, +8.00000000
980-4.00000000, +4.00000000
981-1.00000000, +5.00000000
982-8.00000000, -6.00000000
983
984dim = 1; ndim = getUnifRand(1, 3); nsam = getUnifRand(4, 5)
985[dim, ndim, nsam]
986+1, +2, +4
987tlate = getUnifRand(-5, +5, ndim)
988tlate
989-1.00000000, +1.00000000
990sample = getUnifRand(-9, +9, nsam, ndim)
991sample
992+2.00000000, +1.00000000
993-1.00000000, -1.00000000
994-1.00000000, +6.00000000
995+5.00000000, -4.00000000
996call setResized(affinity, shape(sample, IK))
997call setResized(affinInv, shape(sample, IK))
998tform = getCovRand(1., ndim)
999tform
1000+1.00000000, -0.720322132
1001-0.720322132, +1.00000000
1002call setAffinity(affinity, sample, dim, tform) ! whole sample transformation.
1003affinity
1004+1.27967787, -0.440644264
1005-0.279677868, -0.279677868
1006-5.32193279, +6.72032213
1007+7.88128853, -7.60161066
1008do isam = 1, nsam
1009 call setAffinity(affinity(isam,:), sample(isam,:), tform, genrecmat) ! point-wise transformation
1010end do
1011affinity
1012+1.27967787, -0.440644264
1013-0.279677868, -0.279677868
1014-5.32193279, +6.72032213
1015+7.88128853, -7.60161066
1016call setAffinity(affinInv, affinity, dim, getMatInv(tform), genrecmat) ! inverse transformation.
1017affinInv
1018+2.00000000, +1.00000024
1019-1.00000000, -1.00000000
1020-1.00000095, +5.99999905
1021+5.00000191, -3.99999809
1022sample ! for comparison with affinInv.
1023+2.00000000, +1.00000000
1024-1.00000000, -1.00000000
1025-1.00000000, +6.00000000
1026+5.00000000, -4.00000000
1027
1028call setMatInit(tform(2:,1:ndim-1), lowDia, vlow = 0._RKG, vdia = 0._RKG) ! make tform an upper-diagonal matrix.
1029tform
1030+1.00000000, -0.720322132
1031+0.00000000, +1.00000000
1032call setAffinity(affinity, sample, dim, tform, class = upperDiag) ! whole sample transformation.
1033affinity
1034+1.27967787, +1.00000000
1035-0.279677868, -1.00000000
1036-5.32193279, +6.00000000
1037+7.88128853, -4.00000000
1038do isam = 1, nsam
1039 call setAffinity(affinity(isam,:), sample(isam,:), tform, upperDiag) ! point-wise transformation
1040end do
1041affinity
1042+1.27967787, +1.00000000
1043-0.279677868, -1.00000000
1044-5.32193279, +6.00000000
1045+7.88128853, -4.00000000
1046call setAffinity(affinInv, affinity, dim, getMatInv(tform, upperDiag), class = upperDiag) ! inverse transformation.
1047affinInv
1048+2.00000000, +1.00000000
1049-1.00000000, -1.00000000
1050-1.00000000, +6.00000000
1051+5.00000000, -4.00000000
1052sample ! for comparison with affinInv.
1053+2.00000000, +1.00000000
1054-1.00000000, -1.00000000
1055-1.00000000, +6.00000000
1056+5.00000000, -4.00000000
1057
1058call setAffinity(affinity, sample, dim, tform, class = upperUnit) ! whole sample transformation.
1059affinity
1060+1.27967787, +1.00000000
1061-0.279677868, -1.00000000
1062-5.32193279, +6.00000000
1063+7.88128853, -4.00000000
1064do isam = 1, nsam
1065 call setAffinity(affinity(isam,:), sample(isam,:), tform, upperUnit) ! point-wise transformation
1066end do
1067affinity
1068+1.27967787, +1.00000000
1069-0.279677868, -1.00000000
1070-5.32193279, +6.00000000
1071+7.88128853, -4.00000000
1072call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = upperUnit) ! inverse transformation.
1073affinInv
1074+2.00000000, +1.00000000
1075-1.00000000, -1.00000000
1076-1.00000000, +6.00000000
1077+5.00000000, -4.00000000
1078sample ! for comparison with affinInv.
1079+2.00000000, +1.00000000
1080-1.00000000, -1.00000000
1081-1.00000000, +6.00000000
1082+5.00000000, -4.00000000
1083
1084call setAffinity(affinity, sample, dim, tform, class = lowerDiag) ! whole sample transformation.
1085affinity
1086+2.00000000, +1.00000000
1087-1.00000000, -1.00000000
1088-1.00000000, +6.00000000
1089+5.00000000, -4.00000000
1090do isam = 1, nsam
1091 call setAffinity(affinity(isam,:), sample(isam,:), tform, lowerDiag) ! point-wise transformation
1092end do
1093affinity
1094+2.00000000, +1.00000000
1095-1.00000000, -1.00000000
1096-1.00000000, +6.00000000
1097+5.00000000, -4.00000000
1098call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerDiag) ! inverse transformation.
1099affinInv
1100+2.00000000, +1.00000000
1101-1.00000000, -1.00000000
1102-1.00000000, +6.00000000
1103+5.00000000, -4.00000000
1104sample ! for comparison with affinInv.
1105+2.00000000, +1.00000000
1106-1.00000000, -1.00000000
1107-1.00000000, +6.00000000
1108+5.00000000, -4.00000000
1109
1110call setAffinity(affinity, sample, dim, tform, class = lowerUnit) ! whole sample transformation.
1111affinity
1112+2.00000000, +1.00000000
1113-1.00000000, -1.00000000
1114-1.00000000, +6.00000000
1115+5.00000000, -4.00000000
1116do isam = 1, nsam
1117 call setAffinity(affinity(isam,:), sample(isam,:), tform, lowerUnit) ! point-wise transformation
1118end do
1119affinity
1120+2.00000000, +1.00000000
1121-1.00000000, -1.00000000
1122-1.00000000, +6.00000000
1123+5.00000000, -4.00000000
1124call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerUnit) ! inverse transformation.
1125affinInv
1126+2.00000000, +1.00000000
1127-1.00000000, -1.00000000
1128-1.00000000, +6.00000000
1129+5.00000000, -4.00000000
1130sample ! for comparison with affinInv.
1131+2.00000000, +1.00000000
1132-1.00000000, -1.00000000
1133-1.00000000, +6.00000000
1134+5.00000000, -4.00000000
1135
1136dim = 1; ndim = getUnifRand(1, 3); nsam = getUnifRand(4, 5)
1137[dim, ndim, nsam]
1138+1, +2, +4
1139tlate = getUnifRand(-5, +5, ndim)
1140tlate
1141+5.00000000, -5.00000000
1142sample = getUnifRand(-9, +9, nsam, ndim)
1143sample
1144+1.00000000, -3.00000000
1145+6.00000000, -6.00000000
1146+6.00000000, -9.00000000
1147-6.00000000, +3.00000000
1148call setResized(affinity, shape(sample, IK))
1149call setResized(affinInv, shape(sample, IK))
1150tform = getCovRand(1., ndim)
1151tform
1152+1.00000000, -0.942091107
1153-0.942091107, +1.00000000
1154call setAffinity(affinity, sample, dim, tform) ! whole sample transformation.
1155affinity
1156+3.82627344, -3.94209099
1157+11.6525469, -11.6525469
1158+14.4788198, -14.6525469
1159-8.82627296, +8.65254688
1160do isam = 1, nsam
1161 call setAffinity(affinity(isam,:), sample(isam,:), tform, genrecmat) ! point-wise transformation
1162end do
1163affinity
1164+3.82627344, -3.94209099
1165+11.6525469, -11.6525469
1166+14.4788198, -14.6525469
1167-8.82627296, +8.65254688
1168call setAffinity(affinInv, affinity, dim, getMatInv(tform), genrecmat) ! inverse transformation.
1169affinInv
1170+1.00000381, -2.99999619
1171+6.00000000, -6.00000000
1172+6.00000000, -9.00000000
1173-5.99999237, +3.00000000
1174sample ! for comparison with affinInv.
1175+1.00000000, -3.00000000
1176+6.00000000, -6.00000000
1177+6.00000000, -9.00000000
1178-6.00000000, +3.00000000
1179
1180call setMatInit(tform(2:,1:ndim-1), lowDia, vlow = 0._RKG, vdia = 0._RKG) ! make tform an upper-diagonal matrix.
1181tform
1182+1.00000000, -0.942091107
1183+0.00000000, +1.00000000
1184call setAffinity(affinity, sample, dim, tform, class = upperDiag) ! whole sample transformation.
1185affinity
1186+3.82627344, -3.00000000
1187+11.6525469, -6.00000000
1188+14.4788198, -9.00000000
1189-8.82627296, +3.00000000
1190do isam = 1, nsam
1191 call setAffinity(affinity(isam,:), sample(isam,:), tform, upperDiag) ! point-wise transformation
1192end do
1193affinity
1194+3.82627344, -3.00000000
1195+11.6525469, -6.00000000
1196+14.4788198, -9.00000000
1197-8.82627296, +3.00000000
1198call setAffinity(affinInv, affinity, dim, getMatInv(tform, upperDiag), class = upperDiag) ! inverse transformation.
1199affinInv
1200+1.00000000, -3.00000000
1201+6.00000000, -6.00000000
1202+6.00000000, -9.00000000
1203-5.99999952, +3.00000000
1204sample ! for comparison with affinInv.
1205+1.00000000, -3.00000000
1206+6.00000000, -6.00000000
1207+6.00000000, -9.00000000
1208-6.00000000, +3.00000000
1209
1210call setAffinity(affinity, sample, dim, tform, class = upperUnit) ! whole sample transformation.
1211affinity
1212+3.82627344, -3.00000000
1213+11.6525469, -6.00000000
1214+14.4788198, -9.00000000
1215-8.82627296, +3.00000000
1216do isam = 1, nsam
1217 call setAffinity(affinity(isam,:), sample(isam,:), tform, upperUnit) ! point-wise transformation
1218end do
1219affinity
1220+3.82627344, -3.00000000
1221+11.6525469, -6.00000000
1222+14.4788198, -9.00000000
1223-8.82627296, +3.00000000
1224call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = upperUnit) ! inverse transformation.
1225affinInv
1226+1.00000000, -3.00000000
1227+6.00000000, -6.00000000
1228+6.00000000, -9.00000000
1229-5.99999952, +3.00000000
1230sample ! for comparison with affinInv.
1231+1.00000000, -3.00000000
1232+6.00000000, -6.00000000
1233+6.00000000, -9.00000000
1234-6.00000000, +3.00000000
1235
1236call setAffinity(affinity, sample, dim, tform, class = lowerDiag) ! whole sample transformation.
1237affinity
1238+1.00000000, -3.00000000
1239+6.00000000, -6.00000000
1240+6.00000000, -9.00000000
1241-6.00000000, +3.00000000
1242do isam = 1, nsam
1243 call setAffinity(affinity(isam,:), sample(isam,:), tform, lowerDiag) ! point-wise transformation
1244end do
1245affinity
1246+1.00000000, -3.00000000
1247+6.00000000, -6.00000000
1248+6.00000000, -9.00000000
1249-6.00000000, +3.00000000
1250call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerDiag) ! inverse transformation.
1251affinInv
1252+1.00000000, -3.00000000
1253+6.00000000, -6.00000000
1254+6.00000000, -9.00000000
1255-6.00000000, +3.00000000
1256sample ! for comparison with affinInv.
1257+1.00000000, -3.00000000
1258+6.00000000, -6.00000000
1259+6.00000000, -9.00000000
1260-6.00000000, +3.00000000
1261
1262call setAffinity(affinity, sample, dim, tform, class = lowerUnit) ! whole sample transformation.
1263affinity
1264+1.00000000, -3.00000000
1265+6.00000000, -6.00000000
1266+6.00000000, -9.00000000
1267-6.00000000, +3.00000000
1268do isam = 1, nsam
1269 call setAffinity(affinity(isam,:), sample(isam,:), tform, lowerUnit) ! point-wise transformation
1270end do
1271affinity
1272+1.00000000, -3.00000000
1273+6.00000000, -6.00000000
1274+6.00000000, -9.00000000
1275-6.00000000, +3.00000000
1276call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerUnit) ! inverse transformation.
1277affinInv
1278+1.00000000, -3.00000000
1279+6.00000000, -6.00000000
1280+6.00000000, -9.00000000
1281-6.00000000, +3.00000000
1282sample ! for comparison with affinInv.
1283+1.00000000, -3.00000000
1284+6.00000000, -6.00000000
1285+6.00000000, -9.00000000
1286-6.00000000, +3.00000000
1287
1288dim = 1; ndim = getUnifRand(1, 3); nsam = getUnifRand(4, 5)
1289[dim, ndim, nsam]
1290+1, +3, +4
1291tlate = getUnifRand(-5, +5, ndim)
1292tlate
1293-3.00000000, -2.00000000, -5.00000000
1294sample = getUnifRand(-9, +9, nsam, ndim)
1295sample
1296-6.00000000, +3.00000000, +9.00000000
1297-2.00000000, -8.00000000, -5.00000000
1298-3.00000000, +0.00000000, +0.00000000
1299-2.00000000, -8.00000000, +6.00000000
1300call setResized(affinity, shape(sample, IK))
1301call setResized(affinInv, shape(sample, IK))
1302tform = getCovRand(1., ndim)
1303tform
1304+1.00000000, -0.989137650, -0.713971853
1305-0.989137650, +1.00000000, +0.771194875
1306-0.713971853, +0.771194875, +1.00000000
1307call setAffinity(affinity, sample, dim, tform) ! whole sample transformation.
1308affinity
1309-15.3931599, +15.8755798, +15.5974159
1310+9.48296070, -9.87769890, -9.74161530
1311-3.00000000, +2.96741295, +2.14191556
1312+1.62927008, -1.39455557, +1.25838470
1313do isam = 1, nsam
1314 call setAffinity(affinity(isam,:), sample(isam,:), tform, genrecmat) ! point-wise transformation
1315end do
1316affinity
1317-15.3931599, +15.8755798, +15.5974159
1318+9.48296070, -9.87769890, -9.74161530
1319-3.00000000, +2.96741295, +2.14191556
1320+1.62927008, -1.39455557, +1.25838470
1321call setAffinity(affinInv, affinity, dim, getMatInv(tform), genrecmat) ! inverse transformation.
1322affinInv
1323-5.99999237, +3.00003052, +9.00000381
1324-2.00000000, -8.00000000, -5.00001144
1325-2.99999809, -0.572204590E-5, +0.476837158E-6
1326-2.00002384, -8.00001431, +6.00000143
1327sample ! for comparison with affinInv.
1328-6.00000000, +3.00000000, +9.00000000
1329-2.00000000, -8.00000000, -5.00000000
1330-3.00000000, +0.00000000, +0.00000000
1331-2.00000000, -8.00000000, +6.00000000
1332
1333call setMatInit(tform(2:,1:ndim-1), lowDia, vlow = 0._RKG, vdia = 0._RKG) ! make tform an upper-diagonal matrix.
1334tform
1335+1.00000000, -0.989137650, -0.713971853
1336+0.00000000, +1.00000000, +0.771194875
1337+0.00000000, +0.00000000, +1.00000000
1338call setAffinity(affinity, sample, dim, tform, class = upperDiag) ! whole sample transformation.
1339affinity
1340-15.3931599, +9.94075394, +9.00000000
1341+9.48296070, -11.8559742, -5.00000000
1342-3.00000000, +0.00000000, +0.00000000
1343+1.62927008, -3.37283087, +6.00000000
1344do isam = 1, nsam
1345 call setAffinity(affinity(isam,:), sample(isam,:), tform, upperDiag) ! point-wise transformation
1346end do
1347affinity
1348-15.3931599, +9.94075394, +9.00000000
1349+9.48296070, -11.8559742, -5.00000000
1350-3.00000000, +0.00000000, +0.00000000
1351+1.62927008, -3.37283087, +6.00000000
1352call setAffinity(affinInv, affinity, dim, getMatInv(tform, upperDiag), class = upperDiag) ! inverse transformation.
1353affinInv
1354-6.00000000, +3.00000000, +9.00000000
1355-1.99999905, -8.00000000, -5.00000000
1356-3.00000000, +0.00000000, +0.00000000
1357-2.00000000, -8.00000000, +6.00000000
1358sample ! for comparison with affinInv.
1359-6.00000000, +3.00000000, +9.00000000
1360-2.00000000, -8.00000000, -5.00000000
1361-3.00000000, +0.00000000, +0.00000000
1362-2.00000000, -8.00000000, +6.00000000
1363
1364call setAffinity(affinity, sample, dim, tform, class = upperUnit) ! whole sample transformation.
1365affinity
1366-15.3931599, +9.94075394, +9.00000000
1367+9.48296070, -11.8559742, -5.00000000
1368-3.00000000, +0.00000000, +0.00000000
1369+1.62927008, -3.37283087, +6.00000000
1370do isam = 1, nsam
1371 call setAffinity(affinity(isam,:), sample(isam,:), tform, upperUnit) ! point-wise transformation
1372end do
1373affinity
1374-15.3931599, +9.94075394, +9.00000000
1375+9.48296070, -11.8559742, -5.00000000
1376-3.00000000, +0.00000000, +0.00000000
1377+1.62927008, -3.37283087, +6.00000000
1378call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = upperUnit) ! inverse transformation.
1379affinInv
1380-6.00000000, +3.00000000, +9.00000000
1381-1.99999905, -8.00000000, -5.00000000
1382-3.00000000, +0.00000000, +0.00000000
1383-2.00000000, -8.00000000, +6.00000000
1384sample ! for comparison with affinInv.
1385-6.00000000, +3.00000000, +9.00000000
1386-2.00000000, -8.00000000, -5.00000000
1387-3.00000000, +0.00000000, +0.00000000
1388-2.00000000, -8.00000000, +6.00000000
1389
1390call setAffinity(affinity, sample, dim, tform, class = lowerDiag) ! whole sample transformation.
1391affinity
1392-6.00000000, +3.00000000, +9.00000000
1393-2.00000000, -8.00000000, -5.00000000
1394-3.00000000, +0.00000000, +0.00000000
1395-2.00000000, -8.00000000, +6.00000000
1396do isam = 1, nsam
1397 call setAffinity(affinity(isam,:), sample(isam,:), tform, lowerDiag) ! point-wise transformation
1398end do
1399affinity
1400-6.00000000, +3.00000000, +9.00000000
1401-2.00000000, -8.00000000, -5.00000000
1402-3.00000000, +0.00000000, +0.00000000
1403-2.00000000, -8.00000000, +6.00000000
1404call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerDiag) ! inverse transformation.
1405affinInv
1406-6.00000000, +3.00000000, +9.00000000
1407-2.00000000, -8.00000000, -5.00000000
1408-3.00000000, +0.00000000, +0.00000000
1409-2.00000000, -8.00000000, +6.00000000
1410sample ! for comparison with affinInv.
1411-6.00000000, +3.00000000, +9.00000000
1412-2.00000000, -8.00000000, -5.00000000
1413-3.00000000, +0.00000000, +0.00000000
1414-2.00000000, -8.00000000, +6.00000000
1415
1416call setAffinity(affinity, sample, dim, tform, class = lowerUnit) ! whole sample transformation.
1417affinity
1418-6.00000000, +3.00000000, +9.00000000
1419-2.00000000, -8.00000000, -5.00000000
1420-3.00000000, +0.00000000, +0.00000000
1421-2.00000000, -8.00000000, +6.00000000
1422do isam = 1, nsam
1423 call setAffinity(affinity(isam,:), sample(isam,:), tform, lowerUnit) ! point-wise transformation
1424end do
1425affinity
1426-6.00000000, +3.00000000, +9.00000000
1427-2.00000000, -8.00000000, -5.00000000
1428-3.00000000, +0.00000000, +0.00000000
1429-2.00000000, -8.00000000, +6.00000000
1430call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerUnit) ! inverse transformation.
1431affinInv
1432-6.00000000, +3.00000000, +9.00000000
1433-2.00000000, -8.00000000, -5.00000000
1434-3.00000000, +0.00000000, +0.00000000
1435-2.00000000, -8.00000000, +6.00000000
1436sample ! for comparison with affinInv.
1437-6.00000000, +3.00000000, +9.00000000
1438-2.00000000, -8.00000000, -5.00000000
1439-3.00000000, +0.00000000, +0.00000000
1440-2.00000000, -8.00000000, +6.00000000
1441
1442dim = 1; ndim = getUnifRand(1, 3); nsam = getUnifRand(4, 5)
1443[dim, ndim, nsam]
1444+1, +3, +4
1445tlate = getUnifRand(-5, +5, ndim)
1446tlate
1447+4.00000000, +2.00000000, -5.00000000
1448sample = getUnifRand(-9, +9, nsam, ndim)
1449sample
1450+7.00000000, +2.00000000, +7.00000000
1451-1.00000000, +8.00000000, -7.00000000
1452+8.00000000, -2.00000000, +8.00000000
1453+7.00000000, -5.00000000, +8.00000000
1454call setResized(affinity, shape(sample, IK))
1455call setResized(affinInv, shape(sample, IK))
1456tform = getCovRand(1., ndim)
1457tform
1458+1.00000000, +0.996356666, -0.705461860
1459+0.996356666, +1.00000000, -0.645069540
1460-0.705461860, -0.645069540, +1.00000000
1461call setAffinity(affinity, sample, dim, tform) ! whole sample transformation.
1462affinity
1463+4.05448008, +4.45901012, +0.771627903
1464+11.9090862, +11.5191307, -11.4550943
1465+0.363591671, +0.810297012, +3.64644432
1466-3.62547827, -3.18605947, +6.28711510
1467do isam = 1, nsam
1468 call setAffinity(affinity(isam,:), sample(isam,:), tform, genrecmat) ! point-wise transformation
1469end do
1470affinity
1471+4.05448008, +4.45901012, +0.771627903
1472+11.9090862, +11.5191307, -11.4550943
1473+0.363591671, +0.810297012, +3.64644432
1474-3.62547827, -3.18605947, +6.28711510
1475call setAffinity(affinInv, affinity, dim, getMatInv(tform), genrecmat) ! inverse transformation.
1476affinInv
1477+6.99862671, +2.00109863, +6.99986267
1478-1.00341797, +8.00366211, -7.00054932
1479+8.00006104, -2.00012207, +8.00001526
1480+7.00048828, -5.00097656, +8.00004578
1481sample ! for comparison with affinInv.
1482+7.00000000, +2.00000000, +7.00000000
1483-1.00000000, +8.00000000, -7.00000000
1484+8.00000000, -2.00000000, +8.00000000
1485+7.00000000, -5.00000000, +8.00000000
1486
1487call setMatInit(tform(2:,1:ndim-1), lowDia, vlow = 0._RKG, vdia = 0._RKG) ! make tform an upper-diagonal matrix.
1488tform
1489+1.00000000, +0.996356666, -0.705461860
1490+0.00000000, +1.00000000, -0.645069540
1491+0.00000000, +0.00000000, +1.00000000
1492call setAffinity(affinity, sample, dim, tform, class = upperDiag) ! whole sample transformation.
1493affinity
1494+4.05448055, -2.51548672, +7.00000000
1495+11.9090862, +12.5154867, -7.00000000
1496+0.363591671, -7.16055632, +8.00000000
1497-3.62547874, -10.1605568, +8.00000000
1498do isam = 1, nsam
1499 call setAffinity(affinity(isam,:), sample(isam,:), tform, upperDiag) ! point-wise transformation
1500end do
1501affinity
1502+4.05448055, -2.51548672, +7.00000000
1503+11.9090862, +12.5154867, -7.00000000
1504+0.363591671, -7.16055632, +8.00000000
1505-3.62547874, -10.1605568, +8.00000000
1506call setAffinity(affinInv, affinity, dim, getMatInv(tform, upperDiag), class = upperDiag) ! inverse transformation.
1507affinInv
1508+7.00000000, +2.00000000, +7.00000000
1509-1.00000000, +8.00000000, -7.00000000
1510+8.00000000, -2.00000000, +8.00000000
1511+7.00000000, -5.00000048, +8.00000000
1512sample ! for comparison with affinInv.
1513+7.00000000, +2.00000000, +7.00000000
1514-1.00000000, +8.00000000, -7.00000000
1515+8.00000000, -2.00000000, +8.00000000
1516+7.00000000, -5.00000000, +8.00000000
1517
1518call setAffinity(affinity, sample, dim, tform, class = upperUnit) ! whole sample transformation.
1519affinity
1520+4.05448055, -2.51548672, +7.00000000
1521+11.9090862, +12.5154867, -7.00000000
1522+0.363591671, -7.16055632, +8.00000000
1523-3.62547874, -10.1605568, +8.00000000
1524do isam = 1, nsam
1525 call setAffinity(affinity(isam,:), sample(isam,:), tform, upperUnit) ! point-wise transformation
1526end do
1527affinity
1528+4.05448055, -2.51548672, +7.00000000
1529+11.9090862, +12.5154867, -7.00000000
1530+0.363591671, -7.16055632, +8.00000000
1531-3.62547874, -10.1605568, +8.00000000
1532call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = upperUnit) ! inverse transformation.
1533affinInv
1534+7.00000000, +2.00000000, +7.00000000
1535-1.00000000, +8.00000000, -7.00000000
1536+8.00000000, -2.00000000, +8.00000000
1537+7.00000000, -5.00000048, +8.00000000
1538sample ! for comparison with affinInv.
1539+7.00000000, +2.00000000, +7.00000000
1540-1.00000000, +8.00000000, -7.00000000
1541+8.00000000, -2.00000000, +8.00000000
1542+7.00000000, -5.00000000, +8.00000000
1543
1544call setAffinity(affinity, sample, dim, tform, class = lowerDiag) ! whole sample transformation.
1545affinity
1546+7.00000000, +2.00000000, +7.00000000
1547-1.00000000, +8.00000000, -7.00000000
1548+8.00000000, -2.00000000, +8.00000000
1549+7.00000000, -5.00000000, +8.00000000
1550do isam = 1, nsam
1551 call setAffinity(affinity(isam,:), sample(isam,:), tform, lowerDiag) ! point-wise transformation
1552end do
1553affinity
1554+7.00000000, +2.00000000, +7.00000000
1555-1.00000000, +8.00000000, -7.00000000
1556+8.00000000, -2.00000000, +8.00000000
1557+7.00000000, -5.00000000, +8.00000000
1558call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerDiag) ! inverse transformation.
1559affinInv
1560+7.00000000, +2.00000000, +7.00000000
1561-1.00000000, +8.00000000, -7.00000000
1562+8.00000000, -2.00000000, +8.00000000
1563+7.00000000, -5.00000000, +8.00000000
1564sample ! for comparison with affinInv.
1565+7.00000000, +2.00000000, +7.00000000
1566-1.00000000, +8.00000000, -7.00000000
1567+8.00000000, -2.00000000, +8.00000000
1568+7.00000000, -5.00000000, +8.00000000
1569
1570call setAffinity(affinity, sample, dim, tform, class = lowerUnit) ! whole sample transformation.
1571affinity
1572+7.00000000, +2.00000000, +7.00000000
1573-1.00000000, +8.00000000, -7.00000000
1574+8.00000000, -2.00000000, +8.00000000
1575+7.00000000, -5.00000000, +8.00000000
1576do isam = 1, nsam
1577 call setAffinity(affinity(isam,:), sample(isam,:), tform, lowerUnit) ! point-wise transformation
1578end do
1579affinity
1580+7.00000000, +2.00000000, +7.00000000
1581-1.00000000, +8.00000000, -7.00000000
1582+8.00000000, -2.00000000, +8.00000000
1583+7.00000000, -5.00000000, +8.00000000
1584call setAffinity(affinInv, affinity, dim, getMatInv(tform), class = lowerUnit) ! inverse transformation.
1585affinInv
1586+7.00000000, +2.00000000, +7.00000000
1587-1.00000000, +8.00000000, -7.00000000
1588+8.00000000, -2.00000000, +8.00000000
1589+7.00000000, -5.00000000, +8.00000000
1590sample ! for comparison with affinInv.
1591+7.00000000, +2.00000000, +7.00000000
1592-1.00000000, +8.00000000, -7.00000000
1593+8.00000000, -2.00000000, +8.00000000
1594+7.00000000, -5.00000000, +8.00000000
1595

Postprocessing of the example output
1#!/usr/bin/env python
2
3import matplotlib.pyplot as plt
4import pandas as pd
5import numpy as np
6import glob
7import sys
8
9linewidth = 2
10fontsize = 17
11
12for kind in ["warp", "rotation"]:
13
14 pattern = "*."+kind+"*.txt"
15 fileList = glob.glob(pattern)
16
17 for file in fileList:
18
19 df = pd.read_csv(file, delimiter = ",", header = None)
20
21 # definitions for the axes
22 left, width = 0.1, 0.65
23 bottom, height = 0.1, 0.65
24 spacing = 0.015
25
26 # start with a square Figure
27 fig = plt.figure(figsize = (8, 8))
28
29 plt.rcParams.update({'font.size': fontsize - 2})
30 ax = fig.add_axes([left, bottom, width, height]) # scatter plot
31 ax_histx = fig.add_axes([left, bottom + height + spacing, width, 0.2], sharex = ax) # histx
32 ax_histy = fig.add_axes([left + width + spacing, bottom, 0.2, height], sharey = ax) # histy
33
34 for axes in [ax, ax_histx, ax_histy]:
35 axes.grid(visible = True, which = "both", axis = "both", color = "0.85", linestyle = "-")
36 axes.tick_params(axis = "y", which = "minor")
37 axes.tick_params(axis = "x", which = "minor")
38
39 # no labels
40 ax_histy.tick_params(axis = "y", labelleft = False)
41 ax_histx.tick_params(axis = "x", labelbottom = False)
42
43 # the scatter plot:
44 ax.scatter ( df.values[:, 0]
45 , df.values[:, 1]
46 , s = 8
47 , zorder = 1000
48 )
49
50 ax_histx.hist(df.values[:, 0], bins = 50, zorder = 1000)
51 ax_histy.hist(df.values[:, 1], bins = 50, orientation = "horizontal", zorder = 1000)
52
53 ax.set_xlabel("X", fontsize = 17)
54 ax.set_ylabel("Y", fontsize = 17)
55 ax.legend([file.split(".")[-2]], fontsize = fontsize)
56
57 #plt.tight_layout()
58 plt.savefig(file.replace(".txt",".png"))

Visualization of the example output
Test:
test_pm_sampleAffinity
Todo:
High Priority: The performance of the algorithm for the case dim = 1 can be further improved by looping over the sample in the innermost layer.


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, Thursday 2:45 AM, August 19, 2021, Dallas, TX
Fatemeh Bagheri, Wednesday 00:01 AM, August 25, 2021, Dallas, TX

Definition at line 502 of file pm_sampleAffinity.F90.


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