use of ffx.numerics.ComplexNumber in project ffx by mjschnie.
the class SigmaAMinimize method setWEstimate.
private void setWEstimate() {
// generate initial w estimate
ReflectionSpline spline = new ReflectionSpline(reflectionList, refinementData.nbins);
int[] nmean = new int[refinementData.nbins];
for (int i = 0; i < refinementData.nbins; i++) {
nmean[i] = 0;
}
double mean = 0.0, tot = 0.0;
double[][] fc = refinementData.fctot;
double[][] fo = refinementData.fsigf;
for (HKL ih : reflectionList.hkllist) {
int i = ih.index();
if (ih.allowed() == 0.0 || Double.isNaN(fc[i][0]) || Double.isNaN(fo[i][0])) {
continue;
}
double s2 = Crystal.invressq(crystal, ih);
double epsc = ih.epsilonc();
ComplexNumber fct = new ComplexNumber(fc[i][0], fc[i][1]);
double ecscale = spline.f(s2, refinementData.fcesq);
double eoscale = spline.f(s2, refinementData.foesq);
double ec = fct.times(sqrt(ecscale)).abs();
double eo = fo[i][0] * sqrt(eoscale);
double wi = pow(eo - ec, 2.0) / epsc;
nmean[spline.i1()]++;
tot++;
x[spline.i1() + refinementData.nbins] += (wi - x[spline.i1() + refinementData.nbins]) / nmean[spline.i1()];
mean += (wi - mean) / tot;
}
logger.info(String.format(" Starting mean w: %8.3f", mean));
logger.info(String.format(" Starting w scaling: %8.3f", 1.0 / mean));
for (int i = 0; i < refinementData.nbins; i++) {
x[i] -= x[i + refinementData.nbins];
x[i] *= scaling[i];
scaling[i + refinementData.nbins] = 1.0 / mean;
x[i + refinementData.nbins] *= scaling[i + refinementData.nbins];
}
sigmaAEnergy.setScaling(scaling);
}
use of ffx.numerics.ComplexNumber in project ffx by mjschnie.
the class CrystalStats method getRFree.
/**
* simply return the current Rfree value
*
* @return rfree value as a percent
*/
public double getRFree() {
double sum = 0.0;
double sumfo = 0.0;
for (HKL ih : reflectionlist.hkllist) {
int i = ih.index();
// ignored cases
if (Double.isNaN(fc[i][0]) || Double.isNaN(fo[i][0]) || fo[i][1] <= 0.0) {
continue;
}
// spline setup
double ss = Crystal.invressq(crystal, ih);
double fh = spline.f(ss, refinementdata.spline);
ComplexNumber c = new ComplexNumber(fc[i][0], fc[i][1]);
if (refinementdata.isFreeR(i)) {
sum += abs(abs(fo[i][0]) - fh * abs(c.abs()));
sumfo += abs(fo[i][0]);
}
}
rfree = (sum / sumfo) * 100.0;
return rfree;
}
use of ffx.numerics.ComplexNumber in project ffx by mjschnie.
the class CrystalStats method getR.
/**
* simply return the current R value
*
* @return r value as a percent
*/
public double getR() {
double numer = 0.0;
double denom = 0.0;
double sum = 0.0;
double sumfo = 0.0;
double sumall = 0.0;
double sumfoall = 0.0;
for (HKL ih : reflectionlist.hkllist) {
int i = ih.index();
// ignored cases
if (Double.isNaN(fc[i][0]) || Double.isNaN(fo[i][0]) || fo[i][1] <= 0.0) {
continue;
}
// spline setup
double ss = Crystal.invressq(crystal, ih);
double fh = spline.f(ss, refinementdata.spline);
ComplexNumber c = new ComplexNumber(fc[i][0], fc[i][1]);
numer = abs(abs(fo[i][0]) - fh * abs(c.abs()));
denom = abs(fo[i][0]);
sumall += numer;
sumfoall += denom;
if (!refinementdata.isFreeR(i)) {
sum += numer;
sumfo += denom;
}
}
rall = (sumall / sumfoall) * 100.0;
r = (sum / sumfo) * 100.0;
return r;
}
use of ffx.numerics.ComplexNumber in project ffx by mjschnie.
the class MTZFilter method readFcs.
/**
* Read the structure factors.
*
* @param mtzFile
* @param reflectionList
* @param fcData
* @param properties
* @return
*/
public boolean readFcs(File mtzFile, ReflectionList reflectionList, DiffractionRefinementData fcData, CompositeConfiguration properties) {
int nRead, nIgnore, nRes, nFriedel, nCut;
ByteOrder byteOrder = ByteOrder.nativeOrder();
FileInputStream fileInputStream;
DataInputStream dataInputStream;
StringBuilder sb = new StringBuilder();
try {
fileInputStream = new FileInputStream(mtzFile);
dataInputStream = new DataInputStream(fileInputStream);
byte[] headerOffset = new byte[4];
byte[] bytes = new byte[80];
int offset = 0;
// Eat "MTZ" title.
dataInputStream.read(bytes, offset, 4);
String mtzString = null;
// Header offset.
dataInputStream.read(headerOffset, offset, 4);
// Machine stamp.
dataInputStream.read(bytes, offset, 4);
ByteBuffer byteBuffer = ByteBuffer.wrap(bytes);
int stamp = byteBuffer.order(ByteOrder.BIG_ENDIAN).getInt();
String stampString = Integer.toHexString(stamp);
switch(stampString.charAt(0)) {
case '1':
case '3':
if (byteOrder.equals(ByteOrder.LITTLE_ENDIAN)) {
byteOrder = ByteOrder.BIG_ENDIAN;
}
break;
case '4':
if (byteOrder.equals(ByteOrder.BIG_ENDIAN)) {
byteOrder = ByteOrder.LITTLE_ENDIAN;
}
break;
}
byteBuffer = ByteBuffer.wrap(headerOffset);
int headerOffsetI = byteBuffer.order(byteOrder).getInt();
// Skip to header and parse.
dataInputStream.skipBytes((headerOffsetI - 4) * 4);
for (Boolean parsing = true; parsing; dataInputStream.read(bytes, offset, 80)) {
mtzString = new String(bytes);
parsing = parseHeader(mtzString);
}
// Column identifiers.
fc = phiC = fs = phiS = -1;
boolean print = true;
parseFcColumns(print);
if (h < 0 || k < 0 || l < 0) {
String message = " Fatal error in MTZ file - no H K L indexes?\n";
logger.log(Level.SEVERE, message);
return false;
}
// Reopen to start at beginning.
fileInputStream = new FileInputStream(mtzFile);
dataInputStream = new DataInputStream(fileInputStream);
// Skip initial header.
dataInputStream.skipBytes(80);
float[] data = new float[nColumns];
HKL mate = new HKL();
// Read in data.
ComplexNumber complexNumber = new ComplexNumber();
nRead = nIgnore = nRes = nFriedel = nCut = 0;
for (int i = 0; i < nReflections; i++) {
for (int j = 0; j < nColumns; j++) {
dataInputStream.read(bytes, offset, 4);
byteBuffer = ByteBuffer.wrap(bytes);
data[j] = byteBuffer.order(byteOrder).getFloat();
}
int ih = (int) data[h];
int ik = (int) data[k];
int il = (int) data[l];
boolean friedel = reflectionList.findSymHKL(ih, ik, il, mate, false);
HKL hkl = reflectionList.getHKL(mate);
if (hkl != null) {
if (fc > 0 && phiC > 0) {
complexNumber.re(data[fc] * cos(toRadians(data[phiC])));
complexNumber.im(data[fc] * sin(toRadians(data[phiC])));
fcData.setFc(hkl.index(), complexNumber);
}
if (fs > 0 && phiS > 0) {
complexNumber.re(data[fs] * cos(toRadians(data[phiS])));
complexNumber.im(data[fs] * sin(toRadians(data[phiS])));
fcData.setFs(hkl.index(), complexNumber);
}
nRead++;
} else {
HKL tmp = new HKL(ih, ik, il);
if (!reflectionList.resolution.inInverseResSqRange(Crystal.invressq(reflectionList.crystal, tmp))) {
nRes++;
} else {
nIgnore++;
}
}
}
if (logger.isLoggable(Level.INFO)) {
sb.append(format(" MTZ file type (machine stamp): %s\n", stampString));
sb.append(format(" Fc HKL read in: %d\n", nRead));
sb.append(format(" Fc HKL read as friedel mates: %d\n", nFriedel));
sb.append(format(" Fc HKL NOT read in (too high resolution): %d\n", nRes));
sb.append(format(" Fc HKL NOT read in (not in internal list?): %d\n", nIgnore));
sb.append(format(" Fc HKL NOT read in (F/sigF cutoff): %d\n", nCut));
sb.append(format(" HKL in internal list: %d\n", reflectionList.hkllist.size()));
logger.info(sb.toString());
}
} catch (EOFException e) {
String message = " MTZ end of file reached.";
logger.log(Level.WARNING, message, e);
return false;
} catch (IOException e) {
String message = " MTZ IO Exception.";
logger.log(Level.WARNING, message, e);
return false;
}
return true;
}
use of ffx.numerics.ComplexNumber in project ffx by mjschnie.
the class CrystalReciprocalSpaceTest method test1N7SPermanent.
/**
* Test of permanent method, of class CrystalReciprocalSpace.
*/
@Test
public void test1N7SPermanent() {
String filename = "ffx/xray/structures/1N7S.pdb";
int index = filename.lastIndexOf(".");
String name = filename.substring(0, index);
// load the structure
ClassLoader cl = this.getClass().getClassLoader();
File structure = new File(cl.getResource(filename).getPath());
PotentialsUtils potutil = new PotentialsUtils();
MolecularAssembly mola = potutil.open(structure);
CompositeConfiguration properties = mola.getProperties();
Crystal crystal = new Crystal(39.767, 51.750, 132.938, 90.00, 90.00, 90.00, "P212121");
Resolution resolution = new Resolution(1.45);
ReflectionList reflectionList = new ReflectionList(crystal, resolution);
DiffractionRefinementData refinementData = new DiffractionRefinementData(properties, reflectionList);
mola.finalize(true, mola.getForceField());
ForceFieldEnergy energy = mola.getPotentialEnergy();
List<Atom> atomList = mola.getAtomList();
Atom[] atomArray = atomList.toArray(new Atom[atomList.size()]);
// set up FFT and run it
ParallelTeam parallelTeam = new ParallelTeam();
CrystalReciprocalSpace crs = new CrystalReciprocalSpace(reflectionList, atomArray, parallelTeam, parallelTeam);
crs.computeAtomicDensity(refinementData.fc);
// tests
ComplexNumber b = new ComplexNumber(-828.584, -922.704);
HKL hkl = reflectionList.getHKL(1, 1, 4);
ComplexNumber a = refinementData.getFc(hkl.index());
System.out.println("1 1 4: " + a.toString() + " | " + b.toString() + " | " + a.divides(b).toString());
assertEquals("1 1 4 reflection should be correct", -753.4722104328416, a.re(), 0.0001);
assertEquals("1 1 4 reflection should be correct", -1012.1341308707799, a.im(), 0.0001);
b.re(-70.4582);
b.im(-486.142);
hkl = reflectionList.getHKL(2, 1, 10);
a = refinementData.getFc(hkl.index());
System.out.println("2 1 10: " + a.toString() + " | " + b.toString() + " | " + a.divides(b).toString());
assertEquals("2 1 10 reflection should be correct", -69.39660884054359, a.re(), 0.0001);
assertEquals("2 1 10 reflection should be correct", -412.0147625765328, a.im(), 0.0001);
}
Aggregations