Page Menu
Home
c4science
Search
Configure Global Search
Log In
Files
F68904919
surface_statistics.hh
No One
Temporary
Actions
Download File
Edit File
Delete File
View Transforms
Subscribe
Mute Notifications
Award Token
Subscribers
None
File Metadata
Details
File Info
Storage
Attached
Created
Sat, Jun 29, 10:21
Size
6 KB
Mime Type
text/x-c++
Expires
Mon, Jul 1, 10:21 (2 d)
Engine
blob
Format
Raw Data
Handle
18519443
Attached To
rTAMAAS tamaas
surface_statistics.hh
View Options
/**
* @file
* @section LICENSE
*
* Copyright (©) 2016-2020 EPFL (École Polytechnique Fédérale de Lausanne),
* Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
*
* This program is free software: you can redistribute it and/or modify
* it under the terms of the GNU Affero General Public License as published
* by the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* This program is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU Affero General Public License for more details.
*
* You should have received a copy of the GNU Affero General Public License
* along with this program. If not, see <https://www.gnu.org/licenses/>.
*
*/
/* -------------------------------------------------------------------------- */
#ifndef SURFACE_STATISTICS_HH
#define SURFACE_STATISTICS_HH
/* -------------------------------------------------------------------------- */
#include "surface.hh"
/* -------------------------------------------------------------------------- */
namespace
tamaas
{
/**
* @deprecated should move code to Statistics class instead
* @brief Statistics calculations on surfaces
*/
class
SurfaceStatistics
{
/* ------------------------------------------------------------------------ */
/* Methods */
/* ------------------------------------------------------------------------ */
public
:
static
Real
computeMaximum
(
const
Surface
<
Real
>&
s
);
static
Real
computeMinimum
(
const
Surface
<
Real
>&
s
);
static
Real
computeAverage
(
const
Surface
<
Real
>&
s
,
Real
epsilon
=
std
::
numeric_limits
<
Real
>::
lowest
())
{
return
computeAverage
<
Real
>
(
s
,
epsilon
);
}
template
<
typename
T
>
static
T
computeAverage
(
const
Surface
<
T
>&
s
,
T
epsilon
=
std
::
numeric_limits
<
T
>::
lowest
());
static
Real
computeStdev
(
const
Grid
<
Real
,
2
>&
s
,
Real
epsilon
=
std
::
numeric_limits
<
Real
>::
lowest
());
static
Real
computeSkewness
(
const
Surface
<
Real
>&
surface
,
Real
avg
,
Real
std
);
static
Real
computeKurtosis
(
const
Surface
<
Real
>&
surface
,
Real
avg
,
Real
std
);
static
Real
computeRMSSlope
(
const
Surface
<
Real
>&
surface
);
static
Real
computeSpectralRMSSlope
(
Surface
<
Real
>&
surface
);
static
Real
computeSpectralMeanCurvature
(
Surface
<
Real
>&
surface
);
static
Real
computeSpectralStdev
(
Surface
<
Real
>&
surface
);
// static void computeStatistics(Surface<Real> & surface);
static
std
::
vector
<
Real
>
computeMoments
(
Surface
<
Real
>&
surface
);
static
Real
computeAnalyticFractalMoment
(
UInt
order
,
UInt
k1
,
UInt
k2
,
Real
hurst
,
Real
rms
,
Real
L
);
//! compute perimeter
static
Real
computePerimeter
(
const
Surface
<
Real
>&
surface
,
Real
eps
=
1e-10
);
//! get real contact area
static
Real
computeContactArea
(
const
Surface
<
Real
>&
surface
);
//! get real contact area
static
Real
computeContactAreaRatio
(
const
Surface
<
Real
>&
surface
);
template
<
bool
consider_zeros
=
false
>
static
std
::
vector
<
int
>
computeDistribution
(
const
Surface
<
Real
>&
surface
,
unsigned
int
n_bins
,
Real
min_value
=
std
::
numeric_limits
<
Real
>::
max
(),
Real
max_value
=
std
::
numeric_limits
<
Real
>::
lowest
());
static
std
::
vector
<
Complex
>
computeSpectralDistribution
(
const
Surface
<
Real
>&
surface
,
unsigned
int
n_bins
,
UInt
cutoff
=
std
::
numeric_limits
<
UInt
>::
max
());
//! compute average amplitude
static
Real
computeAverageAmplitude
(
const
Surface
<
Real
>&
surface
);
static
Complex
computeAverageAmplitude
(
const
SurfaceComplex
<
Real
>&
surface
);
//! sum all heights
static
Real
computeSum
(
const
Surface
<
Real
>&
surface
);
//! compute auto correlation
static
Surface
<
Real
>
computeAutocorrelation
(
Surface
<
Real
>&
surface
);
//! compute powerspectrum
static
SurfaceComplex
<
Real
>
computePowerSpectrum
(
Surface
<
Real
>&
surface
);
};
/* -------------------------------------------------------------------------- */
template
<
bool
consider_zeros
>
std
::
vector
<
int
>
SurfaceStatistics
::
computeDistribution
(
const
Surface
<
Real
>&
surface
,
unsigned
int
n_bins
,
Real
min_value
,
Real
max_value
)
{
if
(
n_bins
==
0
)
SURFACE_FATAL
(
"badly set number of bins"
);
Real
min_heights
;
Real
max_heights
;
std
::
vector
<
int
>
bins
;
bins
.
resize
(
n_bins
);
bins
.
assign
(
n_bins
,
0
);
min_heights
=
min_value
;
max_heights
=
max_value
;
if
(
max_value
==
-
1.
*
std
::
numeric_limits
<
Real
>::
max
())
max_heights
=
computeMaximum
(
surface
);
if
(
min_value
==
std
::
numeric_limits
<
Real
>::
max
())
min_heights
=
computeMinimum
(
surface
);
int
counter
=
0
;
const
int
n
=
surface
.
size
();
for
(
int
i
=
0
;
i
<
n
*
n
;
++
i
)
{
Real
value
=
surface
(
i
);
if
(
!
consider_zeros
&&
value
<
1e-14
)
continue
;
if
(
value
>
max_heights
||
value
<
min_heights
)
SURFACE_FATAL
(
"out of bound value: "
<<
min_heights
<<
" <? "
<<
value
<<
" <? "
<<
max_heights
<<
" : bounds probably badly chosen"
);
Real
ind
=
(
value
-
min_heights
)
/
(
max_heights
-
min_heights
);
unsigned
int
index
=
ind
*
n_bins
;
if
(
index
==
n_bins
)
--
index
;
++
bins
[
index
];
++
counter
;
}
if
(
consider_zeros
&&
counter
!=
n
*
n
)
SURFACE_FATAL
(
"problem in counting points"
);
// Real offset = (max_heights - min_heights)/n_bins;
// Real x = 0;
// Real f1 = 0, f2 = 0;
// for (unsigned int i = 0 ; i < n_bins ; ++i)
// {
// x = i*offset+min_heights;
// f1 = bins[i]/Real(n*n);
// f2 = bins[i]/Real(counter);
//}
Real
min_non_zero
=
std
::
numeric_limits
<
Real
>::
max
();
#pragma omp parallel for reduction(min : min_non_zero)
for
(
int
i
=
0
;
i
<
n
*
n
;
++
i
)
{
Real
value
=
surface
(
i
);
if
(
value
<
min_non_zero
&&
value
>
1e-14
)
min_non_zero
=
value
;
}
return
bins
;
}
/* -------------------------------------------------------------------------- */
template
<
typename
T
>
T
SurfaceStatistics
::
computeAverage
(
const
Surface
<
T
>&
surface
,
T
epsilon
)
{
Real
avg
=
0.0
;
UInt
cpt
=
0
;
UInt
n
=
surface
.
size
();
#pragma omp parallel for reduction(+ : avg, cpt)
for
(
UInt
i
=
0
;
i
<
n
*
n
;
i
++
)
{
Real
val
=
surface
(
i
);
if
(
std
::
abs
(
val
)
>
epsilon
)
{
avg
+=
val
;
++
cpt
;
}
}
avg
/=
cpt
;
return
avg
;
}
}
// namespace tamaas
#endif
/* __SURFACE_STATISTICS_HH__ */
Event Timeline
Log In to Comment