| Title: | Retinal Reconstruction Program |
|---|---|
| Description: | Reconstructs retinae by morphing a flat surface with cuts (a dissected flat-mount retina) onto a curvilinear surface (the standard retinal shape). It can estimate the position of a point on the intact adult retina to within 8 degrees of arc (3.6% of nasotemporal axis). The coordinates in reconstructed retinae can be transformed to visuotopic coordinates. For more details see Sterratt, D. C., Lyngholm, D., Willshaw, D. J. and Thompson, I. D. (2013) <doi:10.1371/journal.pcbi.1002921>. |
| Authors: | David C. Sterratt [aut, cre, cph], Daniel Lyngholm [aut, cph], Jan Okul [aut, cph] |
| Maintainer: | David C. Sterratt <[email protected]> |
| License: | CC BY-NC-SA 4.0 |
| Version: | 0.8.2 |
| Built: | 2026-05-21 09:23:47 UTC |
| Source: | https://github.com/davidcsterratt/retistruct |
An AnnotatedOutline contains a function to annotate tears on the outline.
AnnotatedOutline object, with extra fields for tears, full
cuts, the latitude of rim phi0 and the index of fixed
point i0.
retistruct::OutlineCommon -> retistruct::Outline -> retistruct::PathOutline -> AnnotatedOutline
tearsMatrix in which each row represents a cut by the
indices into the outline points of the apex (V0) and
backward (VB) and forward (VF) points
fullcutsMatrix in which each row represents a full cut by the
indices into the outline points of the forward cut (VF0 and VF1)and
backward cut (VB0 and VB1).
phi0rim angle in radians
lambda0longitude of fixed point
i0index of fixed point
retistruct::OutlineCommon$clearFeatureSets()retistruct::OutlineCommon$getFeatureSet()retistruct::OutlineCommon$getFeatureSetTypes()retistruct::OutlineCommon$getFeatureSets()retistruct::OutlineCommon$getIDs()retistruct::Outline$addFeatureSet()retistruct::Outline$getDepth()retistruct::Outline$getFragment()retistruct::Outline$getFragmentIDs()retistruct::Outline$getFragmentIDsFromPointIDs()retistruct::Outline$getFragmentPointIDs()retistruct::Outline$getFragmentPoints()retistruct::Outline$getImage()retistruct::Outline$getOutlineLengths()retistruct::Outline$getOutlineSet()retistruct::Outline$getPoints()retistruct::Outline$getPointsScaled()retistruct::Outline$getPointsXY()retistruct::Outline$mapFragment()retistruct::Outline$mapPids()retistruct::Outline$replaceImage()retistruct::PathOutline$insertPoint()retistruct::PathOutline$nextPoint()retistruct::PathOutline$stitchSubpaths()new()
Constructor
AnnotatedOutline$new(...)
...Parameters to PathOutline
labelTearPoints()
Label a set of three unlabelled points supposed
to refer to the apex and vertices of a tear with the V0
(Apex), VF (forward vertex) and VB (backward vertex) labels.
AnnotatedOutline$labelTearPoints(pids)
pidsthe vector of three indices
Vector of indices labelled with V0, VF and VB
whichTear()
Return index of tear in an AnnotatedOutline in which a point appears
AnnotatedOutline$whichTear(pid)
pidID of point
ID of tear
getTear()
Return indices of tear in AnnotatedOutline
AnnotatedOutline$getTear(tid)
tidTear ID, which can be returned from whichTear()
Vector of three point IDs, labelled with V0,
VF and VB
getTears()
Get tears
AnnotatedOutline$getTears()
Matrix of tears
computeTearRelationships()
Compute the parent relationships for a potential set of tears. The function throws an error if tears overlap.
AnnotatedOutline$computeTearRelationships(tears = NULL)
tearsMatrix containing columns V0 (Apices of tears)
VB (Backward vertices of tears) and VF (Forward
vertices of tears)
List containing
Rsetthe set of points on the rim
TFsetlist containing indices of points in each forward tear
TBsetlist containing indices of points in each backward tear
hcorrespondence mapping
hfcorrespondence mapping in forward direction for points on boundary
hbcorrespondence mapping in backward direction for points on boundary
addTear()
Add tear to an AnnotatedOutline
AnnotatedOutline$addTear(pids)
pidsVector of three point IDs to be added
removeTear()
Remove tear from an AnnotatedOutline
AnnotatedOutline$removeTear(tid)
tidTear ID, which can be returned from whichTear()
checkTears()
Check that all tears are correct.
AnnotatedOutline$checkTears()
If all is OK, returns empty vector. If not, returns indices of problematic tears.
setFixedPoint()
Set fixed point
AnnotatedOutline$setFixedPoint(i0, name)
i0Index of fixed point
nameName of fixed point
getFixedPoint()
Get point ID of fixed point
AnnotatedOutline$getFixedPoint()
Point ID of fixed point
getRimSet()
Get point IDs of points on rim
AnnotatedOutline$getRimSet()
Point IDs of points on rim. If the outline has been
stitched (see StitchedOutline), the point IDs
will be ordered in the direction of the
forward pointer.
getBoundarySets()
Get point IDs of points on boundaries
AnnotatedOutline$getBoundarySets()
List of Point IDs of points on the boundaries.
If the outline has been stitched (see StitchedOutline),
the point IDs in each
element of the list will be ordered in the direction of the
forward pointer, and the boundary that is longest will be
named as Rim. If the outline has not been stitched,
the list will have one element named Rim.
ensureFixedPointInRim()
Ensure that the fixed point i0 is in the rim, not a tear.
Alters object in which i0 may have been changed.
AnnotatedOutline$ensureFixedPointInRim()
labelFullCutPoints()
Label a set of four unlabelled points supposed to refer to a cut.
AnnotatedOutline$labelFullCutPoints(pids)
pidsthe vector of point indices
addFullCut()
Add cut to an AnnotatedOutline
AnnotatedOutline$addFullCut(pids)
pidsVector of three point IDs to be added
whichFullCut()
Return index of cut in an AnnotatedOutline in which a point appears
AnnotatedOutline$whichFullCut(pid)
pidID of point
ID of cut
removeFullCut()
Remove cut from an AnnotatedOutline
AnnotatedOutline$removeFullCut(cid)
cidFullCut ID, which can be returned from
whichFullCut
computeFullCutRelationships()
Compute the cut relationships between the points
AnnotatedOutline$computeFullCutRelationships(fullcuts)
fullcutsMatrix containing columns VB0,
and VB1 (Backward vertices of fullcuts) and VF0 and VF1 (Forward
vertices of fullcuts)
List containing
Rsetthe set of points on the rim
TFsetlist containing indices of points in each forward cut
TBsetlist containing indices of points in each backward cut
hcorrespondence mapping
hfcorrespondence mapping in forward direction for points on boundary
hbcorrespondence mapping in backward direction for points on boundary
getFullCut()
Return indices of fullcuts in AnnotatedOutline
AnnotatedOutline$getFullCut(cid)
cidFullCut ID, which can be returned from whichFullCut
Vector of four point IDs, labelled with VF1,
VF1, VB0 and VB1
getFullCuts()
Return indices of fullcuts in AnnotatedOutline
AnnotatedOutline$getFullCuts()
Matrix in which each row contains point IDs, for the forward and backward
sides of the cut: VF0, VF1, VB0 and VB1
addPoints()
Add points to the outline register of points
AnnotatedOutline$addPoints(P, fid)
P2 column matrix of points to add
fidfragment id of the points
The ID of each added point in the register. If points already exist a point will not be created in the register, but an ID will be returned
getRimLengths()
Get lengths of edges on rim
AnnotatedOutline$getRimLengths()
Vector of rim lengths
clone()
The objects of this class are cloneable with this method.
AnnotatedOutline$clone(deep = FALSE)
deepWhether to make a deep clone.
David Sterratt
P <- rbind(c(1,1), c(2,1), c(2,-1), c(1,-1), c(1,-2), c(-1,-2), c(-1,-1), c(-2,-1),c(-2,1), c(-1,1), c(-1,2), c(1,2)) o <- TriangulatedOutline$new(P) o$addTear(c(3, 4, 5)) o$addTear(c(6, 7, 8)) o$addTear(c(9, 10, 11)) o$addTear(c(12, 1, 2)) flatplot(o) P <- list(rbind(c(1,1), c(2,1), c(2.5,2), c(3,1), c(4,1), c(1,4)), rbind(c(-1,1), c(-1,4), c(-2,3), c(-2,2), c(-3,2), c(-4,1)), rbind(c(-4,-1), c(-1,-1), c(-1,-4)), rbind(c(1,-1), c(2,-1), c(2.5,-2), c(3,-1), c(4,-1), c(1,-4))) o <- AnnotatedOutline$new(P) o$addTear(c(2, 3, 4)) o$addTear(c(17, 18, 19)) o$addTear(c(9, 10, 11)) o$addFullCut(c(1, 5, 16, 20)) flatplot(o)P <- rbind(c(1,1), c(2,1), c(2,-1), c(1,-1), c(1,-2), c(-1,-2), c(-1,-1), c(-2,-1),c(-2,1), c(-1,1), c(-1,2), c(1,2)) o <- TriangulatedOutline$new(P) o$addTear(c(3, 4, 5)) o$addTear(c(6, 7, 8)) o$addTear(c(9, 10, 11)) o$addTear(c(12, 1, 2)) flatplot(o) P <- list(rbind(c(1,1), c(2,1), c(2.5,2), c(3,1), c(4,1), c(1,4)), rbind(c(-1,1), c(-1,4), c(-2,3), c(-2,2), c(-3,2), c(-4,1)), rbind(c(-4,-1), c(-1,-1), c(-1,-4)), rbind(c(1,-1), c(2,-1), c(2.5,-2), c(3,-1), c(4,-1), c(1,-4))) o <- AnnotatedOutline$new(P) o$addTear(c(2, 3, 4)) o$addTear(c(17, 18, 19)) o$addTear(c(9, 10, 11)) o$addFullCut(c(1, 5, 16, 20)) flatplot(o)
Convert azimuth-elevation coordinates to spherical coordinates
azel.to.sphere.colatitude(r, r0)azel.to.sphere.colatitude(r, r0)
r |
Coordinates of points in azimuth-elevation coordinates
represented as 2 column matrix with column names |
r0 |
Direction of the axis of the sphere on which to project
represented as a 2 column matrix of with column names |
2-column matrix of spherical coordinates of points with
column names psi (colatitude) and lambda (longitude).
David Sterratt
r0 <- cbind(alpha=0, theta=0) r <- rbind(r0, r0+c(1,0), r0-c(1,0), r0+c(0,1), r0-c(0,1)) azel.to.sphere.colatitude(r, r0)r0 <- cbind(alpha=0, theta=0) r <- rbind(r0, r0+c(1,0), r0-c(1,0), r0+c(0,1), r0-c(0,1)) azel.to.sphere.colatitude(r, r0)
Azimuthal conformal or stereographic or Wulff projection
azimuthal.conformal(r, ...)azimuthal.conformal(r, ...)
r |
2-column Matrix of spherical coordinates of points on
sphere. Column names are |
... |
Arguments not used by this projection. |
2-column Matrix of Cartesian coordinates of points on polar
projection. Column names should be x and y.
This is a special case with the point centred on the projection being the South Pole. The MathWorld equations are for the more general case.
David Sterratt
https://en.wikipedia.org/wiki/Map_projection, http://mathworld.wolfram.com/StereographicProjection.html Fisher, N. I., Lewis, T., and Embleton, B. J. J. (1987). Statistical analysis of spherical data. Cambridge University Press, Cambridge, UK.
Lambert azimuthal equal area projection
azimuthal.equalarea(r, ...)azimuthal.equalarea(r, ...)
r |
2-column Matrix of spherical coordinates of points on
sphere. Column names are |
... |
Arguments not used by this projection. |
2-column Matrix of Cartesian coordinates of points on polar
projection. Column names should be x and y.
This is a special case with the point centred on the projection being the South Pole. The MathWorld equations are for the more general case.
David Sterratt
https://en.wikipedia.org/wiki/Map_projection, http://mathworld.wolfram.com/LambertAzimuthalEqual-AreaProjection.html Fisher, N. I., Lewis, T., and Embleton, B. J. J. (1987). Statistical analysis of spherical data. Cambridge University Press, Cambridge, UK.
Azimuthal equidistant projection
azimuthal.equidistant(r, ...)azimuthal.equidistant(r, ...)
r |
2-column Matrix of spherical coordinates of points on
sphere. Column names are |
... |
Arguments not used by this projection. |
2-column Matrix of Cartesian coordinates of points on polar
projection. Column names should be x and y.
This is a special case with the point centred on the projection being the South Pole. The MathWorld equations are for the more general case.
David Sterratt
https://en.wikipedia.org/wiki/Map_projection, http://mathworld.wolfram.com/AzimuthalEquidistantProjection.html
Given a triangular mesh on a sphere described by mesh locations
(phi, lambda), a radius R and a triangulation
Trt, determine the Cartesian coordinates of points cb
given in barycentric coordinates with respect to the mesh.
bary.to.sphere.cart(phi, lambda, R, Trt, cb)bary.to.sphere.cart(phi, lambda, R, Trt, cb)
phi |
Latitudes of mesh points |
lambda |
Longitudes of mesh points |
R |
Radius of sphere |
Trt |
Triangulation |
cb |
Object returned by tsearch containing information on the triangle in which a point occurs and the barycentric coordinates within that triangle |
An N-by-3 matrix of the Cartesian coordinates of the points
David Sterratt
On a sphere the central angle between two points is defined as the
angle whose vertex is the centre of the sphere and that subtends
the arc formed by the great circle between the points. This
function computes the central angle for two points and .
central.angle(phi1, lambda1, phi2, lambda2)central.angle(phi1, lambda1, phi2, lambda2)
phi1 |
Latitude of first point |
lambda1 |
Longitude of first point |
phi2 |
Latitude of second point |
lambda2 |
Longitude of second point |
Central angle
David Sterratt
Wikipedia https://en.wikipedia.org/wiki/Central_angle
Check the whether directory contains valid data
checkDatadir(dir = NULL)checkDatadir(dir = NULL)
dir |
Directory to check. |
TRUE if dir contains valid data;
FALSE otherwise.
David Sterratt
Return points on the unit circle in an anti-clockwise
direction. If L is not specified n points are
returned. If L is specified, the same number of points are
returned as there are elements in L, the interval between
successive points being proportional to L.
circle(n = 12, L = NULL)circle(n = 12, L = NULL)
n |
Number of points |
L |
Intervals between points |
The cartesian coordinates of the points
David Sterratt
Find the intersections of the plane defined by the normal n and the
distance d expressed as a fractional distance along the side of
each triangle.
compute.intersections.sphere(phi, lambda, T, n, d)compute.intersections.sphere(phi, lambda, T, n, d)
phi |
Latitude of grid points on sphere centred on origin. |
lambda |
Longitude of grid points on sphere centred on origin. |
T |
Triangulation |
n |
Normal of plane |
d |
Distance of plane along normal from origin. |
Matrix with same dimensions as T. Each row gives
the intersection of the plane with the corresponding triangle in
T. Column 1 gives the fractional distance from vertex 2 to
vertex 3. Column 2 gives the fractional distance from vertex 3 to
vertex 1. Column 2 gives the fractional distance from vertex 1 to
vertex 2. A value of NaN indicates that the corresponding
edge lies in the plane. A value of Inf indicates that the
edge lies parallel to the plane but outside it.
David Sterratt
Compute a kernel estimate over a grid and do a contour analysis
of this estimate. The contour heights the determined by finding
heights that exclude a certain fraction of the probability. For
example, the 95
and it should enclose about 5
are specified by the contour.levels option; by default
they are c(5, 25, 50, 75, 95).
compute.kernel.estimate(Dss, phi0, fhat, compute.conc)compute.kernel.estimate(Dss, phi0, fhat, compute.conc)
Dss |
List of datasets. The first two columns of each datasets
are coordinates of points on the sphere in spherical polar
(latitude, |
phi0 |
Rim angle in radians |
fhat |
Function such as |
compute.conc |
Function to return the optimal value of the concentration parameter kappa given the data. |
A list containing
kappa |
The concentration parameter |
h |
A pseudo-bandwidth parameter, the inverse of the square root of |
flevels |
Contour levels. |
labels |
Labels of the contours. |
g |
Raw density estimate drawn on non-area-preserving projection. Comprises locations of gridlines in Cartesian coordinates ( |
gpa |
Raw density estimate drawn on area-preserving projection. Comprises same elements as above. |
contour.areas |
Area of each individual contour. One level may have more than one contour; this shows the areas of all such contours. |
tot.contour.areas |
Data frame containing the total area within the contours at each level. |
David Sterratt
FeatureSet to represent counts centred
on pointsA CountSet contains information about points located
on Outlines. Each CountSet contains a list of
matrices, each of which has columns labelled X and
Y describing the cartesian coordinates (in the unscaled
coordinate frame) of the centres of boxes in the Outline, and a
column C representing the counts in those boxes.
retistruct::FeatureSetCommon -> retistruct::FeatureSet -> CountSet
new()
Constructor
CountSet$new(data = NULL, cols = NULL)
dataList of matrices describing data. Each matrix
should have columns named X, Y and C
colsNamed vector of colours for each data set. The name is
used as the ID (label) for the data set. The colours should be names
present in the output of the colors function
reconstruct()
Map the CountSet to a ReconstructedOutline
CountSet$reconstruct(ro)
roclone()
The objects of this class are cloneable with this method.
CountSet$clone(deep = FALSE)
deepWhether to make a deep clone.
David Sterratt
Create grid on projection of hemisphere onto plane
create.polar.cart.grid(pa, res, phi0)create.polar.cart.grid(pa, res, phi0)
pa |
If |
res |
Resolution of grid |
phi0 |
Value of |
List containing:
s |
Grid locations in spherical coordinates |
c |
Grid locations in Cartesian coordinates on plane |
xs |
X grid line locations in Cartesian coordinates on plane |
ys |
Y grid line locations in Cartesian coordinates on plane |
David Sterratt
Read a retinal dataset in CSV format. Each dataset is a folder
containing a file called outline.csv that specifies the outline in
X-Y coordinates. It may also contain a file datapoints.csv,
containing the locations of data points and a file
datacounts.csv, containing the locations of data counts;
see read.datapoints and
read.datacounts for the formats of these files. The
folder may also contain a file od.csv specifying the
coordinates of the optic disc.
csv.read.dataset(dataset, report = message)csv.read.dataset(dataset, report = message)
dataset |
Path to directory containing |
report |
Function to report progress |
A RetinalOutline object
David Sterratt
The function that computes the gradient of the energy (or error)
of the deformation of the mesh from the flat outline to the
sphere. This depends on the locations of the points given in
spherical coordinates. The function is designed to take these as a
vector that is received from the optim function.
dE( p, Cu, C, L, B, Tr, A, R, Rset, i0, phi0, lambda0, Nphi, N, alpha = 1, x0, nu = 1, verbose = FALSE )dE( p, Cu, C, L, B, Tr, A, R, Rset, i0, phi0, lambda0, Nphi, N, alpha = 1, x0, nu = 1, verbose = FALSE )
p |
Parameter vector of |
Cu |
The upper part of the connectivity matrix |
C |
The connectivity matrix |
L |
Length of each edge in the flattened outline |
B |
Connectivity matrix |
Tr |
Triangulation in the flattened outline |
A |
Area of each triangle in the flattened outline |
R |
Radius of the sphere |
Rset |
Indices of points on the rim |
i0 |
Index of fixed point on rim |
phi0 |
Latitude at which sphere curtailed |
lambda0 |
Longitude of fixed points |
Nphi |
Number of free values of |
N |
Number of points in sphere |
alpha |
Area penalty scaling coefficient |
x0 |
Area penalty cut-off coefficient |
nu |
Power to which to raise area |
verbose |
How much information to report |
A vector representing the derivative of the energy of this particular configuration with respect to the parameter vector
David Sterratt
Draw the "flat" outline in 3D with depth information
depthplot3D(r, ...)depthplot3D(r, ...)
r |
|
... |
Parameters depending on class of |
David Sterratt
File system directories used by shinyFiles
directories()directories()
The function that computes the energy (or error) of the
deformation of the mesh from the flat outline to the sphere. This
depends on the locations of the points given in spherical
coordinates. The function is designed to take these as a vector
that is received from the optim function.
E( p, Cu, C, L, B, Tr, A, R, Rset, i0, phi0, lambda0, Nphi, N, alpha = 1, x0, nu = 1, verbose = FALSE )E( p, Cu, C, L, B, Tr, A, R, Rset, i0, phi0, lambda0, Nphi, N, alpha = 1, x0, nu = 1, verbose = FALSE )
p |
Parameter vector of |
Cu |
The upper part of the connectivity matrix |
C |
The connectivity matrix |
L |
Length of each edge in the flattened outline |
B |
Connectivity matrix |
Tr |
Triangulation in the flattened outline |
A |
Area of each triangle in the flattened outline |
R |
Radius of the sphere |
Rset |
Indices of points on the rim |
i0 |
Index of fixed point on rim |
phi0 |
Latitude at which sphere curtailed |
lambda0 |
Longitude of fixed points |
Nphi |
Number of free values of |
N |
Number of points in sphere |
alpha |
Area scaling coefficient |
x0 |
Area cut-off coefficient |
nu |
Power to which to raise area |
verbose |
How much information to report |
A single value, representing the energy of this particular configuration
David Sterratt
The function that computes the energy (or error) of the
deformation of the mesh from the flat outline to the sphere. This
depends on the locations of the points given in spherical
coordinates. The function is designed to take these as a vector
that is received from the optim function.
Ecart(P, Cu, L, Tr, A, R, alpha = 1, x0, nu = 1, verbose = FALSE)Ecart(P, Cu, L, Tr, A, R, alpha = 1, x0, nu = 1, verbose = FALSE)
P |
N-by-3 matrix of point coordinates |
Cu |
The upper part of the connectivity matrix |
L |
Length of each edge in the flattened outline |
Tr |
Triangulation in the flattened outline |
A |
Area of each triangle in the flattened outline |
R |
Radius of sphere |
alpha |
Area penalty scaling coefficient |
x0 |
Area penalty cut-off coefficient |
nu |
Power to which to raise area |
verbose |
How much information to report |
A single value, representing the energy of this particular configuration
David Sterratt
Piecewise, smooth function that increases linearly with negative arguments.
f(x, x0)f(x, x0)
x |
Main argument |
x0 |
The cut-off parameter. Above this value the function is zero. |
The value of the function.
David Sterratt
The function that computes the gradient of the energy (or error)
of the deformation of the mesh from the flat outline to the
sphere. This depends on the locations of the points given in
spherical coordinates. The function is designed to take these as a
vector that is received from the optim function.
Fcart(P, C, L, Tr, A, R, alpha = 1, x0, nu = 1, verbose = FALSE)Fcart(P, C, L, Tr, A, R, alpha = 1, x0, nu = 1, verbose = FALSE)
P |
N-by-3 matrix of point coordinates |
C |
The connectivity matrix |
L |
Length of each edge in the flattened outline |
Tr |
Triangulation in the flattened outline |
A |
Area of each triangle in the flattened outline |
R |
Radius of sphere |
alpha |
Area penalty scaling coefficient |
x0 |
Area penalty cut-off coefficient |
nu |
Power to which to raise area |
verbose |
How much information to report |
A vector representing the derivative of the energy of this particular configuration with respect to the parameter vector
David Sterratt
OutlinesA FeatureSet contains information about features
located on Outlines. Each FeatureSet contains a
list of matrices, each of which has columns labelled X
and Y describing the cartesian coordinates of points on
the Outline, in the unscaled coordinate frame. Derived classes,
e.g. a CountSet, may have extra columns. Each matrix
in the list has an associated label and colour, which is used by
plotting functions.
retistruct::FeatureSetCommon -> FeatureSet
new()
Constructor
FeatureSet$new(data = NULL, cols = NULL, type = NULL)
dataList of matrices describing data. Each matrix
should have columns named X and Y
colsNamed vector of colours for each data set. The name is
used as the ID (label) for the data set. The colours should be names
present in the output of the colors function
typeString
clone()
The objects of this class are cloneable with this method.
FeatureSet$clone(deep = FALSE)
deepWhether to make a deep clone.
David Sterratt
FeatureSets and
ReconstructedFeatureSetsAn FeatureSetCommon has functionality for retrieving sets of features (e.g. points or landmarks associated with an outline)
dataList of matrices describing data
colsVector of colours for each data set
typeString giving type of feature set
getIndex()
Get numeric index of features
FeatureSetCommon$getIndex(fid)
fidFeature ID (string)
getIDs()
Get IDs of features
FeatureSetCommon$getIDs()
Vector of IDs of features
setID()
Set name
FeatureSetCommon$setID(i, fid)
iNumeric index of feature
fidFeature ID (string)
getFeature()
Get feature by feature ID
FeatureSetCommon$getFeature(fid)
fidFeature ID string
Matrix describing feature
getFeatures()
Get all features
FeatureSetCommon$getFeatures()
getCol()
Get colour in which to plot feature ID
FeatureSetCommon$getCol(fid)
fidFeature ID string
clone()
The objects of this class are cloneable with this method.
FeatureSetCommon$clone(deep = FALSE)
deepWhether to make a deep clone.
David Sterratt
This is an implementation of the FIRE algorithm for structural relaxation put forward by Bitzek et al. (2006)
fire( r, force, restraint, m = 1, dt = 0.1, maxmove = 100, dtmax = 1, Nmin = 5, finc = 1.1, fdec = 0.5, astart = 0.1, fa = 0.99, a = 0.1, nstep = 100, tol = 1e-05, verbose = FALSE, report = message )fire( r, force, restraint, m = 1, dt = 0.1, maxmove = 100, dtmax = 1, Nmin = 5, finc = 1.1, fdec = 0.5, astart = 0.1, fa = 0.99, a = 0.1, nstep = 100, tol = 1e-05, verbose = FALSE, report = message )
r |
Initial locations of particles |
force |
Force function |
restraint |
Restraint function |
m |
Masses of points |
dt |
Initial time step |
maxmove |
Maximum distance to move in any time step |
dtmax |
Maximum time step |
Nmin |
Number of steps after which to start increasing |
finc |
Fractional increase in |
fdec |
Fractional decrease in |
astart |
Starting value of |
fa |
Fraction of |
a |
Initial value of |
nstep |
Maximum number of steps |
tol |
Tolerance - if RMS force is below this value, stop and report convergence |
verbose |
If |
report |
Function to report progress when |
List containing x, the positions of the points,
conv, which is 0 if convergence as occurred and 1 otherwise,
and frms, the root mean square of the forces on the
particles.
David Sterratt
Bitzek, E., Koskinen, P., Gähler, F., Moseler, M., and Gumbsch, P. (2006). Structural relaxation made simple. Phys. Rev. Lett., 97:170201.
Plot "flat" (unreconstructed) representation of outline
flatplot(x, ...)flatplot(x, ...)
x |
|
... |
Other plotting parameters |
David Sterratt
Plot flat AnnotatedOutline. The user markup is
displayed by default.
## S3 method for class 'AnnotatedOutline' flatplot(x, axt = "n", xlim = NULL, ylim = NULL, markup = TRUE, ...)## S3 method for class 'AnnotatedOutline' flatplot(x, axt = "n", xlim = NULL, ylim = NULL, markup = TRUE, ...)
x |
|
axt |
whether to plot axes |
xlim |
x-limits |
ylim |
y-limits |
markup |
If |
... |
Other plotting parameters |
David Sterratt
Plot flat Outline.
## S3 method for class 'Outline' flatplot( x, axt = "n", xlim = NULL, ylim = NULL, add = FALSE, image = TRUE, scalebar = 1, rimset = FALSE, pids = FALSE, pid.joggle = 0, lwd.outline = 1, ... )## S3 method for class 'Outline' flatplot( x, axt = "n", xlim = NULL, ylim = NULL, add = FALSE, image = TRUE, scalebar = 1, rimset = FALSE, pids = FALSE, pid.joggle = 0, lwd.outline = 1, ... )
x |
|
axt |
whether to plot axes |
xlim |
x limits |
ylim |
y limits |
add |
If |
image |
If |
scalebar |
If numeric and if the Outline has a |
rimset |
If |
pids |
If |
pid.joggle |
Amount to joggle point IDs by randomly |
lwd.outline |
Line width of outline |
... |
Other plotting parameters |
David Sterratt
Plot ReconstructedOutline object. This adds a mesh
of gridlines from the spherical retina (described by points
phi, lambda and triangulation Trt and cut-off
point phi0) onto a flattened retina (described by points
P and triangulation T).
## S3 method for class 'ReconstructedOutline' flatplot( x, axt = "n", xlim = NULL, ylim = NULL, grid = TRUE, strain = FALSE, ... )## S3 method for class 'ReconstructedOutline' flatplot( x, axt = "n", xlim = NULL, ylim = NULL, grid = TRUE, strain = FALSE, ... )
x |
|
axt |
whether to plot axes |
xlim |
x-limits |
ylim |
y-limits |
grid |
Whether or not to show the grid lines of latitude and longitude |
strain |
Whether or not to show the strain |
... |
Other plotting parameters |
David Sterratt
Plot flat StitchedOutline. If the optional argument
stitch is TRUE the user markup is displayed.
## S3 method for class 'StitchedOutline' flatplot(x, axt = "n", xlim = NULL, ylim = NULL, stitch = TRUE, lwd = 1, ...)## S3 method for class 'StitchedOutline' flatplot(x, axt = "n", xlim = NULL, ylim = NULL, stitch = TRUE, lwd = 1, ...)
x |
|
axt |
whether to plot axes |
xlim |
x-limits |
ylim |
y-limits |
stitch |
If |
lwd |
Line width |
... |
Other parameters |
David Sterratt
TriangulatedOutline.Plot flat TriangulatedOutline.
## S3 method for class 'TriangulatedOutline' flatplot(x, axt = "n", xlim = NULL, ylim = NULL, mesh = TRUE, ...)## S3 method for class 'TriangulatedOutline' flatplot(x, axt = "n", xlim = NULL, ylim = NULL, mesh = TRUE, ...)
x |
|
axt |
whether to plot axes |
xlim |
x-limits |
ylim |
y-limits |
mesh |
If |
... |
Other plotting parameters |
David Sterratt
In the projection of points onto the sphere, some triangles maybe flipped, i.e. in the wrong orientation. This functions determines which triangles are flipped by computing the vector pointing to the centre of each triangle and comparing this direction to vector product of two sides of the triangle.
flipped.triangles(Ps, Trt, R = 1)flipped.triangles(Ps, Trt, R = 1)
Ps |
N-by-2 matrix with columns containing latitudes
( |
Trt |
Triangulation of points |
R |
Radius of sphere |
List containing:
flipped |
Indices of in rows of |
cents |
Vectors of centres. |
areas |
Areas of triangles. |
David Sterratt
In the projection of points onto the sphere, some triangles maybe flipped, i.e. in the wrong orientation. This function determines which triangles are flipped by computing the vector pointing to the centre of each triangle and comparing this direction to vector product of two sides of the triangle.
flipped.triangles.cart(P, Trt, R)flipped.triangles.cart(P, Trt, R)
P |
Points in Cartesian coordinates |
Trt |
Triangulation of points |
R |
Radius of sphere |
List containing:
flipped |
Indices of in rows of |
cents |
Vectors of centres. |
areas |
Areas of triangles. |
David Sterratt
Derivative of f
fp(x, x0)fp(x, x0)
x |
Main argument |
x0 |
The cut-off parameter. Above this value the function is zero. |
The value of the function.
David Sterratt
P, as described below.Construct an outline object. This sanitises the input points
P, as described below.
Construct an outline object. This sanitises the input points
P, as described below.
PA N-by-2 matrix of points of the Outline
arranged in anticlockwise order
gfFor each row of P, the index of P that
is next in the outline travelling anticlockwise (forwards)
gbFor each row of P, the index of P that
is next in the outline travelling clockwise (backwards)
hFor each row of P, the cut of that
point (which will be to itself initially)
A.totTotal area of the Fragment
initializeFromPoints()
Initialise a Fragment from a set of points
Fragment$initializeFromPoints(P)
PAn N-by-2 matrix of points of the Outline
clone()
The objects of this class are cloneable with this method.
Fragment$clone(deep = FALSE)
deepWhether to make a deep clone.
David Sterratt
The identity transformation
identity.transform(r, ...)identity.transform(r, ...)
r |
Coordinates of points in spherical coordinates
represented as 2 column matrix with column names |
... |
Other arguments |
Identical matrix
David Sterratt
Read one of the Thompson lab's retinal datasets. Each dataset is a folder containing a SYS file in SYSTAT format and a MAP file in text format. The SYS file specifies the locations of the data points and the MAP file specifies the outline.
idt.read.dataset(dataset, report = message, d.close = 0.25)idt.read.dataset(dataset, report = message, d.close = 0.25)
dataset |
Path to directory containing as SYS and MAP file |
report |
Function to report progress |
d.close |
Maximum distance between points for them to count as the same point. This is expressed as a fraction of the width of the outline. |
The function returns the outline of the retina. In order to do so,
it has to join up the segments of the MAP file. The tracings are
not always precise; sometimes there are gaps between points that
are actually the same point. The parameter d.close specifies
how close points must be to count as the same point.
dataset |
The path to the directory given as an argument |
raw |
List containing
|
P |
The points of the outline |
gf |
Forward pointers along the outline |
gb |
Backward pointers along the outline |
Ds |
List of datapoints |
Ss |
List of landmark lines |
David Sterratt
Read a retinal dataset in IJROI format. Each dataset is a folder
containing a file called outline.roi that specifies the
outline in X-Y coordinates. It may also contain a file
datapoints.csv, containing the locations of data points;
see read.datapoints for the format of this file. The
folder may also contain a file od.roi specifying the
coordinates of the optic disc.
ijroi.read.dataset(dataset, report = report)ijroi.read.dataset(dataset, report = report)
dataset |
Path to directory containing |
report |
Function to report progress |
A RetinalOutline object
David Sterratt
Read a retinal dataset in IJROI format. Each dataset is a folder
containing a file called outline.roi that specifies the
outline in X-Y coordinates. It may also contain a file
datapoints.csv, containing the locations of data points;
see read.datapoints for the format of this file. The
folder may also contain a file od.roi specifying the
coordinates of the optic disc.
ijroimulti.read.dataset(dataset)ijroimulti.read.dataset(dataset)
dataset |
Path to directory containing |
A RetinalOutline object
David Sterratt
Interpolate values in image
interpolate.image(im, P, invert.y = FALSE, wmin = 10, wmax = 100)interpolate.image(im, P, invert.y = FALSE, wmin = 10, wmax = 100)
im |
image to interpolate |
P |
N by 2 matrix of x, y values at which to interpolate. x
is in range |
invert.y |
If |
wmin |
minimum window size for inferring NA values |
wmax |
maximum window size for inferring NA values |
Vector of N interpolated values
David Sterratt
Invert sphere about its centre
invert.sphere(r, ...)invert.sphere(r, ...)
r |
Coordinates of points in spherical coordinates
represented as 2 column matrix with column names |
... |
Other arguments |
Matrix in same format, but with pi added to lambda
and phi negated.
David Sterratt
Invert image of a partial sphere and scale the longitude so that
points at latitude phi0 is projected onto a longitude of 0
degrees (the equator).
invert.sphere.to.hemisphere(r, phi0, ...)invert.sphere.to.hemisphere(r, phi0, ...)
r |
Coordinates of points in spherical coordinates
represented as 2 column matrix with column names |
phi0 |
The latitude to map onto the equator |
... |
Other arguments |
Matrix in same format, but with pi added to lambda
and phi negated and scaled so that the longitude
phi0 is projected to 0 degrees (the equator)
David Sterratt
The Karcher mean of a set of points on a manifold is defined as the point whose sum of squared Riemann distances to the points is minimal. On a sphere using spherical coordinates this distance can be computed using the formula for central angle.
karcher.mean.sphere(x, na.rm = FALSE, var = FALSE)karcher.mean.sphere(x, na.rm = FALSE, var = FALSE)
x |
Matrix of points on sphere as N-by-2 matrix with labelled
columns |
na.rm |
logical value indicating whether |
var |
logical value indicating whether variance should be returned too. |
Vector of means with components named phi and
lambda. If var is TRUE, a list containing
mean and variance in elements mean and var.
David Sterratt
Heo, G. and Small, C. G. (2006). Form representations and means for landmarks: A survey and comparative study. Computer Vision and Image Understanding, 102:188-203.
Find the optimal concentration for a set of data
kde.compute.concentration(mu)kde.compute.concentration(mu)
mu |
Data in spherical coordinates |
The optimal concentration
David Sterratt
Kernel density estimate on sphere using Fisherian density with polar coordinates
kde.fhat(r, mu, kappa)kde.fhat(r, mu, kappa)
r |
Locations at which to estimate density in polar coordinates |
mu |
Locations of data points in polar coordinates |
kappa |
Concentration parameter |
Vector of density estimates
David Sterratt
Kernel density estimate on sphere using Fisherian density with Cartesian coordinates
kde.fhat.cart(r, mu, kappa)kde.fhat.cart(r, mu, kappa)
r |
Locations at which to estimate density in Cartesian coordinates on unit sphere |
mu |
Locations of data points in Cartesian coordinates on unit sphere |
kappa |
Concentration parameter |
Vector of density estimates
David Sterratt
Estimate of the log likelihood of the points mu given a particular value of the concentration kappa
kde.L(mu, kappa)kde.L(mu, kappa)
mu |
Locations of data points in Cartesian coordinates on unit sphere |
kappa |
Concentration parameter |
Log likelihood of data
David Sterratt
Find the optimal concentration for a set of data
kr.compute.concentration(mu, y)kr.compute.concentration(mu, y)
mu |
Locations in Cartesian coordinates (independent variables) |
y |
Values at locations (dependent variables) |
The optimal concentration
David Sterratt
Cross validation estimate of the least squares error of the points mu given a particular value of the concentration kappa
kr.sscv(mu, y, kappa)kr.sscv(mu, y, kappa)
mu |
Locations in Cartesian coordinates (independent variables) |
y |
Values at locations (dependent variables) |
kappa |
Concentration parameter |
Least squares error
David Sterratt
Kernel regression on sphere using Fisherian density with polar coordinates
kr.yhat(r, mu, y, kappa)kr.yhat(r, mu, y, kappa)
r |
Locations at which to estimate dependent variables in polar coordinates |
mu |
Locations in polar coordinates (independent variables) |
y |
Values at data points (dependent variables) |
kappa |
Concentration parameter |
Estimates of dependent variables at locations r
David Sterratt
Kernel regression on sphere using Fisherian density with Cartesian coordinates
kr.yhat.cart(r, mu, y, kappa)kr.yhat.cart(r, mu, y, kappa)
r |
Locations at which to estimate dependent variables in Cartesian coordinates |
mu |
Locations in Cartesian coordinates (independent variables) |
y |
Values at locations (dependent variables) |
kappa |
Concentration parameter |
Estimates of dependent variables at locations r
David Sterratt
FeatureSet to represent pointsA LandmarkSet contains information about points
located on Outlines. Each LandmarkSet contains a
list of matrices, each of which has columns labelled X
and Y describing the cartesian coordinates (in the
unscaled coordinate frame) of points in landmarks on the
Outline.
retistruct::FeatureSetCommon -> retistruct::FeatureSet -> LandmarkSet
new()
Constructor
LandmarkSet$new(data = NULL, cols = NULL)
dataList of matrices describing data. Each matrix
should have columns named X and Y
colsNamed vector of colours for each data set. The name is
used as the ID (label) for the data set. The colours should be names
present in the output of the colors function
reconstruct()
Map the LandmarkSet to a ReconstructedOutline
LandmarkSet$reconstruct(ro)
roclone()
The objects of this class are cloneable with this method.
LandmarkSet$clone(deep = FALSE)
deepWhether to make a deep clone.
David Sterratt
Determine the intersection of two lines L1 and L2 in two dimensions, using the formula described by Weisstein.
line.line.intersection(P1, P2, P3, P4, interior.only = FALSE)line.line.intersection(P1, P2, P3, P4, interior.only = FALSE)
P1 |
vector containing x,y coordinates of one end of L1 |
P2 |
vector containing x,y coordinates of other end of L1 |
P3 |
vector containing x,y coordinates of one end of L2 |
P4 |
vector containing x,y coordinates of other end of L2 |
interior.only |
boolean flag indicating whether only intersections inside L1 and L2 should be returned. |
Vector containing x,y coordinates of intersection of L1
and L2. If L1 and L2 are parallel, this is infinite-valued. If
interior.only is TRUE, then when the intersection
does not occur between P1 and P2 and P3 and P4, a vector
containing NAs is returned.
David Sterratt
Weisstein, Eric W. "Line-Line Intersection." From MathWorld–A Wolfram Web Resource. http://mathworld.wolfram.com/Line-LineIntersection.html
## Intersection of two intersecting lines line.line.intersection(c(0, 0), c(1, 1), c(0, 1), c(1, 0)) ## Two lines that don't intersect line.line.intersection(c(0, 0), c(0, 1), c(1, 0), c(1, 1))## Intersection of two intersecting lines line.line.intersection(c(0, 0), c(1, 1), c(0, 1), c(1, 0)) ## Two lines that don't intersect line.line.intersection(c(0, 0), c(0, 1), c(1, 0), c(1, 1))
Convert an list created by R6_to_list() into an R6 object.
list_to_R6(l)list_to_R6(l)
l |
list created by R6_to_list() |
R6 object or list list
David Sterratt
List valid datasets underneath a directory. This reports all directories that appear to be valid.
list.datasets(path = ".", verbose = FALSE)list.datasets(path = ".", verbose = FALSE)
path |
Directory path to start searching from |
verbose |
If |
A vector of directories containing datasets
David Sterratt
Plot the fractional change in length of mesh edges. The length of each edge in the mesh in the reconstructed object is plotted against each edge in the spherical object. The points are colour-coded according to the amount of log strain each edge is under.
lvsLplot(r, ...)lvsLplot(r, ...)
r |
|
... |
Other plotting parameters |
David Sterratt
Morph a flat dataset to a parabola for testing purposes
morph.dataset.to.parabola( orig.dataset = file.path(system.file(package = "retistruct"), "extdata", "smi32-csv"), morphed.dataset = NA, f = 100 )morph.dataset.to.parabola( orig.dataset = file.path(system.file(package = "retistruct"), "extdata", "smi32-csv"), morphed.dataset = NA, f = 100 )
orig.dataset |
Directory of original dataset |
morphed.dataset |
Directory to write morphed dataset to. If NA, a temporary directory will be created |
f |
Focus of parabola |
Path to morphed.dataset
David Sterratt
Return a new version of the list in which any unnamed elements have been given standardised names
name.list(l)name.list(l)
l |
the list with unnamed elements |
The list with standardised names
David Sterratt
Bring angle into range
normalise.angle(theta)normalise.angle(theta)
theta |
Angle to bring into range |
Normalised angle
David Sterratt
Orthographic projection
orthographic(r, proj.centre = cbind(phi = 0, lambda = 0), ...)orthographic(r, proj.centre = cbind(phi = 0, lambda = 0), ...)
r |
Latitude-longitude coordinates in a matrix with columns
labelled |
proj.centre |
Location of centre of projection as matrix with
column names |
... |
Arguments not used by this projection.n |
Two-column matrix with columns labelled x and
y of locations of projection of coordinates on plane
David Sterratt
https://en.wikipedia.org/wiki/Map_projection, http://mathworld.wolfram.com/OrthographicProjection.html
An Outline has contains the polygon describing the outline and an image associated with the outline.
retistruct::OutlineCommon -> Outline
PA N-by-2 matrix of points of the Outline
arranged in anticlockwise order
scaleThe length of one unit of P in the X-Y plane in arbitrary units
scalezThe length of one unit of P in the Z-direction in arbitrary units
unitsString giving units of scaled P, e.g. “um”
gfFor each row of P, the index of P that
is next in the outline travelling anticlockwise (forwards)
gbFor each row of P, the index of P that
is next in the outline travelling clockwise (backwards)
hFor each row of P, the cut of that
point (which will be to itself initially)
imAn image as a raster object
dmDepthmap, with same dimensions as im, which indicates
height of each pixel in Z-direction
A.fragmentsAreas of fragments
dm.inferna.window.minMinimum window size (in pixels) for inferring missing values in depthmaps
dm.inferna.window.maxMinimum window size (in pixels) for inferring missing values in depthmaps
new()
Construct an outline object. This sanitises the
input points P.
Outline$new( fragments = list(), scale = NA, im = NULL, scalez = NA, dm = NULL, units = NA )
fragmentsA list of N-by-2 matrix of points for each fragment of the Outline
scaleThe length of one unit of P in arbitrary units
imThe image as a raster object
scalezThe length of one unit of P in the Z-direction in arbitrary units. If NA, the depthmap is ignored.
dmDepthmap, with same dimensions as im, which indicates
height of each pixel in Z-direction
unitsString giving units of scaled P, e.g. “um”
getImage()
Image accessor
Outline$getImage()
An image as a raster object
replaceImage()
Image setter
Outline$replaceImage(im)
imAn image as a raster object
mapFragment()
Map the point IDs of a Fragment on the point IDs of this Outline
Outline$mapFragment(fragment, pids)
mapPids()
Map references to points
Outline$mapPids(x, y, pids)
xReferences to point indices in source
yReferences to existing point indices in target
pidsIDs of points in point register
New references to point indices in target
getDepth()
Get depth of points P
Outline$getDepth(P)
Pmatrix containing unscaled X-Y coordinates of points
Vector of unscaled Z coordinates of points P
addPoints()
Add points to the outline register of points
Outline$addPoints(P, fid)
P2 or 3 column matrix of points to add
fidID of fragment to which to add the points
The ID of each added point in the register. If points already exist a point will not be created in the register, but an ID will be returned
getFragmentPointIDs()
Get the point IDs in a fragment
Outline$getFragmentPointIDs(fid)
fidfragment id of the points
Vector of point IDs, i.e. indices of the rows in
the matrices returned by getPoints and
getPointsScaled
getFragmentPoints()
Get the points in a fragment
Outline$getFragmentPoints(fid)
fidfragment id of the points
Vector of points
getFragment()
Get fragment
Outline$getFragment(fid)
fidFragment ID
The Fragment object with ID fid
getFragmentIDsFromPointIDs()
Get fragment IDs from point IDS
Outline$getFragmentIDsFromPointIDs(pids)
pidsVector of point IDs
The Fragment ID to which each point belongs
getFragmentIDs()
Get fragment IDs
Outline$getFragmentIDs()
IDs of all fragments
getPoints()
Get unscaled mesh points
Outline$getPoints(pids = NULL)
pidsIDs of point to return
Matrix with columns X, Y and Z
getPointsXY()
Get X-Y coordinates of unscaled mesh points
Outline$getPointsXY(pids = NULL)
pidsIDs of point to return
Matrix with columns X and Y
getPointsScaled()
Get scaled mesh points
Outline$getPointsScaled()
Matrix with columns X and Y which is
exactly scale times the matrix returned by getPoints
getRimSet()
Get set of points on rim
Outline$getRimSet()
Vector of point IDs, i.e. indices of the rows in
the matrices returned by getPoints and
getPointsScaled
getOutlineSet()
Get points on the edge of the outline
Outline$getOutlineSet()
Vector of points IDs on outline
getOutlineLengths()
Get lengths of edges of the outline
Outline$getOutlineLengths()
Vector of lengths of edges connecting neighbouring points
addFeatureSet()
Add a FeatureSet, e.g. a PointSet or LandmarkSet
Outline$addFeatureSet(fs)
fsFeatureSet to add
clone()
The objects of this class are cloneable with this method.
Outline$clone(deep = FALSE)
deepWhether to make a deep clone.
David Sterratt
An OutlineCommon has functionality for retrieving sets of features (e.g. points or landmarks associated with an outline)
versionVersion of reconstruction file data format
featureSetsList of feature sets associated with the outline, which may be of various types, e.g. a PointSet or LandmarkSet
getFeatureSets()
Get all the feature sets
OutlineCommon$getFeatureSets()
List of FeatureSets associated with the outline
getFeatureSet()
Get all feature sets of a particular type, e.g. PointSet or LandmarkSet
OutlineCommon$getFeatureSet(type)
typeThe type of the feature set as a string
All FeatureSets of that type
clearFeatureSets()
Clear all feature sets from the outline
OutlineCommon$clearFeatureSets()
getIDs()
Get all the distinct IDs contained in the FeatureSets
OutlineCommon$getIDs()
Vector of IDs
getFeatureSetTypes()
Get all the distinct types of FeatureSets
OutlineCommon$getFeatureSetTypes()
Vector of types as strings, e.g. PointSet, LandmarkSet
clone()
The objects of this class are cloneable with this method.
OutlineCommon$clone(deep = FALSE)
deepWhether to make a deep clone.
Ancillary function to place labels
panlabel(panlabel, line = -0.7)panlabel(panlabel, line = -0.7)
panlabel |
Label text |
line |
Line on which to appear |
David Sterratt
Arc length of a parabola y=x^2/4f
parabola.arclength(x1, x2, f)parabola.arclength(x1, x2, f)
x1 |
x co-ordinate of start of arc |
x2 |
x co-ordinate of end of arc |
f |
focal length of parabola |
length of parabola arc
David Sterratt
Inverse arc length of a parabola y=x^2/4f
parabola.invarclength(x1, s, f)parabola.invarclength(x1, s, f)
x1 |
co-ordinate of start of arc |
s |
length of parabola arc to follow |
f |
focal length of parabola |
x co-ordinate of end of arc
David Sterratt
Parse dependencies
parse.dependencies(deps)parse.dependencies(deps)
deps |
Text produced by, e.g., |
Table with package column, relationship column and version number
David Sterratt
Add point fullcuts to the outline
Add point fullcuts to the outline
The member function stitchSubpaths() stitches together two
subpaths of the outline. One subpath is stitched in the forward
direction from the point indexed by VF0 to the point
indexed by VF1. The other is stitched in the backward
direction from VB0 to VB1. Each point in the subpath
is linked to points in the opposing pathway at an equal or
near-equal fraction along. If a point exists in the opposing
pathway within a distance epsilon of the projection, this
point is connected. If no point exists within this tolerance, a
new point is created.
To the Outline object this adds
hf |
point cut mapping in forward direction for points on boundary |
hb |
point cut mapping in backward direction for points on boundary |
retistruct::OutlineCommon -> retistruct::Outline -> PathOutline
hfForward fullcuts
hbBackward fullcuts
retistruct::OutlineCommon$clearFeatureSets()retistruct::OutlineCommon$getFeatureSet()retistruct::OutlineCommon$getFeatureSetTypes()retistruct::OutlineCommon$getFeatureSets()retistruct::OutlineCommon$getIDs()retistruct::Outline$addFeatureSet()retistruct::Outline$getDepth()retistruct::Outline$getFragment()retistruct::Outline$getFragmentIDs()retistruct::Outline$getFragmentIDsFromPointIDs()retistruct::Outline$getFragmentPointIDs()retistruct::Outline$getFragmentPoints()retistruct::Outline$getImage()retistruct::Outline$getOutlineLengths()retistruct::Outline$getOutlineSet()retistruct::Outline$getPoints()retistruct::Outline$getPointsScaled()retistruct::Outline$getPointsXY()retistruct::Outline$getRimSet()retistruct::Outline$initialize()retistruct::Outline$mapFragment()retistruct::Outline$mapPids()retistruct::Outline$replaceImage()addPoints()
Add points to the outline register of points
PathOutline$addPoints(P, fid)
P2 column matrix of points to add
fidfragment id of the points
The ID of each added point in the register. If points already exist a point will not be created in the register, but an ID will be returned
nextPoint()
Get next point in path for
PathOutline$nextPoint(pids)
pidsPoint IDs of points to get next position
insertPoint()
Insert point at a fractional distance between points
PathOutline$insertPoint(i0, i1, f)
i0Point ID of first point
i1Point ID of second point
fFraction of distance between points i0 and
i1 at which to insert point
stitchSubpaths()
Stitch subpaths
PathOutline$stitchSubpaths(VF0, VF1, VB0, VB1, epsilon)
VF0First vertex of “forward” subpath
VF1Second vertex of “forward” subpath
VB0First vertex of “backward” subpath
VB1Second vertex of “backward” subpath
epsilonMinimum distance between points
clone()
The objects of this class are cloneable with this method.
PathOutline$clone(deep = FALSE)
deepWhether to make a deep clone.
FeatureSet to represent pointsA PointSet contains information about points located
on Outlines. Each PointSet contains a list of
matrices, each of which has columns labelled X and
Y describing the cartesian coordinates (in the unscaled
coordinate frame) of points on the Outline.
retistruct::FeatureSetCommon -> retistruct::FeatureSet -> PointSet
new()
Constructor
PointSet$new(data = NULL, cols = NULL)
dataList of matrices describing data. Each matrix
should have columns named X and Y
colsNamed vector of colours for each data set. The name is
used as the ID (label) for the data set. The colours should be names
present in the output of the colors function
reconstruct()
Map the PointSet to a ReconstructedOutline
PointSet$reconstruct(ro)
roclone()
The objects of this class are cloneable with this method.
PointSet$clone(deep = FALSE)
deepWhether to make a deep clone.
David Sterratt
This is the inverse of sphere.spherical.to.polar.cart
polar.cart.to.sphere.spherical(r, pa = FALSE, preserve = "latitude")polar.cart.to.sphere.spherical(r, pa = FALSE, preserve = "latitude")
r |
2-column Matrix of Cartesian coordinates of points on
polar projection. Column names should be |
pa |
If |
preserve |
Quantity to preserve locally in the
projection. Options are |
2-column Matrix of spherical coordinates of points on
sphere. Column names are phi and lambda.
David Sterratt
Place text at bottom right of projection
polartext(text)polartext(text)
text |
Test to place |
David Sterratt
Plot projection of a reconstructed outline
projection(r, ...)projection(r, ...)
r |
Object such as a |
... |
Other plotting parameters |
David Sterratt
Draw a projection of a ReconstructedOutline. This method sets up
the grid lines and the angular labels and draws the image.
## S3 method for class 'ReconstructedOutline' projection( r, transform = identity.transform, axisdir = cbind(phi = 90, lambda = 0), projection = azimuthal.equalarea, proj.centre = cbind(phi = 0, lambda = 0), lambdalim = c(-180, 180), philim = c(-90, 90), labels = c(0, 90, 180, 270), mesh = FALSE, grid = TRUE, grid.bg = "transparent", grid.int.minor = 15, grid.int.major = 45, colatitude = TRUE, pole = FALSE, image = TRUE, markup = TRUE, add = FALSE, max.proj.dim = getOption("max.proj.dim"), ... )## S3 method for class 'ReconstructedOutline' projection( r, transform = identity.transform, axisdir = cbind(phi = 90, lambda = 0), projection = azimuthal.equalarea, proj.centre = cbind(phi = 0, lambda = 0), lambdalim = c(-180, 180), philim = c(-90, 90), labels = c(0, 90, 180, 270), mesh = FALSE, grid = TRUE, grid.bg = "transparent", grid.int.minor = 15, grid.int.major = 45, colatitude = TRUE, pole = FALSE, image = TRUE, markup = TRUE, add = FALSE, max.proj.dim = getOption("max.proj.dim"), ... )
r |
|
transform |
Transform function to apply to spherical coordinates before rotation |
axisdir |
Direction of axis (North pole) of sphere in
external space as matrix with column names |
projection |
Projection in which to display object,
e.g. |
proj.centre |
Location of centre of projection as matrix with
column names |
lambdalim |
Limits of longitude (in degrees) to display |
philim |
Limits of latitude (in degrees) to display |
labels |
Vector of 4 labels to plot at 0, 90, 180 and 270 degrees |
mesh |
If |
grid |
Whether or not to show the grid lines of latitude and longitude |
grid.bg |
Background colour of the grid |
grid.int.minor |
Interval between minor grid lines in degrees |
grid.int.major |
Interval between major grid lines in degrees |
colatitude |
If |
pole |
If |
image |
If |
markup |
If |
add |
If |
max.proj.dim |
Maximum width of the image created in pixels |
... |
Graphical parameters to pass to plotting functions |
Plot projection of reconstructed dataset
## S3 method for class 'RetinalReconstructedOutline' projection( r, transform = identity.transform, projection = azimuthal.equalarea, axisdir = cbind(phi = 90, lambda = 0), proj.centre = cbind(phi = 0, lambda = 0), lambdalim = c(-180, 180), datapoints = TRUE, datapoint.means = TRUE, datapoint.contours = FALSE, grouped = FALSE, grouped.contours = FALSE, landmarks = TRUE, mesh = FALSE, grid = TRUE, image = TRUE, ids = r$getIDs(), ... )## S3 method for class 'RetinalReconstructedOutline' projection( r, transform = identity.transform, projection = azimuthal.equalarea, axisdir = cbind(phi = 90, lambda = 0), proj.centre = cbind(phi = 0, lambda = 0), lambdalim = c(-180, 180), datapoints = TRUE, datapoint.means = TRUE, datapoint.contours = FALSE, grouped = FALSE, grouped.contours = FALSE, landmarks = TRUE, mesh = FALSE, grid = TRUE, image = TRUE, ids = r$getIDs(), ... )
r |
|
transform |
Transform function to apply to spherical coordinates before rotation |
projection |
Projection in which to display object,
e.g. |
axisdir |
Direction of axis (North pole) of sphere in external space |
proj.centre |
Location of centre of projection as matrix with
column names |
lambdalim |
Limits of longitude (in degrees) to display |
datapoints |
If |
datapoint.means |
If |
datapoint.contours |
If |
grouped |
If |
grouped.contours |
If |
landmarks |
If |
mesh |
If |
grid |
If |
image |
If |
ids |
IDs of groups of data within a dataset, returned using
|
... |
Graphical parameters to pass to plotting functions |
Convert an R6 object into a list, ignoring functions and environments
R6_to_list(r, path = "", envs = list())R6_to_list(r, path = "", envs = list())
r |
R6 object or list |
path |
root of the path to the list - no need to supply. Not used but could be developed for pretty-printing |
envs |
list of environments already encountered - do not set |
List with structure mirroring the R6 object.
David Sterratt
Restore points to spherical manifold after an update of the Lagrange integration rule
Rcart(P, R, Rset, i0, phi0, lambda0)Rcart(P, R, Rset, i0, phi0, lambda0)
P |
Point positions as N-by-3 matrix |
R |
Radius of sphere |
Rset |
Indices of points on rim |
i0 |
Index of fixed point |
phi0 |
FullCut-off of curtailed sphere in radians |
lambda0 |
Longitude of fixed point on rim |
Points projected back onto sphere
David Sterratt
Read data counts from a file ‘datacounts.csv’ in the
directory dataset. The CSV file should contain two columns
for every dataset. Each pair of columns must contain a unique name
in the first cell of the first row and a valid colour in the
second cell of the first row. In the remaining rows, the X
coordinates of data counts should be in the first column and the Y
coordinates should be in the second column.
read.datacounts(dataset)read.datacounts(dataset)
dataset |
Path to directory containing |
List containing
Ds |
List of sets of data counts. Each set comprises a 2-column matrix and each set is named. |
cols |
List of colours for each dataset. There is one element that corresponds to each element of |
David Sterratt
Read data points from a file dataponts.csv in the directory
dataset. The CSV should contain two columns for every
dataset. Each pair of columns must contain a unique name in the
first cell of the first row and a valid colour in the second
cell of the first row. In the remaining rows, the X coordinates of
data points should be in the first column and the Y coordinates
should be in the second column.
read.datapoints(dataset)read.datapoints(dataset)
dataset |
Path to directory containing |
List containing
Ds |
List of sets of datapoints. Each set comprises a 2-column matrix and each set is named. |
cols |
List of colours for each dataset. There is one element that corresponds to each element of |
David Sterratt
A ReconstructedCountSet contains information about
features located on ReconstructedOutlines. Each
ReconstructedCountSet contains a list of matrices, each of which
has columns labelled phi (latitude) and lambda
(longitude) describing the spherical coordinates of points on
the ReconstructedOutline, and a column C representing the
counts at these points.
retistruct::FeatureSetCommon -> retistruct::ReconstructedFeatureSet -> ReconstructedCountSet
KRKernel regression
new()
Constructor
ReconstructedCountSet$new(fs = NULL, ro = NULL)
fsFeatureSet to reconstruct
roReconstructedOutline to which feature
set should be mapped
getKR()
Get kernel regression estimate of grouped data points
ReconstructedCountSet$getKR()
Kernel regression computed using
compute.kernel.estimate
clone()
The objects of this class are cloneable with this method.
ReconstructedCountSet$clone(deep = FALSE)
deepWhether to make a deep clone.
David Sterratt
A ReconstructedFeatureSet contains information about
features located on ReconstructedOutlines. Each
ReconstructedFeatureSet contains a list of matrices, each of
which has columns labelled phi (latitude) and
lambda (longitude) describing the spherical coordinates
of points on the ReconstructedOutline. Derived classes, e.g. a
ReconstructedCountSet, may have extra columns.
Each matrix in the list has an associated label and colour,
which is used by plotting functions.
retistruct::FeatureSetCommon -> ReconstructedFeatureSet
new()
Constructor
ReconstructedFeatureSet$new(fs = NULL, ro = NULL)
fsFeatureSet to reconstruct
roReconstructedOutline to which feature
set should be mapped
clone()
The objects of this class are cloneable with this method.
ReconstructedFeatureSet$clone(deep = FALSE)
deepWhether to make a deep clone.
David Sterratt
A ReconstructedLandmarkSet contains information about
features located on ReconstructedOutlines. Each
ReconstructedLandmarkSet contains a list of matrices, each of
which has columns labelled phi (latitude) and
lambda (longitude) describing the spherical coordinates
of points on the ReconstructedOutline.
retistruct::FeatureSetCommon -> retistruct::ReconstructedFeatureSet -> ReconstructedLandmarkSet
clone()
The objects of this class are cloneable with this method.
ReconstructedLandmarkSet$clone(deep = FALSE)
deepWhether to make a deep clone.
David Sterratt
The function reconstruct reconstructs outline
into spherical surface Reconstruct outline into spherical
surface.
retistruct::OutlineCommon -> ReconstructedOutline
olAnnotated outline
ol0Original Annotated outline
PtTransformed cartesian mesh points
TrtTransformed triangulation
CtTransformed links
CutTransformed links
BtTransformed binary vector representation of edge indices onto a binary vector representation of the indices of the points linked by the edge
LtTransformed lengths
htTransformed correspondences
uIndices of unique points in untransformed space
UTransformed indices of unique points in untransformed space
RsettTransformed rim set
i0tTransformed marker
Hmapping from edges onto corresponding edges
HtTransformed mapping from edges onto corresponding edges
phi0Rim angle
RRadius of spherical template
lambda0Longitude of pole on rim
lambdaLongitudes of transformed mesh points
phiLatitudes of transformed mesh points
PsLocation of mesh point on sphere in spherical coordinates
nNumber of mesh points
alphaWeighting of areas in energy function
x0Area cut-off coefficient
nflip0Initial number flipped triangles
nflipFinal number flipped triangles
optOptimisation object
E.totEnergy function including area
E.lEnergy function based on lengths alone
mean.strainMean strain
mean.logstrainMean log strain
titrationTitrated data structure, saved by titrate
debugDebug function
loadOutline()
Load AnnotatedOutline into ReconstructedOutline object
ReconstructedOutline$loadOutline( ol, n = 500, alpha = 8, x0 = 0.5, plot.3d = FALSE, dev.flat = NA, dev.polar = NA, report = retistruct::report, debug = FALSE )
olAnnotatedOutline object, containing the following information
nNumber of points in triangulation.
alphaArea scaling coefficient
x0Area cut-off coefficient
plot.3dWhether to show 3D picture during optimisation.
dev.flatDevice to plot grid onto. Value of NA (default)
means no plotting.
dev.polarDevice display projection. Value of NA (default) means no plotting.
reportFunction to report progress.
debugIf TRUE print extra debugging output
reconstruct()
Reconstruct Reconstruction proceeds in a number of stages:
The flat object is triangulated with at least n
triangles. This can introduce new vertices in the rim.
The triangulated object is stitched.
The stitched object is triangulated again, but this time it is not permitted to add extra vertices to the rim.
The corresponding points determined by the stitching process are merged to form a new set of merged points and a new triangulation.
The merged points are projected roughly to a sphere.
The locations of the points on the sphere are moved so as to minimise the energy function.
ReconstructedOutline$reconstruct(
plot.3d = FALSE,
dev.flat = NA,
dev.polar = NA,
shinyOutput = NA,
report = getOption("retistruct.report")
)plot.3dIf TRUE make a 3D plot in an RGL window
dev.flatDevice handle for plotting flatplot updates to. If
NA don't make any flat plots
dev.polarDevice handle for plotting polar plot updates
to. If NA don't make any polar plots.
shinyOutputA Shiny output element used to render and display a
plot in the application. If NA or NULL don't output to Shiny.
reportFunction to report progress.
Controlargument to pass to optim
mergePointsEdges()
Merge stitched points and edges.
Create merged and transformed versions (all suffixed with t)
of a number of existing variables, as well
as a matrix Bt, which maps a binary vector representation
of edge indices onto a binary vector representation of the
indices of the points linked by the edge.
Sets following fields
PtTransformed point locations
TrtTransformed triangulation
CtTransformed connection set
CutTransformed symmetric connection set
BtTransformed binary vector representation of edge indices onto a binary vector representation of the indices of the points linked by the edge
LtTransformed edge lengths
htTransformed correspondences
uIndices of unique points in untransformed space
UTransformed indices of unique points in untransformed space
RsetThe set of points on the rim (which has been reordered)
RsettTransformed set of points on rim
i0tTransformed index of the landmark
Hmapping from edges onto corresponding edges
HtTransformed mapping from edges onto corresponding edges
ReconstructedOutline$mergePointsEdges()
projectToSphere()
Project mesh points in the flat outline onto a sphere
This takes the mesh points from the flat outline and maps them to
the curtailed sphere. It uses the area of the flat outline and
phi0 to determine the radius R of the sphere. It
tries to get a good first approximation by using the function
stretchMesh.
The following fields are set:
phiLatitude of mesh points.
lmabdaLongitude of mesh points.
RRadius of sphere.
ReconstructedOutline$projectToSphere()
getStrains()
Return strains edges are under in spherical retina Set information about how edges on the sphere have been deformed from their flat state.
ReconstructedOutline$getStrains()
A list containing two data frames flat and spherical.
Each data frame contains for each edge in the flat or spherical meshes:
LLength of the edge in the flat outline
lLength of the corresponding edge on the sphere
strainThe strain of each connection
logstrainThe logarithmic strain of each connection
optimiseMapping()
Optimise the mapping from the flat outline to the sphere
ReconstructedOutline$optimiseMapping( alpha = 4, x0 = 0.5, nu = 1, optim.method = "BFGS", plot.3d = FALSE, dev.flat = NA, dev.polar = NA, shinyOutput = NULL, control = list() )
alphaArea penalty scaling coefficient
x0Area penalty cut-off coefficient
nuPower to which to raise area
optim.methodMethod to pass to optim
plot.3dIf TRUE make a 3D plot in an RGL window
dev.flatDevice handle for plotting flatplot updates to. If
dev.polarDevice handle for plotting polar plot updates
to. If NA don't make any polar plots.
shinyOutputA Shiny output element used to render and display a
plot in the application.
NA don't make any flat plots
controlControl argument to pass to optim
optimiseMappingCart()
Optimise the mapping from the flat outline to the sphere
ReconstructedOutline$optimiseMappingCart( alpha = 4, x0 = 0.5, nu = 1, method = "BFGS", plot.3d = FALSE, dev.flat = NA, dev.polar = NA, shinyOutput = NULL, ... )
alphaArea penalty scaling coefficient
x0Area penalty cut-off coefficient
nuPower to which to raise area
methodMethod to pass to optim
plot.3dIf TRUE make a 3D plot in an RGL window
dev.flatDevice handle for plotting grid to
dev.polarDevice handle for plotting polar plot to
shinyOutputA Shiny output element used to render and display a plot in the application.
...Extra arguments to pass to fire
transformImage()
Transform an image into the reconstructed space
Transform an image into the reconstructed space. The four corner
coordinates of each pixel are transformed into spherical
coordinates and a mask matrix with the same dimensions as
im is created. This has TRUE for pixels that should
be displayed and FALSE for ones that should not.
Sets the field
imsCoordinates of corners of pixels in spherical coordinates
ReconstructedOutline$transformImage()
getIms()
Get coordinates of corners of pixels of image in spherical coordinates
ReconstructedOutline$getIms()
Coordinates of corners of pixels in spherical coordinates
getTearCoords()
Get locations of tears in spherical coordinates
ReconstructedOutline$getTearCoords()
List containing locations of tears in spherical coordinates
getFullCutCoords()
Get locations of fullcuts in spherical coordinates
ReconstructedOutline$getFullCutCoords()
List containing locations of fullcuts in spherical coordinates
getNonRimBoundaryCoords()
Get location of non-rim boundaries in spherical coordinates
ReconstructedOutline$getNonRimBoundaryCoords()
List containing locations of non-rim boundaries in spherical coordinates
getFeatureSet()
ReconstructedOutline$getFeatureSet(type)
typeBase type of FeatureSet as string.
E.g. PointSet returns a ReconstructedPointSet
reconstructFeatureSets()
Reconstruct any attached feature sets.
ReconstructedOutline$reconstructFeatureSets()
getPoints()
Get mesh points in spherical coordinates
ReconstructedOutline$getPoints()
Matrix with columns phi (latitude) and lambda
(longitude)
mapFlatToSpherical()
Return location of point on sphere corresponding to point on the flat outline
ReconstructedOutline$mapFlatToSpherical(P, dropna = TRUE)
PCartesian coordinates on flat outline as a matrix
with X and Y columns
dropnaWhether to drop points outwith the outline
titrate()
Try a range of values of phi0s in the reconstruction, recording the energy of the mapping in each case.
ReconstructedOutline$titrate( alpha = 8, x0 = 0.5, byd = 1, len.up = 5, len.down = 20 )
alphaArea penalty scaling coefficient
x0Area cut-off coefficient
bydIncrements in degrees
len.upHow many increments to go up from starting value of
phi0 in r.
len.downHow many increments to go up from starting value
of phi0 in r.
clone()
The objects of this class are cloneable with this method.
ReconstructedOutline$clone(deep = FALSE)
deepWhether to make a deep clone.
David Sterratt
A ReconstructedPointSet contains information about
features located on ReconstructedOutlines. Each
ReconstructedPointSet contains a list of matrices, each of
which has columns labelled phi (latitude) and
lambda (longitude) describing the spherical coordinates
of points on the ReconstructedOutline.
retistruct::FeatureSetCommon -> retistruct::ReconstructedFeatureSet -> ReconstructedPointSet
KDEKernel density estimate, computed using
compute.kernel.estimate in getKDE
getMean()
Get Karcher mean of datapoints in spherical coordinates
ReconstructedPointSet$getMean()
Karcher mean of datapoints in spherical coordinates
getHullarea()
Get area of convex hull around data points on sphere
ReconstructedPointSet$getHullarea()
Area in degrees squared
getKDE()
Get kernel density estimate of data points
ReconstructedPointSet$getKDE()
clone()
The objects of this class are cloneable with this method.
ReconstructedPointSet$clone(deep = FALSE)
deepWhether to make a deep clone.
David Sterratt
This is similar to unique(), but spares rows which are duplicated, but at different points in the matrix
remove.identical.consecutive.rows(P)remove.identical.consecutive.rows(P)
P |
Source matrix |
Matrix with identical consecutive rows removed.
David Sterratt
Suppose segments AB and CD intersect. Point B is replaced by the intersection point, defined B'. Point C is replaced by a point C' on the line B'D. The maximum distance of B'C' is given by the parameter d. If the distance l B'D is less than 2d, the distance B'C' is l/2.
remove.intersections(P, d = 50)remove.intersections(P, d = 50)
P |
The points, as a 2-column matrix |
d |
Criterion for maximum distance when points are inserted |
A new closed path without intersections
David Sterratt
Calls function specified by option retistruct.report
report(x, ...)report(x, ...)
x |
First arguments to reporting function |
... |
Arguments to reporting function |
David Sterratt
In addition to fields inherited from
StitchedOutline, a RetinalOutline contains a
dataset field, describing the system path to dataset
directory and metadata specific to retinae and some formats of
retinae
An retinalOutline object. This contains the following fields:
retistruct::OutlineCommon -> retistruct::Outline -> retistruct::PathOutline -> retistruct::AnnotatedOutline -> retistruct::TriangulatedOutline -> retistruct::StitchedOutline -> RetinalOutline
DVflipTRUE if the raw data is flipped in
the dorsoventral direction
sideThe side of the eye (“Left” or “Right”)
datasetFile system path to dataset directory
retistruct::OutlineCommon$clearFeatureSets()retistruct::OutlineCommon$getFeatureSet()retistruct::OutlineCommon$getFeatureSetTypes()retistruct::OutlineCommon$getFeatureSets()retistruct::OutlineCommon$getIDs()retistruct::Outline$addFeatureSet()retistruct::Outline$getDepth()retistruct::Outline$getFragment()retistruct::Outline$getFragmentIDs()retistruct::Outline$getFragmentIDsFromPointIDs()retistruct::Outline$getFragmentPointIDs()retistruct::Outline$getFragmentPoints()retistruct::Outline$getImage()retistruct::Outline$getOutlineLengths()retistruct::Outline$getOutlineSet()retistruct::Outline$getPoints()retistruct::Outline$getPointsScaled()retistruct::Outline$getPointsXY()retistruct::Outline$mapFragment()retistruct::Outline$mapPids()retistruct::Outline$replaceImage()retistruct::PathOutline$insertPoint()retistruct::PathOutline$nextPoint()retistruct::PathOutline$stitchSubpaths()retistruct::AnnotatedOutline$addFullCut()retistruct::AnnotatedOutline$addPoints()retistruct::AnnotatedOutline$addTear()retistruct::AnnotatedOutline$checkTears()retistruct::AnnotatedOutline$computeFullCutRelationships()retistruct::AnnotatedOutline$computeTearRelationships()retistruct::AnnotatedOutline$ensureFixedPointInRim()retistruct::AnnotatedOutline$getFixedPoint()retistruct::AnnotatedOutline$getFullCut()retistruct::AnnotatedOutline$getFullCuts()retistruct::AnnotatedOutline$getRimLengths()retistruct::AnnotatedOutline$getRimSet()retistruct::AnnotatedOutline$getTear()retistruct::AnnotatedOutline$getTears()retistruct::AnnotatedOutline$labelFullCutPoints()retistruct::AnnotatedOutline$labelTearPoints()retistruct::AnnotatedOutline$removeFullCut()retistruct::AnnotatedOutline$removeTear()retistruct::AnnotatedOutline$setFixedPoint()retistruct::AnnotatedOutline$whichFullCut()retistruct::AnnotatedOutline$whichTear()retistruct::TriangulatedOutline$mapTriangulatedFragment()retistruct::TriangulatedOutline$triangulate()retistruct::StitchedOutline$getBoundarySets()retistruct::StitchedOutline$isStitched()retistruct::StitchedOutline$stitchFullCuts()retistruct::StitchedOutline$stitchTears()new()
Constructor
RetinalOutline$new(..., dataset = NULL)
...Parameters to superclass constructors
datasetFile system path to dataset directory
clone()
The objects of this class are cloneable with this method.
RetinalOutline$clone(deep = FALSE)
deepWhether to make a deep clone.
David Sterratt
A RetinalReconstructedOutline overrides methods of
ReconstructedOutline so that they return data point and
landmark coordinates that have been transformed according to the
values of DVflip and side. When reconstructing, it
also computes the “Optic disc displacement”, i.e. the
number of degrees subtended between the optic disc and the pole.
retistruct::OutlineCommon -> retistruct::ReconstructedOutline -> RetinalReconstructedOutline
EODOptic disc displacement in degrees
retistruct::OutlineCommon$clearFeatureSets()retistruct::OutlineCommon$getFeatureSetTypes()retistruct::OutlineCommon$getFeatureSets()retistruct::OutlineCommon$getIDs()retistruct::ReconstructedOutline$getFullCutCoords()retistruct::ReconstructedOutline$getNonRimBoundaryCoords()retistruct::ReconstructedOutline$getPoints()retistruct::ReconstructedOutline$getStrains()retistruct::ReconstructedOutline$loadOutline()retistruct::ReconstructedOutline$mapFlatToSpherical()retistruct::ReconstructedOutline$mergePointsEdges()retistruct::ReconstructedOutline$optimiseMapping()retistruct::ReconstructedOutline$optimiseMappingCart()retistruct::ReconstructedOutline$projectToSphere()retistruct::ReconstructedOutline$reconstructFeatureSets()retistruct::ReconstructedOutline$titrate()retistruct::ReconstructedOutline$transformImage()getIms()
Get coordinates of corners of pixels of image in spherical
coordinates, transformed according to the value of DVflip
RetinalReconstructedOutline$getIms()
Coordinates of corners of pixels in spherical coordinates
getTearCoords()
Get location of tear coordinates in spherical coordinates,
transformed according to the value of DVflip
RetinalReconstructedOutline$getTearCoords()
Location of tear coordinates in spherical coordinates
reconstruct()
RetinalReconstructedOutline$reconstruct(...)
...Parameters to ReconstructedOutline
getFeatureSet()
Get ReconstructedFeatureSet, transformed
according to the value of DVflip
RetinalReconstructedOutline$getFeatureSet(type)
typeBase type of FeatureSet as string.
E.g. PointSet returns a ReconstructedPointSet
clone()
The objects of this class are cloneable with this method.
RetinalReconstructedOutline$clone(deep = FALSE)
deepWhether to make a deep clone.
David Sterratt
Start the Retistruct GUI
retistruct()retistruct()
Object with getData() method to return
reconstructed retina data and environment this which
contains variables in object.
This function reconstructs a number of datasets, using the R
parallel package to distribute the reconstruction of
multiple datasets across CPUs. If datasets is not specified
the function recurses through a directory tree starting at
tldir, determining whether the directory contains valid raw
data and markup, and performing the reconstruction if it does.
retistruct.batch( tldir = ".", outputdir = tldir, datasets = NULL, device = "pdf", titrate = FALSE, cpu.time.limit = 3600, mc.cores = getOption("mc.cores", 2L) )retistruct.batch( tldir = ".", outputdir = tldir, datasets = NULL, device = "pdf", titrate = FALSE, cpu.time.limit = 3600, mc.cores = getOption("mc.cores", 2L) )
tldir |
If datasets is not specified, the top level of the directory tree through which to recurse in order to find datasets. |
outputdir |
directory in which to dump a log file and images |
datasets |
Vector of dataset directories to reconstruct |
device |
string indicating what type of graphics output
required. Options are |
titrate |
Whether to "titrate" the reconstruction for
different values of |
cpu.time.limit |
amount of CPU after which to terminate the process |
mc.cores |
The number of cores to use. Defaults to the value
given by the option |
David Sterratt
Recurse through a directory tree, determining whether the directory contains valid derived data and converting ‘r.rData’ files to files in MATLAB format named ‘r.mat’
retistruct.batch.export.matlab(tldir = ".")retistruct.batch.export.matlab(tldir = ".")
tldir |
The top level of the directory tree through which to recurse |
David Sterratt
Recurse through a directory tree, determining whether the directory contains valid derived data and plotting graphs if it does.
retistruct.batch.figures(tldir = ".", outputdir = tldir, ...)retistruct.batch.figures(tldir = ".", outputdir = tldir, ...)
tldir |
The top level directory of the tree through which to recurse. |
outputdir |
Directory in which to dump a log file and images |
... |
Parameters passed to plotting functions |
David Sterratt
Get titrations from a directory of reconstructions
retistruct.batch.get.titrations(tldir = ".")retistruct.batch.get.titrations(tldir = ".")
tldir |
The top level directory of the tree through which to
recurse. The files have to have been reconstructed with the
|
Plot titrations
retistruct.batch.plot.titrations(tdat)retistruct.batch.plot.titrations(tdat)
tdat |
Output of |
Recurse through a directory tree, determining whether the directory contains valid derived data and extracting summary data if it does.
retistruct.batch.summary(tldir = ".", cache = TRUE)retistruct.batch.summary(tldir = ".", cache = TRUE)
tldir |
The top level directory of the tree through which to recurse. |
cache |
If |
Data frame containing summary data
David Sterratt
Check that markup such as tears and the nasal or dorsal points are present.
retistruct.check.markup(o)retistruct.check.markup(o)
o |
Outline object |
If all markup is present, return TRUE. Otherwise
return FALSE.
David Sterratt
This calls retistruct.cli.process with a time limit
specified by cpu.time.limit.
retistruct.cli( dataset, cpu.time.limit = Inf, outputdir = NA, device = "pdf", ... )retistruct.cli( dataset, cpu.time.limit = Inf, outputdir = NA, device = "pdf", ... )
dataset |
Path to dataset to process |
cpu.time.limit |
Time limit in seconds |
outputdir |
Directory in which to save any figures |
device |
String representing device to print figures to |
... |
Other arguments to pass to |
A list comprising
status |
0 for success, 1 for reaching
|
time |
The time take in seconds |
mess |
Any error message |
David Sterratt
Print a figure to file
retistruct.cli.figure( dataset, outputdir, device = "pdf", width = 6, height = 6, res = 100 )retistruct.cli.figure( dataset, outputdir, device = "pdf", width = 6, height = 6, res = 100 )
dataset |
Path to dataset to process |
outputdir |
Directory in which to save any figures |
device |
String representing device to print figures to |
width |
Width of figures in inches |
height |
Height of figures in inches |
res |
Resolution of figures in dpi (only applies to bitmap devices) |
David Sterratt
This function processes a dataset, saving the
reconstruction data and MATLAB export data to the dataset
directory and printing figures to outputdir.
retistruct.cli.process( dataset, outputdir = NA, device = "pdf", titrate = FALSE, matlab = TRUE )retistruct.cli.process( dataset, outputdir = NA, device = "pdf", titrate = FALSE, matlab = TRUE )
dataset |
Path to dataset to process |
outputdir |
Directory in which to save any figures |
device |
String representing device to print figures to |
titrate |
Whether to titrate or not |
matlab |
Whether to save to MATLAB or not |
David Sterratt
Save as a MATLAB object certain fields of an object r of
classRetinalReconstructedOutline to a file called
r.mat in the directory r$dataset.
retistruct.export.matlab(r, filename = NULL)retistruct.export.matlab(r, filename = NULL)
r |
|
filename |
Filename of output file. If not specified, is
|
David Sterratt
Read a retinal dataset in one of three formats; for information on
formats see idt.read.dataset,
csv.read.dataset and
ijroi.read.dataset. The format is autodetected from
the files in the directory.
retistruct.read.dataset(dataset, report = message, ...)retistruct.read.dataset(dataset, report = message, ...)
dataset |
Path to directory containing the files corresponding to each format. |
report |
Function to report progress. Set to |
... |
Parameters passed to the format-specific functions. |
A RetinalOutline object
David Sterratt
Read the markup data contained in the files ‘markup.csv’,
‘P.csv’ and ‘T.csv’ in the directory ‘dataset’,
which is specified in the reconstruction object r.
retistruct.read.markup(a, error = stop)retistruct.read.markup(a, error = stop)
a |
Dataset object, containing |
error |
Function to run on error, by default |
The tear information is contained in the files ‘P.csv’ and
‘T.csv’. The first file contains the locations of outline
points that the tears were marked up on. The second file contains
the indices of the apex and backward and forward vertices of each
tear. It is necessary to have the file of points just in case the
algorithm that determines P in
retistruct.read.dataset has changed since the markup
of the tears.
The remaining information is contained in the file ‘markup.csv’.
If DVflip is specified, the locations of points P
flipped in the -direction. This operation also requires the
swapping of gf and gb and VF and VB.
o RetinalDataset object
V0 |
Indices in |
VB |
Indices in |
VF |
Indices in |
iN |
Index in |
iD |
Index in |
iOD |
Index in |
phi0 |
Angle of rim in degrees |
DVflip |
Boolean variable indicating if dorsoventral (DV) axis has been flipped |
David Sterratt
Given an outline object with a dataset field, read the
reconstruction data from the file ‘dataset/r.Rdata’.
retistruct.read.recdata(o, check = TRUE)retistruct.read.recdata(o, check = TRUE)
o |
Outline object containing |
check |
If |
If the reconstruction data exists, return a reconstruction
object, else return the outline object o.
David Sterratt
Reconstruct a retina
retistruct.reconstruct( a, report = NULL, plot.3d = FALSE, dev.flat = NA, dev.polar = NA, shinyOutput = NULL, debug = FALSE, ... )retistruct.reconstruct( a, report = NULL, plot.3d = FALSE, dev.flat = NA, dev.polar = NA, shinyOutput = NULL, debug = FALSE, ... )
a |
|
report |
Function to report progress. Set to |
plot.3d |
If |
dev.flat |
The ID of the device to which to plot the flat representation |
dev.polar |
The ID of the device to which to plot the polar representation |
shinyOutput |
A Shiny output element used to render and display a plot in the application. |
debug |
If |
... |
Parameters to be passed to
|
A RetinalReconstructedOutline object
David Sterratt
Save the markup in the RetinalOutline a to a
file called markup.csv in the directory a$dataset.
retistruct.save.markup(a)retistruct.save.markup(a)
a |
|
David Sterratt
Save the reconstruction data in an object r of class
RetinalReconstructedOutline to a file called
r.Rdata in the directory r$dataset.
retistruct.save.recdata(r)retistruct.save.recdata(r)
r |
|
David Sterratt
This rotates points on sphere by specifying the direction its polar axis, i.e. the axis going through (90, 0), should point after (a) a rotation about an axis through the points (0, 0) and (0, 180) and (b) rotation about the original polar axis.
rotate.axis(r, r0)rotate.axis(r, r0)
r |
Coordinates of points in spherical coordinates
represented as 2 column matrix with column names |
r0 |
Direction of the polar axis of the sphere on which to project
represented as a 2 column matrix of with column names |
2-column matrix of spherical coordinates of points with
column names phi (latitude) and lambda (longitude).
David Sterratt
r0 <- cbind(phi=0, lambda=-pi/2) r <- rbind(r0, r0+c(1,0), r0-c(1,0), r0+c(0,1), r0-c(0,1)) r <- cbind(phi=pi/2, lambda=0) rotate.axis(r, r0)r0 <- cbind(phi=0, lambda=-pi/2) r <- rbind(r0, r0+c(1,0), r0-c(1,0), r0+c(0,1), r0-c(0,1)) r <- cbind(phi=pi/2, lambda=0) rotate.axis(r, r0)
The R shiny server responsible for storing a state for each session, handling inputs from the UI to the server, and plotting outputs to the UI. The arguments are all handled by the shiny package and this function should not be instantiated manually.
server(input, output, session)server(input, output, session)
input |
object that holds the UI state (Managed automatically by shiny) |
output |
sends new outputs to the UI (Managed automatically by shiny) |
session |
controls each open instance (Managed automatically by shiny) |
Jan Okul
Simplify a fragment object by removing vertices bordering short edges while not encroaching on any of the outline. At present, this is done by finding concave vertices. It is safe to remove these, at the expense of increasing the area a bit.
simplifyFragment(P, min.frac.length = 0.001, plot = FALSE)simplifyFragment(P, min.frac.length = 0.001, plot = FALSE)
P |
points to simplify |
min.frac.length |
the minimum length as a fraction of the total length of the outline. |
plot |
whether to display plotting or not during simplification |
Simplified outline object
David Sterratt
Simplify a outline object by removing vertices bordering short edges while not encroaching on any of the outline. At present, this is done by finding concave vertices. It is safe to remove these, at the expense of increasing the area a bit.
simplifyOutline(P, min.frac.length = 0.001, plot = FALSE)simplifyOutline(P, min.frac.length = 0.001, plot = FALSE)
P |
points to simplify |
min.frac.length |
the minimum length as a fraction of the total length of the outline. |
plot |
whether to display plotting or not during simplification |
Simplified outline object
David Sterratt
Sinusoidal projection
sinusoidal( r, proj.centre = cbind(phi = 0, lambda = 0), lambdalim = NULL, lines = FALSE, ... )sinusoidal( r, proj.centre = cbind(phi = 0, lambda = 0), lambdalim = NULL, lines = FALSE, ... )
r |
Latitude-longitude coordinates in a matrix with columns
labelled |
proj.centre |
Location of centre of projection as matrix with
column names |
lambdalim |
Limits of longitude to plot |
lines |
If this is |
... |
Arguments not used by this projection. |
Two-column matrix with columns labelled x and
y of locations of projection of coordinates on plane
David Sterratt
https://en.wikipedia.org/wiki/Map_projection, http://mathworld.wolfram.com/SinusoidalProjection.html
Convert points in 3D cartesian space to locations of points on sphere in ‘dual-wedge’ coordinates (fx, fy). Wedges are defined by planes inclined at angle running through a line between poles on the rim above the x axis or the y-axis. fx and fy are the fractional distances along the circle defined by the intersection of this plane and the curtailed sphere.
sphere.cart.to.sphere.dualwedge(P, phi0, R = 1)sphere.cart.to.sphere.dualwedge(P, phi0, R = 1)
P |
locations of points on sphere as N-by-3 matrix with
labelled columns |
phi0 |
rim angle as colatitude |
R |
radius of sphere |
2-column Matrix of ‘wedge’ coordinates of points on
sphere. Column names are phi and lambda.
David Sterratt
Convert locations on the surface of a sphere in cartesian (X, Y, Z) coordinates to spherical (phi, lambda) coordinates.
sphere.cart.to.sphere.spherical(P, R = 1)sphere.cart.to.sphere.spherical(P, R = 1)
P |
locations of points on sphere as N-by-3 matrix with labelled columns "X", "Y" and "Z" |
R |
radius of sphere |
It is assumed that all points are lying on the surface of a sphere of radius R.
N-by-2 Matrix with columns ("phi" and "lambda") of locations of points in spherical coordinates
David Sterratt
Convert points in 3D cartesian space to locations of points on sphere in 'wedge' coordinates (psi, f). Wedges are defined by planes inclined at an angle psi running through a line between poles on the rim above the x axis. f is the fractional distance along the circle defined by the intersection of this plane and the curtailed sphere.
sphere.cart.to.sphere.wedge(P, phi0, R = 1)sphere.cart.to.sphere.wedge(P, phi0, R = 1)
P |
locations of points on sphere as N-by-3 matrix with labelled columns "X", "Y" and "Z" |
phi0 |
rim angle as colatitude |
R |
radius of sphere |
2-column Matrix of 'wedge' coordinates of points on
sphere. Column names are phi and lambda.
David Sterratt
This is the inverse of polar.cart.to.sphere.spherical
sphere.spherical.to.polar.cart(r, pa = FALSE, preserve = "latitude")sphere.spherical.to.polar.cart(r, pa = FALSE, preserve = "latitude")
r |
2-column Matrix of spherical coordinates of points on
sphere. Column names are |
pa |
If |
preserve |
Quantity to preserve locally in the
projection. Options are |
2-column Matrix of Cartesian coordinates of points on polar
projection. Column names should be x and y
David Sterratt
Convert locations of points on sphere in spherical coordinates to points in 3D cartesian space
sphere.spherical.to.sphere.cart(Ps, R = 1)sphere.spherical.to.sphere.cart(Ps, R = 1)
Ps |
N-by-2 matrix with columns containing latitudes
( |
R |
radius of sphere |
An N-by-3 matrix in which each row is the cartesian (X, Y, Z) coordinates of each point
David Sterratt
This uses L'Hullier's theorem to compute the spherical excess and hence the area of the spherical triangle.
sphere.tri.area(P, Tr)sphere.tri.area(P, Tr)
P |
2-column matrix of vertices of triangles given in
spherical polar coordinates. Columns need to be labelled
|
Tr |
3-column matrix of indices of rows of |
Vectors of areas of triangles in units of steradians
David Sterratt
Wolfram MathWorld http://mathworld.wolfram.com/SphericalTriangle.html and http://mathworld.wolfram.com/SphericalExcess.html
## Something that should be an eighth of a sphere, i.e. pi/2 P <- cbind(phi=c(0, 0, pi/2), lambda=c(0, pi/2, pi/2)) Tr <- cbind(1, 2, 3) ## The result of this should be 0.5 print(sphere.tri.area(P, Tr)/pi) ## Now a small triangle P1 <- cbind(phi=c(0, 0, 0.01), lambda=c(0, 0.01, 0.01)) Tr1 <- cbind(1, 2, 3) ## The result of this should approximately 0.01^2/2 print(sphere.tri.area(P, Tr)/(0.01^2/2)) ## Now check that it works for both P <- rbind(P, P1) Tr <- rbind(1:3, 4:6) ## Should have two components print(sphere.tri.area(P, Tr))## Something that should be an eighth of a sphere, i.e. pi/2 P <- cbind(phi=c(0, 0, pi/2), lambda=c(0, pi/2, pi/2)) Tr <- cbind(1, 2, 3) ## The result of this should be 0.5 print(sphere.tri.area(P, Tr)/pi) ## Now a small triangle P1 <- cbind(phi=c(0, 0, 0.01), lambda=c(0, 0.01, 0.01)) Tr1 <- cbind(1, 2, 3) ## The result of this should approximately 0.01^2/2 print(sphere.tri.area(P, Tr)/(0.01^2/2)) ## Now check that it works for both P <- rbind(P, P1) Tr <- rbind(1:3, 4:6) ## Should have two components print(sphere.tri.area(P, Tr))
This in the inverse of sphere.cart.to.sphere.wedge
sphere.wedge.to.sphere.cart(psi, f, phi0, R = 1)sphere.wedge.to.sphere.cart(psi, f, phi0, R = 1)
psi |
vector of slice angles of N points |
f |
vector of fractional distances of N points |
phi0 |
rim angle as colatitude |
R |
radius of sphere |
An N-by-3 matrix in which each row is the cartesian (X, Y, Z) coordinates of each point
David Sterratt
Project spherical coordinate system to a polar
coordinate system such that the area of each
small region is preserved.
spherical.to.polar.area(phi, R = 1)spherical.to.polar.area(phi, R = 1)
phi |
Latitude |
R |
Radius |
This requires
. Hence
. Solving gives
and hence
.
As a check, consider that total area needs to be preserved. If
is maximum value of new variable then
. So
, which agrees with the formula
above.
Coordinate rho that has the dimensions of length
David Sterratt
Spherical plot of reconstructed outline
sphericalplot(r, ...)sphericalplot(r, ...)
r |
Object inheriting |
... |
Parameters depending on class of |
David Sterratt
Draw a spherical plot of reconstructed outline. This method just draws the mesh.
## S3 method for class 'ReconstructedOutline' sphericalplot(r, strain = FALSE, surf = TRUE, ...)## S3 method for class 'ReconstructedOutline' sphericalplot(r, strain = FALSE, surf = TRUE, ...)
r |
|
strain |
If |
surf |
If |
... |
Other graphics parameters – not used at present |
David Sterratt
A StitchedOutline contains a function to stitch the
tears and fullcuts, setting the correspondences hf, hb and
h
retistruct::OutlineCommon -> retistruct::Outline -> retistruct::PathOutline -> retistruct::AnnotatedOutline -> retistruct::TriangulatedOutline -> StitchedOutline
Rsetthe set of points on the rim
TFsetlist containing indices of points in each forward tear
CFsetlist containing indices of points in each forward cut
epsilonthe minimum distance between points, set automatically
tearsStitchedBoolean indicating if tears have been stitched
fullCutsStitchedBoolean indicating if full cuts have been stitched
retistruct::OutlineCommon$clearFeatureSets()retistruct::OutlineCommon$getFeatureSet()retistruct::OutlineCommon$getFeatureSetTypes()retistruct::OutlineCommon$getFeatureSets()retistruct::OutlineCommon$getIDs()retistruct::Outline$addFeatureSet()retistruct::Outline$getDepth()retistruct::Outline$getFragment()retistruct::Outline$getFragmentIDs()retistruct::Outline$getFragmentIDsFromPointIDs()retistruct::Outline$getFragmentPointIDs()retistruct::Outline$getFragmentPoints()retistruct::Outline$getImage()retistruct::Outline$getOutlineLengths()retistruct::Outline$getOutlineSet()retistruct::Outline$getPoints()retistruct::Outline$getPointsScaled()retistruct::Outline$getPointsXY()retistruct::Outline$mapFragment()retistruct::Outline$mapPids()retistruct::Outline$replaceImage()retistruct::PathOutline$insertPoint()retistruct::PathOutline$nextPoint()retistruct::PathOutline$stitchSubpaths()retistruct::AnnotatedOutline$addFullCut()retistruct::AnnotatedOutline$addPoints()retistruct::AnnotatedOutline$addTear()retistruct::AnnotatedOutline$checkTears()retistruct::AnnotatedOutline$computeFullCutRelationships()retistruct::AnnotatedOutline$computeTearRelationships()retistruct::AnnotatedOutline$ensureFixedPointInRim()retistruct::AnnotatedOutline$getFixedPoint()retistruct::AnnotatedOutline$getFullCut()retistruct::AnnotatedOutline$getFullCuts()retistruct::AnnotatedOutline$getRimLengths()retistruct::AnnotatedOutline$getRimSet()retistruct::AnnotatedOutline$getTear()retistruct::AnnotatedOutline$getTears()retistruct::AnnotatedOutline$labelFullCutPoints()retistruct::AnnotatedOutline$labelTearPoints()retistruct::AnnotatedOutline$removeFullCut()retistruct::AnnotatedOutline$removeTear()retistruct::AnnotatedOutline$setFixedPoint()retistruct::AnnotatedOutline$whichFullCut()retistruct::AnnotatedOutline$whichTear()retistruct::TriangulatedOutline$mapTriangulatedFragment()retistruct::TriangulatedOutline$triangulate()new()
Constructor
StitchedOutline$new(...)
...Parameters to superclass constructors
stitchTears()
Stitch together the incisions and tears by inserting new points in the tears and creating correspondences between new points.
StitchedOutline$stitchTears()
stitchFullCuts()
Stitch together the fullcuts by inserting new points in the tears and creating correspondences between new points.
StitchedOutline$stitchFullCuts()
isStitched()
Test if the outline has been stitched
StitchedOutline$isStitched()
Boolean, indicating if the outline has been stitched or not
getBoundarySets()
Get point IDs of points on boundaries
StitchedOutline$getBoundarySets()
List of Point IDs of points on the boundaries.
If the outline has been stitched,
the point IDs in each
element of the list will be ordered in the direction of the
forward pointer, and the boundary that is longest will be
named as Rim. If the outline has not been stitched,
the list will have one element named Rim.
clone()
The objects of this class are cloneable with this method.
StitchedOutline$clone(deep = FALSE)
deepWhether to make a deep clone.
David Sterratt
Generate colours for strain plots
strain.colours(x)strain.colours(x)
x |
Vector of values of log strain |
Vector of colours corresponding to strains
David Sterratt
Stretch the mesh in the flat retina to a circular outline
stretchMesh(Cu, L, i.fix, P.fix)stretchMesh(Cu, L, i.fix, P.fix)
Cu |
Edge matrix |
L |
Lengths in flat outline |
i.fix |
Indices of fixed points |
P.fix |
Coordinates of fixed points |
New matrix of 2D point locations
David Sterratt
Area of triangles on a plane
tri.area(P, Tr)tri.area(P, Tr)
P |
3-column matrix of vertices of triangles |
Tr |
3-column matrix of indices of rows of |
Vectors of areas of triangles
David Sterratt
"Signed area" of triangles on a plane
tri.area.signed(P, Tr)tri.area.signed(P, Tr)
P |
3-column matrix of vertices of triangles |
Tr |
3-column matrix of indices of rows of |
Vectors of signed areas of triangles. Positive sign indicates points are anticlockwise direction; negative indicates clockwise.
David Sterratt
A TriangulatedFragment contains a function to create a triangulated mesh over an fragment, and fields to hold the mesh information.
retistruct::Fragment -> TriangulatedFragment
Tr3 column matrix in which each row contains IDs of points of each triangle
AArea of each triangle in the mesh - has same number of
elements as there are rows of T
Cu2 column matrix in which each row contains IDs of points of edge in mesh
LLength of each edge in the mesh - has same number of
elements as there are rows of Cu
A.signedSigned area of each triangle generated using
tri.area.signed. Positive sign indicates points are
anticlockwise direction; negative indicates clockwise.
new()
Constructor
TriangulatedFragment$new( fragment, n = 200, suppress.external.steiner = FALSE, report = message )
fragmentFragment to triangulate
nMinimum number of points in the triangulation
suppress.external.steinerIf TRUE prevent the
addition of points in the outline. This happens to maintain
triangle quality.
reportFunction to report progress
clone()
The objects of this class are cloneable with this method.
TriangulatedFragment$clone(deep = FALSE)
deepWhether to make a deep clone.
David Sterratt
A TriangulatedOutline contains a function to create a
triangulated mesh over an outline, and fields to hold the mesh
information. Note that areas and lengths are all scaled using
the value of the scale field.
retistruct::OutlineCommon -> retistruct::Outline -> retistruct::PathOutline -> retistruct::AnnotatedOutline -> TriangulatedOutline
Tr3 column matrix in which each row contains IDs of points of each triangle
AArea of each triangle in the mesh - has same number of
elements as there are rows of T
A.totTotal area of the mesh
Cu2 column matrix in which each row contains IDs of
LLength of each edge in the mesh - has same number of
elements as there are rows of Cu
retistruct::OutlineCommon$clearFeatureSets()retistruct::OutlineCommon$getFeatureSet()retistruct::OutlineCommon$getFeatureSetTypes()retistruct::OutlineCommon$getFeatureSets()retistruct::OutlineCommon$getIDs()retistruct::Outline$addFeatureSet()retistruct::Outline$getDepth()retistruct::Outline$getFragment()retistruct::Outline$getFragmentIDs()retistruct::Outline$getFragmentIDsFromPointIDs()retistruct::Outline$getFragmentPointIDs()retistruct::Outline$getFragmentPoints()retistruct::Outline$getImage()retistruct::Outline$getOutlineLengths()retistruct::Outline$getOutlineSet()retistruct::Outline$getPoints()retistruct::Outline$getPointsScaled()retistruct::Outline$getPointsXY()retistruct::Outline$mapFragment()retistruct::Outline$mapPids()retistruct::Outline$replaceImage()retistruct::PathOutline$insertPoint()retistruct::PathOutline$nextPoint()retistruct::PathOutline$stitchSubpaths()retistruct::AnnotatedOutline$addFullCut()retistruct::AnnotatedOutline$addPoints()retistruct::AnnotatedOutline$addTear()retistruct::AnnotatedOutline$checkTears()retistruct::AnnotatedOutline$computeFullCutRelationships()retistruct::AnnotatedOutline$computeTearRelationships()retistruct::AnnotatedOutline$ensureFixedPointInRim()retistruct::AnnotatedOutline$getBoundarySets()retistruct::AnnotatedOutline$getFixedPoint()retistruct::AnnotatedOutline$getFullCut()retistruct::AnnotatedOutline$getFullCuts()retistruct::AnnotatedOutline$getRimLengths()retistruct::AnnotatedOutline$getRimSet()retistruct::AnnotatedOutline$getTear()retistruct::AnnotatedOutline$getTears()retistruct::AnnotatedOutline$initialize()retistruct::AnnotatedOutline$labelFullCutPoints()retistruct::AnnotatedOutline$labelTearPoints()retistruct::AnnotatedOutline$removeFullCut()retistruct::AnnotatedOutline$removeTear()retistruct::AnnotatedOutline$setFixedPoint()retistruct::AnnotatedOutline$whichFullCut()retistruct::AnnotatedOutline$whichTear()triangulate()
Triangulate (mesh) outline
TriangulatedOutline$triangulate(n = 200, suppress.external.steiner = FALSE)
nDesired number of points in mesh
suppress.external.steinerBoolean variable describing whether to insert external Steiner points - see TriangulatedFragment
mapTriangulatedFragment()
Map the point IDs of a TriangulatedFragment on the point IDs of this Outline
TriangulatedOutline$mapTriangulatedFragment(fragment, pids)
fragmentTriangulatedFragment to map
pidsPoint IDs in TriangulatedOutline of points in TriangulatedFragment
clone()
The objects of this class are cloneable with this method.
TriangulatedOutline$clone(deep = FALSE)
deepWhether to make a deep clone.
David Sterratt
P <- rbind(c(1,1), c(2,1), c(2,-1), c(1,-1), c(1,-2), c(-1,-2), c(-1,-1), c(-2,-1),c(-2,1), c(-1,1), c(-1,2), c(1,2)) o <- TriangulatedOutline$new(P) o$addTear(c(3, 4, 5)) o$addTear(c(6, 7, 8)) o$addTear(c(9, 10, 11)) o$addTear(c(12, 1, 2)) flatplot(o) P <- list(rbind(c(1,1), c(2,1), c(2.5,2), c(3,1), c(4,1), c(1,4)), rbind(c(-1,1), c(-1,4), c(-2,3), c(-2,2), c(-3,2), c(-4,1)), rbind(c(-4,-1), c(-1,-1), c(-1,-4)), rbind(c(1,-1), c(2,-1), c(2.5,-2), c(3,-1), c(4,-1), c(1,-4))) o <- TriangulatedOutline$new(P) ##' o$addTear(c(2, 3, 4)) o$addTear(c(17, 18, 19)) o$addTear(c(9, 10, 11)) o$addFullCut(c(1, 5, 16, 20)) flatplot(o)P <- rbind(c(1,1), c(2,1), c(2,-1), c(1,-1), c(1,-2), c(-1,-2), c(-1,-1), c(-2,-1),c(-2,1), c(-1,1), c(-1,2), c(1,2)) o <- TriangulatedOutline$new(P) o$addTear(c(3, 4, 5)) o$addTear(c(6, 7, 8)) o$addTear(c(9, 10, 11)) o$addTear(c(12, 1, 2)) flatplot(o) P <- list(rbind(c(1,1), c(2,1), c(2.5,2), c(3,1), c(4,1), c(1,4)), rbind(c(-1,1), c(-1,4), c(-2,3), c(-2,2), c(-3,2), c(-4,1)), rbind(c(-4,-1), c(-1,-1), c(-1,-4)), rbind(c(1,-1), c(2,-1), c(2.5,-2), c(3,-1), c(4,-1), c(1,-4))) o <- TriangulatedOutline$new(P) ##' o$addTear(c(2, 3, 4)) o$addTear(c(17, 18, 19)) o$addTear(c(9, 10, 11)) o$addFullCut(c(1, 5, 16, 20)) flatplot(o)
The Shiny UI element, runs on a browser and is similar to HTML, attempted to mimic the original Retistruct UI as closely as possible.
ui()ui()
Jan Okul
Vector norm
vecnorm(X)vecnorm(X)
X |
Vector or matrix. |
If a vector, returns the 2-norm of the vector. If a matrix, returns the 2-norm of each row of the matrix
David Sterratt