no display name and $DISPLAY environment variable

I've been spending a lot of time organizing my code for the cross correlation project to run more efficiently. Part of this is setting up my code so that it can process/run a lot of different runs simultaneously. It's been a lot of fun learning all the cool tricks I can do with python to help with this. I should blog about some of these tricks.... but not right now....

This blog is about an error I was getting when I tried to run my newly optimized code on riemann. I use qsub/PBS to submit jobs, and I was getting the following errors in my qsub.out file:

File "/clusterfs/riemann/software/matplotlib/1.0.0/lib64/python2.7/site-packages/matplotlib/pyplot.py", line 270, in figure **kwargs)
File "/clusterfs/riemann/software/matplotlib/1.0.0/lib64/python2.7/site-packages/matplotlib/backends/backend_tkagg.py", line 83, in new_figure_manager
window = Tk.Tk() File "/clusterfs/riemann/software/Python/2.7.1/lib/python2.7/lib-tk/Tkinter.py", line 1685, in __init__
self.tk = _tkinter.create(screenName, baseName, className, interactive, wantobjects, useTk, sync, use)
_tkinter.TclError: no display name and no $DISPLAY environment variable

I know that when I run jobs on other riemann nodes using qsub, that the X11 forwarding doesn't work properly. Or in other words, I can't view any of the plots I make when I am running an interactive PBS session. However, for this code, I simply wanted to save plots to .pdf files. However, for some reason there isn't a simple command (that I can find) to plot something in python directly to a file. If anyone out there knows how to do this, please help me!

But I did come across the following from matplotlib's faq page:
~~~~~~~~~~~~~~~~~
How To Generate Python/Matplotlib Images Without Having A Window Appear

The easiest way to do this is use a non-interactive backend such as Agg (for PNGs), PDF, SVG or PS. In your figure-generating script, just call the matplotlib.use() directive before importing pylab or pyplot:

import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt
plt.plot([1,2,3])
plt.savefig('myfig')
~~~~~~~~~~~~~~~~~
You need to make sure to call matplotlib.use() before pyplot or you will get the following error

UserWarning: This call to matplotlib.use() has no effect because the the backend has already been chosen;
matplotlib.use() must be called *before* pylab, matplotlib.pyplot,
or matplotlib.backends is imported for the first time.
if warn: warnings.warn(_use_error_msg)

I actually had an alias set up so that I call pyplot when I open up python, so I was getting this error until I changed that.

So now, I don't the the "$DISPLAY environment variable" error, because python knows to use the non-interactive backend Agg instead of trying to display to the screen. Yippee!

qsub deleting

This command will come in useful again I am sure. It deletes all the runs by a single user in qsub:

qdel $(qselect -u username)

or more specifically for me on riemann:
qdel $(qselect -u jessica)

Fast Queue

As suggested by Ben Weaver I have modified my qsub scripts such that the likelihood jobs are submitted to the fast queue:

#PBS -q fast

is added to .../likelihood/likeqsub.script

This should make my scripts go to a higher priority because they are running for under an hour... however so far they still aren't running at all.... I also am not 100% sure that they run in under an hour... so we'll see.

Miscellaneous IDL Commands

qselect -u jessica -s Q | xargs qdel

window,xsize=700,ysize=600

colorobject = qsos.z_tot
xobject = ugcoaddcolor
yobject = grcoaddcolor

thiscolorbin = 0.2
thiscolorstart = 2.0

xtitle = 'u-g magnitude'
ytitle = 'g-r magnitude'
mtitle = 'Color-Color Diagram All Input QSOs'

tempplot, colorobject, xobject, yobject, thiscolorbin, thiscolorstart, xtitle, ytitle, mtitle

data = qsos.z_tot
datamin = 0.0
datamax = 5.0
binsize = 0.1

xtitle = 'Redshift Bin'
ytitle = 'QSOs in Bin'
mtitle = 'Redshift Distribution'

.com histplot.pro
histplot, data, datamin, datamax, binsize, xtitle, ytitle, mtitle

Parallelizing Likelihood

I spent most of yesterday and today parallelizing the likelihood to make it compute faster. This is how it works:

1. Start with a list of potential targets that you want to compute the likelihoods on.
2. Split up the targets into junks of 1000 objects using the function splitTargets (in the likelihood directory).
3. This will create separate fits files for every 1000 objects in your set, and save them to the current directory with the format: targetfile#.fits, where # will range from 0 to (number of targets / 1000).
4. Go to the directory with the target files and run the following script: likelihood.script (also in the likelihood directory)
5. The likelihood script loops through all files in the current directory of the format targetfile*.fits and via qsub runs the likelihood on each fits file and outputs the result into a file likefile#.fits.
6. This took about 2.5 hours to run on 775,000 targets and just depends on how many jobs are in the queue. Running the likelihood on each 1000 object file takes about 10 minutes.
   
7. Then use the mergelikelihoods function (in the likelihood directory) to read all the likefile#.fits and merge them into one likelihood file.
8. I then write this file to the file likelihoods.fits, which has a 1-1 mapping with the pre-split targetfile: targetfile.fits
   

The log file which runs this whole code is at the following location: ../logs/100411log.pro

Changing 3D Correlation Function

One of my "things to think about" from yesterday's blog post was to make changes to the 3D correlation function so that it was truly an auto-correlation (only having one set of randoms which should speed up the calculation by a bit).

I made changes to my 3D correlation function to do this (at least I think I have) and I am running a test to see if I get the same solution as when I have two separate random catalogs (what I was doing) and also to see how much this speeds things up.

This set of runs in in the working directory: ../pythonjess/run201048_1527 and the log file where the code to run this is ../logs/100408log.py.

The new version runs in 7 seconds, whereas before it took 20 seconds! So that is a good improvement.

Here is a plot to show they look the same

Buzz Comments from Erin about qsub and splitting up the jobs:

Erin Sheldon - Have you used the PBS system on riemann? You can submit jobs and they get sent to different machines on the system. You could split your random runs into 25 different jobs, run them in parallel on 25 different machines, and then just add up the counts later.

That's a factor of n_machines speedup, better than you could every hope to get by just writing faster code.

Below is a sample script, let's call it test.pbs. You can submit with "qsub test.pbs". Of course change all the log file paths to make sense for you.
Erin

#PBS -l nodes=1:ppn=1
#PBS -l walltime=24:00:00
#PBS -q batch
#PBS -N this_job_name
#PBS -j oe
#PBS -o /path/to/where/the/pbs/script/is/test.pbs.pbslog
#PBS -m a
#PBS -V
#PBS -r n
#PBS -W umask=0022

# put setups here

setup numpy -r ~esheldon/exports/numpy-trunk
setup scipy -r ~esheldon/exports/scipy-trunk
setup esutil -r ~esheldon/exports/esutil-trunk
setup python -r ~esheldon/exports/python/2.6.4

# this log is different than above, you can
# watch it in real time
logfile=/path/to/where/the/pbs/script/is/test.pbs.log

# send commands to python
python &> "$logfile" <6:49 amDeleteUndo deleteReport spamNot spam

Erin Sheldon - oh, and to see the status of your jobs
qstat -n -1 -u your_user_name

or equivalently:
~esheldon/python/bin/pbswatch6:55 am

Housekeeping

I've restructured the way the code is run to that I can do multiple runs at once and to hopefully make it easier to re-produce runs and look at older runs quickly. Below is how it works.

Whenever I do a new run, the first thing that python does is create a new 'working directory' of the following format:

runYYYYMD_HM

where
YYYY is the year (2010)
M is the month (2)
D is the day (1)
H is the hour (1337)
so in this example the directory would be called: run201021_1337

This in done in python by the following code:

now = dt.datetime.now()
year = now.year
month = now.month
day = now.day
hour = now.hour
minute = now.minute
workingDir = "run%d%d%d_%d%d"%(year,month,day,hour,minute)
command = "mkdir "+ workingDir
os.system(command)

Within this directory all the photometric and spectroscopic data files are placed for each of the redshift bins. These files are in the format:

photo2D.dat - ra and dec of all photometric data
spec2D_#.dat - ra and dec of the spectroscopic data in the # redshift bin
spec3D_#.dat - x, y, z of the spectroscopic data in the # redshift bin
photoCD.dat - the comoving distance of the photometric data (from photo-zs)

I then compute the 2D cross correlation functions between the photo2D.dat and spec2D_#.dat files and these are saved in the wps#.dat files. The input arguments for these correlation functions are in the files wpsInputs#.

I then compute the 3D auto correlation functions on the spec3D_#.dat files and these are saved in the xiss#.dat files. The input arguments for these correlation functions are in the files xiInputs#.

Because the correlation functions in each redshift bin are submitted as qsub jobs, they take a while to run... so at this point you need to wait for all those jobs to finish. The jobs are named as follows:

2D Cross Correlation: JessWps#workingDir
3D Auto Correlation: JessXiss#workingDir

where # is the redshift bin and workingDir is the name of the working directory.

I can monitor the qsub jobs by typing in 'qstat' to see what is running. Once all the jobs for my working directory are completed then I do the following:

load in the constants from the run file (I put this in the working directory) and it is named YYMMDD_#run.py

where YY is the year (10)
MM is the month (02)
and DD is the day (01)
# is the run number (if I've done multiple runs on this day)

I set the variable workingDir to the working directory.

I load in the photometric comoving distances (which are saved in a file called photoCD.dat):

photoCD = readPhotoCDfile(photoCDfile)

I create the wps and xi matrixes (which are used in the reconstruction) and save them as files:

wpsMat.dat and xiMat.dat respectively.

Then I can do the reconstruction of the redshift distribution and compare it to the photoCD.

Parallelization

Adam rocks. He helped me figure out a way to parallelize the calculations of my correlation functions to that they all happen at once. Hopefully this will speed up my code by a factor of 10-20. Here is how it works:

A while ago I posted about how to use qsub to run scripts. I've now change my python code to call qsub within the program which calculates the correlation functions.

I have the following script:

runcrossCorr.script

~~~~~~~~~~~~~~

#PBS -l nodes=1

#PBS -N JessCC$job

#PBS -j oe

#PBS -m ea

#PBS -M jakirkpatrick@lbl.gov

#PBS -o qsub.out

cd /home/jessica/repository/ccpzcalib/Jessica/pythonJess

echo "$corr2dfile $argumentfile > $outfile"

$corr2dfile $argumentfile > $outfile

~~~~~~~~~~~~~~

this script is called using the following command:

qsub -v job=0, corr2dfile=./2dcorrFinal, argumentfile=./run2010129_237/wpsInputs0, outfile=./run2010129_237/wps0.dat runcrossCorr.script

where -v basically sets up environment variables for the qsub command so that within the script (in this example)

$job is replaced with 0

$corr2dfile is replaced with ./2dcorrFinal

$argumentfile is replaced with /run2010129_237/wpsInputs0

$outfile is replaced with ./run2010129_237/wps0.dat runcrossCorr.script

Many of these qsub commands can be run in python and they all are sent to the queue and run in parallel on different cluster computers.

I've done this so that python automatically changes the argumentfile and outfile and runs all the correlation functions at once.

These correspond to updated versions of my python functions runcrossCorr and runautoCorr in the correlationData library.

Look at all my jobs running (type qstat into terminal)!

And it even emails me when it's done!

So, now all I do is make a python script that runs all the code and set it going in the background:

python 100128run.py 100128.out 2> 100128.err &

I'm running with 1,000,000 points and finer redshift binning. Hopefully this will work even better!

R Packages on Riemann

I asked Michael Jennings to install R on riemann so I can play with this SVM-light thing for the target selection. He installed R but asked me to install any add-on packages into my home directory. Here is how you do this:

1. Create a directory where you want the R add-on packages to be stored (i.e. /home/jessica/R/library)
2. Tell R to look for libraries in this location by adding this to your .path file: export R_LIBS=/home/jessica/R/library
3. Source your .bashrc file to update path changes
   
4. Open R
5. Type the command: install.packages("e1071", lib = '/home/jessica/R/library', repos = "http://cran.r-project.org") where "e1071" is the name of the package you want to install, lib is the directory you created to store the packages, and repos is the location to download the packages from.
   
6. Test to make sure it was successful by closing R, re-opening and importing package you just installed: library("e1071")
   

And now you have successfully installed a R package on Riemann. Thanks to Michael for working so quickly. Now, to try this SVM thing on the quasar targeting...

Here's the code (also in like:svmcode.r)

library("e1071")
library("FITSio")
setwd("/Users/jessica/Documents/workspace/newman/Jessica/likelihood")
qso = readFITS("smallqso.fits")
star = readFITS("smallstar.fits")
x2 <- as.matrix(star$col[[1]],ncol=5)
x1<- as.matrix(qso$col[[1]],ncol=5)
x <- rbind(x1,x2)
y1= seq(1,smallsize)*0+1
y2= seq(1,smallsize)*0
y=c(y1,y2)
fy = factor(y)
m = svm(x,fy)
save(list = ls(), file = "./svmcode.Rdata")

Where smallqso.fits and smallstar.fits are subsets of the qso/star templates used in the likelihood calculation. In this run I am using 100,001 objects for these small files, but that run is taking 2+ days so far on my laptop, so I think I'm going to try to do it with less points for this run on riemann (30,000). To load this session image, type: load("./svmcode.Rdata")

To run this file as a script type the following (for some reason the qsub way of running a script isn't working):
R CMD BATCH ./svmcode.r svmcode.out 2> svmcode.err &

Running Scripts (Update)

I want to utilize the fact that I can run the correlation function on multiple nodes of riemann. Eric helped me learn to use qsub to run a script on multiple nodes. I've got it working as follows:

In the file on riemann 091118jess.script there is the following:

#PBS -l nodes=1
#PBS -N jessjob1
#PBS -j oe
#PBS -M kirkpatrick@berkeley.edu
#PBS -V
cd /home/jessica/repository/ccpzcalib/Jessica/pythonJess
./runcode.script


where runcode.script contains the following:
python runReconstructionMock.py

I make both runcode.script and 091118jess.script executables, and then type in the terminal:
qsub 091118jess.script

I then exit from riemann/insure (didn't close it using x button, not sure if this makes a difference).

When I log back in the code is still running. (check using qstat)

Now I need to integrate this into the python code so that qsubs are generated for each correlation function for each redshift slice. This should make things go a lot faster.

The email from Eric:
----------
From: Eric Huff
Date: Wed, Nov 18, 2009 at 3:06 PM
Subject: qsub
To: Jessica Kirkpatrick

Here's a script example:

#PBS -l nodes=1
#PBS -l walltime=24:00:00
#PBS -m bea
#PBS -M emhuff@berkeley.edu
#PBS -V
ulimit -s unlimited
export PHOTO_DATA=/clusterfs/riemann/raid006/bosswork/groups/boss/photo/data
export PHOTO_REDUX=/clusterfs/riemann/raid006/bosswork/groups/boss/photo/redux
cd /home/emhuff/Coadds/Code/Coadd_Data_Root/Scripts
echo `hostname`
/home/emhuff/Coadds/Code/Coadd_Data_Root/Scripts/r.out 380720 6 0 00

And here's some online documentation for qsub:

http://www.doesciencegrid.org/public/pbs/qsub.html
http://njal.physics.sunysb.edu/tutorial-1.html

I'm running the full reconstruction on the mock data, I might have a fully working code tomorrow. Woo hoo!