Showing posts with label gnu octave. Show all posts
Showing posts with label gnu octave. Show all posts

03 March 2013

350. Installing GNU Octave on Windows XP

This is a Windows XP post (my first?), so right-thinking linux people can move on. Nothing to see here.

For the rest of you:
Since I can't force my students to switch to Linux -- but I can force the to use Octave and Gnuplot -- I need to make sure that there are some easy, step-by-step (recipe-based) guides available that describe how to install the tools that they are required to use on their platform of choice. OS X has many flaws, but has macports. Windows has cygwin, in addition to native builds of most GNU tools. I'll cover cygwin in a separate post.

Here's a pretty straightforward approach to installing GNU Octave and I've tested it on Windows XP. It suffers from a fatal weakness: it takes well over a minute (!) to actually start Octave this way -- each timeUse the cygwin approach insteadhttp://verahill.blogspot.com.au/2013/03/353-cygwin-with-octave-and-gnuplot-on.html -- this way it's about as fast as on Linux.

0. 7zip and notepad++
You'll need 7zip: http://www.7-zip.org/ and http://notepad-plus-plus.org/. Download and install.

1.  Download Octave
Go to http://wiki.octave.org/Octave_for_Windows
The most current version is 3.6.2 ("Octave-3.6.2-mingw + octaveforge pkgs"), which takes you here: http://downloads.sourceforge.net/project/octave/Octave%20Windows%20binaries/Octave%203.6.2%20for%20Windows%20MinGW%20installer/Octave3.6.2_gcc4.6.2_20120609.7z

2. Set up
The downloaded file is 169 Mb but extracts to 700 Mb(!) so even the extraction process takes a little while. Extract with 7zip by right-clicking the downloaded file, selecting 7zip, and extract to "Octave...".

Create C:\octave, and put the Octave3.6.2_gcc4.6.2 folder in it.

Right click on C:\octave\Octave3.6.2_gcc4.6.2\bin\octave.exe and select create shortcut. Copy the shortcut to e.g. your desktop. Edit it and click on change icon -- select "C:\Octave\Octave3.6.2_gcc4.6.2\share\octave\3.6.2\imagelib\octave-logo.ico".

3. Running
Simply double-click on the shortcut.
Here's the bad news: it takes 1 m 20 s to start, which is unacceptable.

4. Packages
To install packages without having a full build environment, download from
http://downloads.sourceforge.net/project/octave/Octave%20Windows%20binaries/Octave%203.6.2%20for%20Windows%20MinGW%20installer/Octave3.6.2_gcc4.6.2_pkgs_20120613.7z

Extract the file, and copy the bin, include, lib, share folders to c:\octave\Octave3.6.2_gcc4.6.2\ so that the folders merge with those from step 2.

Start octave and run
pkg rebuild -noauto nan

This way all packages (except nan) will auto-load next time you start octave.

29 September 2012

249. Quick but precise isotopic pattern (isotope envelope) calculator in Octave

UPDATE: Below is an accurate calculator,  but it is impractically slow for large molecules. A practical AND accurate calculator is found here:http://verahill.blogspot.com.au/2012/10/isotopic-pattern-caculator-in-python.html

Use the post below to learn about the fundamental theory, but then look at the other post to understand how to implement it.

Old post:
Getting fast and accurate isotopic patterns can be tricky using tools available online, for download or which form part of commercial packages. A particular problem is that different tools give slightly different values -- so which one to trust?

The answer: the tool for which you know that the algorithm is sound.

The extreme conclusion of that way of thinking is to write your own calculator.
Below is the conceptual process of calculating the isotopic pattern of a molecule using GNU Octave.

You need the linear algebra package:
sudo apt-get install octave octave-linear-algebra

b is the isotopic distribution for an element, and bb are the masses of those isotopes.

Once you've got a computational engine it's not too difficult to expand it for more general cases, account for charge, and instrument resolution.


Molecule: Cl4

b=[0.7578,0.2422];
bb=[34.96885,36.96885];
e=prod(cartprod(b,b,b,b),2);
ee=sum(cartprod(bb,bb,bb,bb),2);
n=4;
g=histc([ee e],linspace(min(ee),max(ee),n*(max(ee)-min(ee)+1)),2);
h=linspace(min(ee),max(ee),n*(max(ee)-min(ee)+1));
distr=e'*g;
plot(h,100.*distr/max(distr))
[h' (100.*distr/max(distr))']
Here's the output for n=1:
   139.87540    78.22048
   140.87540     0.00000
   141.87540   100.00000
   142.87540     0.00000
   143.87540    47.94141
   144.87540     0.00000
   145.87540    10.21502
   146.87540     0.00000
   147.87540     0.81620

And here's the output from Matt Monroe's calculator:
Isotopic Abundances for Cl4
  Mass/Charge Fraction  Intensity
   139.87541 0.3297755   78.22
   140.87541 0.0000000    0.00
   141.87541 0.4215974  100.00
   142.87541 0.0000000    0.00
   143.87541 0.2021197   47.94
   144.87541 0.0000000    0.00
   145.87541 0.0430662   10.22
   146.87541 0.0000000    0.00
   147.87541 0.0034411    0.82


Another molecule: Li2Cl2

Here's the code:
a=[0.0759,0.9241];
aa=[6.01512,7.01512];
b=[0.7578,0.2422];
bb=[34.96885,36.96885];
e=prod(cartprod(a,a,b,b),2);
ee=sum(cartprod(aa,aa,bb,bb),2);
n=1;
g=histc([ee e],linspace(min(ee),max(ee),n*(max(ee)-min(ee)+1)),2);
h=linspace(min(ee),max(ee),n*(max(ee)-min(ee)+1));
distr=e'*g;
plot(h,100.*distr/max(distr))
[h' (100.*distr/max(distr))']

ans =

    81.96794     0.67170
    82.96794    16.35626
    83.96794   100.00000
    84.96794    10.45523
    85.96794    63.71604
    86.96794     1.67079
    87.96794    10.17116

vs Matt Monroe's calculator:
Isotopic Abundances for Li2Cl2
  Mass/Charge Fraction  Intensity
    81.96795 0.0033082    0.67
    82.96795 0.0805564   16.36
    83.96795 0.4925109  100.00
    84.96795 0.0514932   10.46
    85.96795 0.3138084   63.72
    86.96795 0.0082288    1.67
    87.96795 0.0500941   10.17

We can then expand the code to allow for plotting
a=[0.0759,0.9241];
aa=[6.01512,7.01512];
b=[0.7578,0.2422];
bb=[34.96885,36.96885];
e=prod(cartprod(a,a,b,b),2);
ee=sum(cartprod(aa,aa,bb,bb),2);
n=1;

g=histc([ee e],linspace(min(ee),max(ee),n*(max(ee)-min(ee)+1)),2);
h=linspace(min(ee),max(ee),n*(max(ee)-min(ee)+1));
distr=e'*g;
gauss= @(x,c,r,s) r.*1./(s.*sqrt(2*pi)).*exp(-0.5*((x-c)./s).^2);
k=100.*distr/max(distr);

npts=1000;
resolution=0.25;

x=linspace(min(ee)-1,max(ee)+1,npts);
l=cumsum(gauss(x,h',k',resolution));
l=100*l./max(l(rows(l),:));
plot(x,l(rows(l),:))

which gives:

Compare with Matt Monroe's calculator: