aboutsummaryrefslogtreecommitdiff
diff options
context:
space:
mode:
authorThomas White <taw@physics.org>2011-11-10 16:08:06 +0100
committerThomas White <taw@physics.org>2012-02-22 15:27:39 +0100
commitf70d3c1535f672c9bd313a229f44d52c979a6ad7 (patch)
tree98ad6da2edcb8e7ab09a89c4232765e6d17c94bc
parente21d69a9141e3f64147aa2f68f9c3634d66a493f (diff)
Add scripts/gen-sfs for calculating scattering
-rw-r--r--doc/man/pattern_sim.110
-rwxr-xr-xscripts/gen-sfs80
2 files changed, 85 insertions, 5 deletions
diff --git a/doc/man/pattern_sim.1 b/doc/man/pattern_sim.1
index b9c8fa33..275c25cd 100644
--- a/doc/man/pattern_sim.1
+++ b/doc/man/pattern_sim.1
@@ -26,10 +26,10 @@ The result will be written to an HDF5 file in the current directory with the nam
.SH REFLECTION LISTS
-pattern_sim does not know about symmetry, so your input reflection list
-(give with "-i") must be expanded. You can do this with:
+You'll need to create a file containing the intensities of the reflections. The normal way to do this is to use CCP4 via the "gen-sfs" script in CrystFEL's script folder. Run it like this:
-$ get_hkl -i myfile.hkl -o output.hkl -y mypointgroup -e 1
+$ gen-sfs mymodel.pdb "P6" 3
-Please be sure to read the "Note about Unit Cell Settings" in the documentation
-for indexamajig.
+You need to give the PDB model, the symmetry of the output reflections (use the lowest symmetry space group with the right point group), and optionally the maximum resolution in Angstroms. If you don't specify the resolution, it'll use 3 Angstroms.
+
+The reflections will be output as "mymodel.pdb.hkl" ready for input to pattern_sim. You'll need to give the Laue class of the symmetry you gave to gen-sfs, "6/m" in this case, to pattern_sim with the "-y" option. By default, gen-sfs calculates the values for CuKa radiation (8.3 keV, 1.5 A). It will not calculate the anomalous contribution to scattering, i.e. the differences in intensities between Bijoet pairs. Both of these are the default behaviour for "sfall" in CCP4, so read the manual for that for further details. If you need something different, get the "ano_sfall.com" script from James Holton and turn the output into a CrystFEL reflection list in a similar way.
diff --git a/scripts/gen-sfs b/scripts/gen-sfs
new file mode 100755
index 00000000..bf86fb01
--- /dev/null
+++ b/scripts/gen-sfs
@@ -0,0 +1,80 @@
+#!/bin/sh
+
+PDB=$1
+SYMM=$2
+RESOLUTION=$3
+
+if [ "x$PDB" = "x" ]; then
+ echo "Syntax: $0 <PDB file> <space group> [<resolution>]"
+ exit
+fi
+
+if [ "x$SYMM" = "x" ]; then
+ echo "Syntax: $0 <PDB file> <space group> [<resolution>]"
+ exit
+fi
+
+if [ "x$RESOLUTION" = "x" ]; then
+ echo "Resolution not given. Using 3 Angstroms."
+ RESOLUTION=3
+fi
+
+echo "Running sfall to calculate structure factors..."
+sfall XYZIN $PDB HKLOUT ${PDB}.mtz > gen-sfs.html << EOF
+MODE SFCALC XYZIN
+RESOLUTION $RESOLUTION 1000
+FORM NGAUSS 5
+SYMM $SYMM
+END
+EOF
+if [ $? -ne 0 ]; then exit 1; fi
+
+echo "Running cad to get the right asymmetric unit..."
+cad HKLIN1 ${PDB}.mtz HKLOUT ${PDB}-sorted.mtz >> gen-sfs.html <<EOF
+TITLE Sorted blah
+LABIN FILE 1 E1=FC E2=PHIC
+CTYPE FILE 1 E1=F E2=P
+EOF
+if [ $? -ne 0 ]; then exit 1; fi
+
+echo "Converting structure factors to text..."
+mtz2various hklin ${PDB}-sorted.mtz hklout ${PDB}-temp.hkl >> gen-sfs.html <<EOF
+LABIN H=H K=K L=L FC=FC PHIC=PHIC
+OUTPUT USER '(3I4,2F9.1)'
+EOF
+if [ $? -ne 0 ]; then exit 1; fi
+rm -f ${PDB}.mtz
+perl < ${PDB}-temp.hkl > ${PDB}.hkl << WIBBLE
+use strict;
+
+my \$line;
+open(FILE, "${PDB}-temp.hkl");
+
+printf(" h k l I phase sigma(I) 1/d(nm^-1) ".
+ "counts fs/px ss/px\\n");
+
+while ( \$line = <FILE> ) {
+
+ if ( \$line =~ /^\s*([\d\-]+)\s+([\d\-]+)\s+([\d\-]+)\s+([\d\-\.]+)\s+([\d\-\.]+)/ ) {
+
+ my \$h = \$1;
+ my \$k = \$2;
+ my \$l = \$3;
+ my \$intensity = \$4*\$4; # Square to convert F->I
+ my \$phase = \$5;
+
+ printf("%3i %3i %3i %10.2f %8.2f %10.2f %s %7i %6.1f %6.1f\n",
+ \$h, \$k, \$l, \$intensity, \$phase, 0.0,
+ " -", 1, 0.0, 0.0);
+
+ } else {
+ printf(STDERR "Couldn't understand line '%s'\n", \$line);
+ }
+
+}
+close(FILE);
+printf("End of reflections\n");
+WIBBLE
+exit
+
+rm -f ${PDB}-temp.hkl