Skip to content
This repository has been archived by the owner on Mar 27, 2024. It is now read-only.

Commit

Permalink
Finished the port of ConvoluteReflectivityTileFixedLambda. It is off by
Browse files Browse the repository at this point in the history 
at most 4% but on average about 1e-4. There are some problems on the end
points of the range.

Signed-off-by: Jay Jay Billings <billingsjj@ornl.gov>
  • Loading branch information
Jay Jay Billings committed Apr 13, 2015
1 parent e638d7b commit dd192f4
Show file tree
Hide file tree
Showing 2 changed files with 127 additions and 39 deletions.
128 changes: 103 additions & 25 deletions ...org.eclipse.ice.reflectivity/src/org/eclipse/ice/reflectivity/ReflectivityCalculator.java
View file
Expand Up @@ -52,7 +52,7 @@ public class ReflectivityCalculator {
* @param wavelength
* the wavelength of the incident neutrons
* @param tiles
* the list of TIles that contains the physical parameters needed
* the list of Tiles that contains the physical parameters needed
* for the calculation, including the scattering densities,
* absorption parameters and thicknesses.
* @return the squared modulus of the specular reflectivity
Expand Down Expand Up @@ -125,11 +125,11 @@ public double getModSqrdSpecRef(double waveVectorQ, double wavelength,
* @param q
* the wave vector (Q) plus additional space for the convolution.
* This array should have length = numPoints + numLowPoints.
* @param delQ0
* @param deltaQ0
* the zeroth order term of a Taylor expansion of the
* reflectometer resolution function dQ = dQ_0 + (dQ/Q)_1 x Q +
* ...
* @param delQ1oQ
* @param deltaQ1ByQ
* the zeroth order term of the Q resolution Taylor expansion
* @param wavelength
* the wavelength of the incident neutrons
Expand All @@ -147,9 +147,9 @@ public double getModSqrdSpecRef(double waveVectorQ, double wavelength,
* OUTPUT - the specular reflectivity values for each Q in q
* convoluted with instrumental resolution.
*/
public void convolute(double[] waveVector, double delQ0, double delQ1oQ,
double wavelength, int numPoints, int numLowPoints,
int numHighPoints, double[] refFit) {
public void convolute(double[] waveVector, double deltaQ0,
double deltaQ1ByQ, double wavelength, int numPoints,
int numLowPoints, int numHighPoints, double[] refFit) {

double ln2 = Math.log(2.0);
double qEff = 0.0, qRes = 0.0, rExp = 0.0, rNorm = 0.0;
Expand All @@ -165,7 +165,7 @@ public void convolute(double[] waveVector, double delQ0, double delQ1oQ,
} else {
qEff = waveVector[i];
}
double qDel = delQ0 + qEff * delQ1oQ;
double qDel = deltaQ0 + qEff * deltaQ1ByQ;
double twSgSq = 2.0 * qDel * qDel / (8.0 * ln2);
if (twSgSq < 1.0e-10) {
twSgSq = 1.0e-10;
Expand Down Expand Up @@ -229,25 +229,25 @@ public void convolute(double[] waveVector, double delQ0, double delQ1oQ,
* @param q
* the wave vector (Q) plus additional space for the convolution.
* This array should have length = numPoints + numLowPoints.
* @param delQ0
* @param deltaQ0
* the zeroth order term of a Taylor expansion of the
* reflectometer resolution function dQ = dQ_0 + (dQ/Q)_1 x Q +
* ...
* @param delQ1oQ
* @param deltaQ1ByQ
* the zeroth order term of the Q resolution Taylor expansion
* @param numPoints
* the number of points in the wave vector
* @return numLowPoints the number of points in the low-Q extension to q
* used for convolution of the data with the resolution function.
* Returned by ExtResFixedLambda.
*/
public int getLowExtensionLength(double[] waveVector, double delQ0,
double delQ1oQ, int numPoints) {
public int getLowExtensionLength(double[] waveVector, double deltaQ0,
double deltaQ1ByQ, int numPoints) {

double ln2 = Math.log(2.0);

// Determine the loq-Q extension
double qDel = delQ0 + waveVector[0] * delQ1oQ;
double qDel = deltaQ0 + waveVector[0] * deltaQ1ByQ;
double qStep = waveVector[1] - waveVector[0];
double twSgSq = Math.max(2.0 * qDel * qDel / (8.0 * ln2), 1.0e-10);
int numLowPoints = 0;
Expand All @@ -268,25 +268,25 @@ public int getLowExtensionLength(double[] waveVector, double delQ0,
* @param q
* the wave vector (Q) plus additional space for the convolution.
* This array should have length = numPoints + numLowPoints.
* @param delQ0
* @param deltaQ0
* the zeroth order term of a Taylor expansion of the
* reflectometer resolution function dQ = dQ_0 + (dQ/Q)_1 x Q +
* ...
* @param delQ1oQ
* @param deltaQ1ByQ
* the zeroth order term of the Q resolution Taylor expansion
* @param numPoints
* the number of points in the wave vector
* @return numHighPoints the number of points in the high-Q extension to q
* used for convolution of the data with the resolution function.
* Returned by ExtResFixedLambda.
*/
public int getHighExtensionLength(double[] waveVector, double delQ0,
double delQ1oQ, int numPoints) {
public int getHighExtensionLength(double[] waveVector, double deltaQ0,
double deltaQ1ByQ, int numPoints) {

double ln2 = Math.log(2.0);

// Determine the high-Q extension
double qDel = delQ0 + waveVector[numPoints - 1] * delQ1oQ;
double qDel = deltaQ0 + waveVector[numPoints - 1] * deltaQ1ByQ;
double qStep = waveVector[numPoints - 1] - waveVector[numPoints - 2];
double twSgSq = 2.0 * qDel * qDel / (8.0 * ln2);
int numHighPoints = 0;
Expand Down Expand Up @@ -438,8 +438,8 @@ public Tile[] generateTiles(Slab[] slabs, int numRough, double[] zInt,
// negative.
gDMid = refSlab.thickness - 0.5 * totalThickness
* (refSlab.interfaceWidth + secondRefSlab.interfaceWidth);
System.out.println("GDMid[" + i + "] = " + gDMid + ", " + refSlab.thickness
+ ", " + refSlab.interfaceWidth + ", "
System.out.println("GDMid[" + i + "] = " + gDMid + ", "
+ refSlab.thickness + ", " + refSlab.interfaceWidth + ", "
+ secondRefSlab.interfaceWidth);
if (gDMid <= 1.0e-10) {
// The interfaces are overlapping. Step through the entire slab
Expand Down Expand Up @@ -486,8 +486,8 @@ public Tile[] generateTiles(Slab[] slabs, int numRough, double[] zInt,
tmpSlab = generatedSlabs[nGlay + j - numRough / 2 - 1];
updateTileByInterface(tmpSlab, thirdRefSlab, refSlab,
zInt[j], rufInt[j]);
System.out.println("q = " + (nGlay + j - numRough/2 -1) + " "
+ tmpSlab.scatteringLength);
System.out.println("q = " + (nGlay + j - numRough / 2 - 1)
+ " " + tmpSlab.scatteringLength);
}
nGlay += numRough / 2 + 1;
// Central, bulk-like portion
Expand Down Expand Up @@ -520,7 +520,7 @@ public Tile[] generateTiles(Slab[] slabs, int numRough, double[] zInt,
tmpSlab = generatedSlabs[nGlay + i - numRough / 2 - 1];
updateTileByInterface(tmpSlab, secondRefSlab, refSlab, zInt[i],
rufInt[i]);
System.out.println("q = " + (nGlay + i - numRough/2 - 1) + " "
System.out.println("q = " + (nGlay + i - numRough / 2 - 1) + " "
+ tmpSlab.scatteringLength);
}
nGlay += numRough / 2 + 1;
Expand Down Expand Up @@ -624,8 +624,86 @@ private double getTileValue(double xm1, double x, double xp1, double tao,
return 0.5 * (xm1 + x + tao * (x - xm1) + x + xp1 + beta * (xp1 - x))
- x;
}

public void convoluteReflectivity() {


/**
* This operation computes the convolution of the reflectivity with a
* variable Gaussian resolution function.
*
* @param deltaQ0
* - FIXME!
* @param deltaQ1ByQ
* - FIXME!
* @param wavelength
* - FIXME!
* @param getRQ4
* - FIXME! True if the routine should compute rq^4, false
* otherwise.
* @param waveVector
* The wave vector - FIXME!
* @param tiles
* The tiles that define the layered structure of the materials.
*/
public double[] convoluteReflectivity(double deltaQ0, double deltaQ1ByQ,
double wavelength, boolean getRQ4, double[] waveVector, Tile[] tiles) {

// Local Declarations
double qEff = 0.0;
int numPoints = waveVector.length, numLowPoints = 0, numHighPoints = 0;
double[] reflectivity = new double[numPoints];

// Determine the length of the high- and low-Q extensions
numLowPoints = getLowExtensionLength(waveVector, deltaQ0, deltaQ1ByQ,
numPoints);
numHighPoints = getHighExtensionLength(waveVector, deltaQ0, deltaQ1ByQ,
numPoints);

// Extend the wave vector in a temporary array
double[] tempWaveVector = new double[numLowPoints + numHighPoints
+ numPoints];
double waveVecStep = waveVector[1] - waveVector[0];
for (int i = 0; i < numLowPoints; i++) {
tempWaveVector[i] = waveVector[0] - waveVecStep
* ((double) numLowPoints + 1 - i);
}
for (int i = 0; i < numPoints; i++) {
tempWaveVector[numLowPoints + i] = waveVector[i];
}
waveVecStep = waveVector[numPoints - 1] - waveVector[numPoints - 2];
for (int i = 0; i < numHighPoints; i++) {
tempWaveVector[i + numLowPoints + numPoints] = waveVector[numPoints - 1]
+ waveVecStep * ((double) i);
}

// Create a temporary array to hold the extended reflectivity.
double[] tempReflectivity = new double[numPoints+numLowPoints+numHighPoints];
// Generate reflectivity values for convolution.
// Calculate perfect-resolution reflectivity on extended wave vector
for (int i = 0; i < numPoints + numLowPoints + numHighPoints; i++) {
if (tempWaveVector[i] < 1.0e-10) {
qEff = 1.0e-10;
} else {
qEff = tempWaveVector[i];
}
tempReflectivity[i] = getModSqrdSpecRef(qEff, wavelength, tiles);
}

// Convolve with instrumental resolution
convolute(tempWaveVector, deltaQ0, deltaQ1ByQ, wavelength, numPoints,
numLowPoints, numHighPoints, tempReflectivity);

// Transfer the results to the reflectivity array.
for (int i = 0; i < numPoints; i++) {
reflectivity[i] = tempReflectivity[i];
}
// Calculate RQ^4 if needed. FIXME! This is not covered by the current
// tests and should actually be a separate function!
if (getRQ4) {
for (int i = 0; i < numPoints; i++) {
reflectivity[i] = Math.pow(waveVector[i], 4.0)
* reflectivity[i];
}
}

return reflectivity;
}
}
38 changes: 24 additions & 14 deletions ...reflectivity.test/src/org/eclipse/ice/reflectivity/test/ReflectivityCalculatorTester.java
View file
Expand Up @@ -475,30 +475,40 @@ public void testConvoluteReflectivity() {
form = reader.read(project.getFile("conRefTileFixedLambda.csv"));
ListComponent<String[]> refLines = (ListComponent<String[]>) form
.getComponent(1);
assertEquals(405, refLines.size());
assertEquals(403, refLines.size());

// Get the parameters
double q0 = Double.valueOf(refLines.get(0)[0]);
double q1ByQ = Double.valueOf(refLines.get(0)[1]);
double deltaQ0 = Double.valueOf(refLines.get(0)[0]);
double deltaQ1ByQ = Double.valueOf(refLines.get(0)[1]);
double lambda = Double.valueOf(refLines.get(0)[2]);
boolean getRQ4 = Boolean.valueOf(refLines.get(0)[3]);

// Get the reference reflectivity results and the q input array
double[] q = new double[402];
double[] refReflectivity = new double[404];
double[] waveVector = new double[402];
double[] refReflectivity = new double[402];

// Load the reference arrays. The reference reflectivity is longer than
// Q so they should be split up.
for (int i = 1; i < 405; i++) {
String[] line = refLines.get(i);
refReflectivity[i-1] = Double.valueOf(line[0]);
}
// Load q
// Load the reference arrays.
for (int i = 1; i < 403; i++) {
String[] line = refLines.get(i);
q[i-1] = Double.valueOf(line[1]);
refReflectivity[i - 1] = Double.valueOf(line[0]);
waveVector[i - 1] = Double.valueOf(line[1]);
}

// Do the calculation
ReflectivityCalculator calc = new ReflectivityCalculator();
double[] reflectivity = calc.convoluteReflectivity(deltaQ0, deltaQ1ByQ,
lambda, getRQ4, waveVector, refTiles);

// Check the results
assertEquals(refReflectivity.length, reflectivity.length);
// FIXME! - Do value comparisons!
for (int i = 0; i < refReflectivity.length; i++) {
// Three percent error is the best that I can do on this because the
// behavior at the end points does not match. I think it is due to
// an accumulation of errors from other sources.
assertEquals(refReflectivity[i], reflectivity[i], Math.abs(refReflectivity[i])*0.04);
}

return;
}

Expand Down

0 comments on commit dd192f4

Please sign in to comment.