Ran's Thoughts

The following post will describe how to install EQMet INTERNET-READY STRONG MOTION ACCELEROGRAPH (TSA-SMA) and set it up for SeisComp3.

Following the instruction in the manual:
The initial steps are to connect the TSA-SMA to electricity, GPS and a computer via the USB device connection.
Then, connect to the device via serial connection.
On windows: read the manual and follow the orders using putty. Don’t forget to download the driver for the device from EQMet site (http://www.eqmet.com/eng/Linux-cdc-acm.zip)

On Linux:
install cu:
# sudo yum install uucp
get the tty of the connection:
# dmesg | grep 'tty'
on my system it was: ttyACM0 but it can be ttyS[1234] or serial.
next I connect to the device:
# cu -l /dev/ttyACM0 -s 115200
it takes a while to prompt for user and password (root,kmi) press enter if does not prompt.
I need to assign an IP based on the mac address so I check the mac:
# ifconfig
after I got the IP, it is easy to connect via http sftp and ssh.

for static IP, edit /etc/network/interfaces:
replace:
iface eth0 inet dhcp
with:
iface eth0 inet static
address <IP>
netmask <MASK>
gateway <GATEWAY>
dns-nameservers <DNS-SERVER>

connect to the device via http.

Next go to Parameters:hardware menu option.

update to your preferable network-station-location-channel ID in:

- Network ID

- Unit ID and Site ID

- dig1 ChX ID

- dig1 ChX Location code

where X in channel 1, 2 or 3.

press OK.

Next, following the response mail I got for EQmet support regarding seedlink connection:

Sir,

The limited version of Rockhound in the TSA-SMA does not support continuous recording. However, you can approximate continuous recording by adding the SEEDLink mirror module and configuring it to accept data from the RockToSlink module at address 127.0.0.1, port 18000 (or as configured). You should add this module after installing and adding in the RockToSlink module.

This module writes the received MiniSEED packets out to storage and breaks the files at approximately the specified duration (depending on MiniSEED packet contents). If you add this module, you should also consider adding an auto file delete module to avoid filling the storage. A file sender module may also be useful (SSH, FTP, or email).

The RockToSlink module for Rockhound is provided by ISTI, New York. The module and documentation is available from ISTI at http://www.isti2.com/RockToSLink/

Technical support for RockToSLink is available from ISTI at isti-info@isti.com. SEEDlink is a public domain program and protocols. Support for these programs including detailed protocol descriptions or installation and setup of the programs is not provided by Kinemetrics.

Regards,
EQMet, Technical Support
support@eqmet.com

installe the RockToSlink module and add it to the workflow:

1. Download the rocktoslik module from http://www.isti2.com/RockToSLink/dist_pub/:
> wget -v http://www.isti2.com/RockToSLink/dist_pub/RockToSLink_v1.0_dist.zip
2. Unzip it:
> unzip RockToSLink_v1.0_dist.zip
following the instractions in http://www.isti2.com/RockToSLink/current_dist/docs/RockToSLink_ovr.html:
3. sftp to the device:
> sftp root@<IP>
4. Change directory to module location:
> cd /usr/rock/SMARTSDist/injar
5 Upload the module-update jar file:
> put RockToSLink_update.jar
6. connect via ssh:
> ssh root@<IP>
7. Change directory to module location:
> cd /usr/rock/SMARTSDist/injar
8. Rename "RockToSLink_update.jar" to "update.jar"
> mv RockToSLink_update.jar update.jar
9. Reboot the instrument:
> reboot
Now add the RockToLink module to device layout:
10. In the device http menue Parameters: Advanced features, add options for Add/Remove Modules, Replace Modules and Advanced Modules.
11. In Module Add, add the RockToLink module.
12. Press apply changes now.

Now go to system operation and make sure all settings are correct.

enjoy.

Download and install Machinoid following: http://www.machinoid.com/?p=10
apt-get install libtool gfortran
apt-get install libusb-1.0
install xenomai following the instructions in the link above.
sysctl vm.min_free_kbytes=8192
vi /etc/ld.so.conf.d/xenomai.conf
add 2 lines:
# xenomai default configuration
/usr/xenomai/lib
ldconfig
install wifi drivers (after I got problems with rt2800usb driver):
apt-get install firmware-ralink

create new file /etc/modprobe.d/8192cu.conf
add:
# Disable power saving
options 8192cu rtw_power_mgnt=0 rtw_enusbss=1 rtw_ips_mode=1

in /etc/dhcp/dhcp.conf

uncomment the timeout line and edit:

timeout 100;

After getting my new LSM303 accelerometer, I got my great electronics tec. Daviv Shtibelman to wire it up for RPI.
Connecting LSM303 pins 1 (SCL), 2 (SDA), 6 (GND) and 8 (V3V) to pins 5 (SCL), 3 (SDA), 6 (GND) and 1 (V3V) on the Pi board.
On my ArchlinuxARM OS I had to install i2c-tools:
> pacman -S i2c-tools
And add i2c modules:
> echo i2c-dev > /etc/modules-load.d/ i2c.conf
> echo i2c-bcm2708 >> /etc/modules-load.d/ i2c.conf
> reboot

now to check it works:
> i2cdetect -y 1
got a message like this:
0 1 2 3 4 5 6 7 8 9 a b c d e f
00: -- -- -- -- -- -- -- -- -- -- -- -- --
10: -- -- -- -- -- -- -- -- -- 19 -- -- -- -- 1e --
20: -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
30: -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
40: -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
50: -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
60: -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
70: -- -- -- -- -- -- -- --

Now, following this, I test it and it works. Now I'm up for creating my own code.

Made some changes to the code:

main.cpp

#include"LSM303DLHC.h"
#include<stdio.h>
#include <unistd.h>

/*
getHeading() calculates a tilt-compensated heading.
A float between 0 and 360 degrees is returned. You need
to pass this function both a magneto and acceleration array.

Headings are calculated as specified in AN3192:
http://www.sparkfun.com/datasheets/Sensors/Magneto/Tilt%20Compensated%20Compass.pdf

*/

float getHeading(LSM303& lsm303dlhc);

#define PI 3.141592654

int main(void)
{
uint8_t bajt;
const char *fileN = "/dev/i2c-1";
LSM303 lsm303dlhc(fileN);

lsm303dlhc.enalbe();

while(1)
{
lsm303dlhc.readAccelerationRaw();
lsm303dlhc.readMagnetometerRaw();
printf("acc [m/s^2]: \e[27;1;31m %f \e[m \e[27;1;32m %f \e[m \e[27;1;34m %f \e[m mag: \e[27;1;31m %d \e[m \e[27;1;32m %d \e[m \e[27;1;34m %d\e[m %fdeg\n",
(int16_t)lsm303dlhc.acc_x_raw*0.00957,(int16_t)lsm303dlhc.acc_y_raw*0.00957,-(int16_t)lsm303dlhc.acc_z_raw*0.00957,
(int16_t)lsm303dlhc.mag_x_raw, (int16_t)lsm303dlhc.mag_y_raw, (int16_t)lsm303dlhc.mag_z_raw,getHeading(lsm303dlhc));
usleep(10);
}

}

float getHeading(LSM303& lsm303dlhc)
{
float heading,pitch,roll,xh,yh,zh;
// see section 1.2 in app note AN3192
int magValue[3];
float accelValue[3];
magValue[0] = (int16_t)lsm303dlhc.mag_x_raw;
magValue[1] = (int16_t)lsm303dlhc.mag_y_raw;
magValue[2] = (int16_t)lsm303dlhc.mag_z_raw;
accelValue[0] = (int16_t)lsm303dlhc.acc_x_raw*0.000976531;
accelValue[1] = (int16_t)lsm303dlhc.acc_y_raw*0.000976531;
accelValue[2] = -(int16_t)lsm303dlhc.acc_z_raw*0.000976531;

// see appendix A in app note AN3192
pitch = asin(-accelValue[0]);
roll = asin(accelValue[1]/cos(pitch));

xh = magValue[0] * cos(pitch) + magValue[2] * sin(pitch);
yh = magValue[0] * sin(roll) * sin(pitch) + magValue[1] * cos(roll) - magValue[2] * sin(roll) * cos(pitch);
zh = -magValue[0] * cos(roll) * sin(pitch) + magValue[1] * sin(roll) + magValue[2] * cos(roll) * cos(pitch);

heading = 180*atan2(yh,xh)/PI;

if (heading <0)
heading += 360;

return heading;
}

LSM303DLHC.cpp:

#include"LSM303DLHC.h"
#include<math.h>
#include<stdio.h>

/*Conection to Raspberry PI:
LSM303 Raspberry PI
VDD -> 3V3(PIN 1)
SDA -> SDA(PIN 3)
SCL -> SCL(PIN 5)
GND -> GND(PIN 6)
*/

#define LSM303DLHC_MAG_ADDRESS (0x3C >> 1)
#define LSM303DLHC_ACC_ADDRESS (0x32 >> 1)

LSM303::LSM303(const char * i2cDeviceName) : i2c_lsm303(i2cDeviceName)
{

}

uint8_t LSM303::readAccRegister(uint8_t regAddr)
{
i2c_lsm303.addrSet(LSM303DLHC_ACC_ADDRESS);
return i2c_lsm303.readByte(regAddr);
}

uint8_t LSM303::readMagRegister(uint8_t regAddr)
{
i2c_lsm303.addrSet(LSM303DLHC_MAG_ADDRESS);
return i2c_lsm303.readByte(regAddr);
}

void LSM303::writeAccRegister(uint8_t regAddr,uint8_t byte)
{
i2c_lsm303.addrSet(LSM303DLHC_ACC_ADDRESS);
i2c_lsm303.writeByte(regAddr, byte);

}

void LSM303::writeMagRegister(uint8_t regAddr, uint8_t byte)
{
i2c_lsm303.addrSet(LSM303DLHC_MAG_ADDRESS);
i2c_lsm303.writeByte(regAddr, byte);

}

void LSM303::enalbe(void)
{
writeAccRegister(LSM303_CTRL_REG1, 0b10010111);
writeAccRegister(LSM303_CTRL_REG4, 0b00001000);

writeMagRegister(LSM303_MR_REG, 0x00);
}

void LSM303::readAccelerationRaw(void)
{
uint8_t block[6];

i2c_lsm303.addrSet(LSM303DLHC_ACC_ADDRESS);

i2c_lsm303.readBlock(0x80 | LSM303_OUT_X_L_A, sizeof(block), block);
acc_x_raw = (int16_t)(block[0] | (block[1] << 8)) >> 4;
acc_y_raw = (int16_t)(block[2] | block[3] << 8) >> 4;
acc_z_raw = (int16_t)(block[4] | block[5] << 8) >> 4;

}

void LSM303::readMagnetometerRaw(void)
{
uint8_t block[6];

i2c_lsm303.addrSet(LSM303DLHC_MAG_ADDRESS);
i2c_lsm303.readBlock(0x80 | LSM303_OUT_X_H_M, sizeof(block), block);

mag_x_raw = (int16_t)(block[1] | block[0] << 8);
mag_y_raw = (int16_t)(block[5] | block[4] << 8);
mag_z_raw = (int16_t)(block[3] | block[2] << 8);

}

void LSM303::readAcceleration(void)
{
readAccelerationRaw();
}

After installing and running ArchLinuxARM on my RaspberryPi, I tried to compile a real-time kernel.
With the help of moonman, following this ArchLinuxARM forum thread, I downloaded the RT patch:
https://www.kernel.org/pub/linux/kernel/projects/rt/3.6/patch-3.6.11.4-rt36.patch.gz
and the PKGBUILD:
# git clone https://github.com/archlinuxarm/PKGBUILDs.git
following this, I updated pacman:
# pacman -Syy
and installed some packages:
# pacman -S kernel26-headers file base-devel abs
next
# cd PKGBUILDs/core/linux-raspberrypi
and edit the PKGBUILD file to include the patch in the source array:

source=('config'
'change-default-console-loglevel.patch'
'usb-add-reset-resume-quirk-for-several-webcams.patch'
'args-uncompressed.txt'
'boot-uncompressed.txt'
'imagetool-uncompressed.py'
'https://raspy-juice.googlecode.com/svn/trunk/linux-rtc/0001-rtc-pcf8523.patch'
'https://raspy-juice.googlecode.com/svn/trunk/linux-rtc/0002-pcf8523-i2c-register-dt.patch'
'patch-3.6.11.4-rt36.patch')

add a patch line:

patch -Np1 -i "${srcdir}/patch-3.6.11.4-rt36.patch"

and add md5sums at the end of the PKGBUILD file. to get the correct number run:

# md5sum patch-3.6.11.4-rt36.patch

next for the build:

# makepkg --asroot -Acs

The process has failed with some rejections of patches. So we need to manually apply the patches:

first look for rejections:

# find . -name "*.rej"

then edit the file and the .rej file:

# vim XXX.c XXX.c.rej

replacing XXX with the file path and name.

Then use :split to split the screen :next to change one screen to .rej file and ctrl+ww to switch between files.

add the lines marked by + and remove lines marked by - (lines with no +- sign are used for orientation.

after manually patching all files, edit the PKGBUILD file to ignore the kernel git and the patches (adding # at the start of the line)
add a line to src/linux/security/apparmor/sid.c:
"#include <linux/cache.h>"

re-run:

# makepkg --asroot

Now for the installation part:
# pacman -U linux-headers-raspberrypi-3.6.11-11-armv6h.pkg.tar.xz linux-raspberrypi-3.6.11-11-armv6h.pkg.tar.xz

And reboot.

Well, it didn't work.

I initiated a new project involving low-cost micro-computer RaspberryPi.
In order to begin I ordered some basic parts (costs are in USD for Israeli market, should be cheaper elsewhere):
a RaspberryPi (55$)
5v charger with microusb adapter (10$)
8GB SD card (12$)
HDMI cable (8$)

and some extra parts:
3.5" LCD screen (45$)
micro-keyboard (62$)
micro-sd low profile adapter (11.9$)
wifi dongle (18.6$)

Not so cheap but might be worth it.

after I got a couple of RaspberryPi (model B rev2 512MB RAM) devices for a project I'm working on (details about the project will come later if it will work).
I bought a SD 8GB card and hooked it to the my desktop computer.
I downloaded the ARCLinux image unzipped it and followed the instructions:

Replacing sdX with the location of the SD card, run:
dd bs=1M if=/path/to/archlinux-hf-2013-06-15.img of=/dev/sdX

In order to increase the SD space to it's full 8GB I used gparted.

upgrade system:
> pacman -Sy pacman
> pacman-key --init
> pacman -Syu
sync (3 times):
# sync
reboot:
# reboot
Now I can add packages using pacman. to look for a package use:
> pacman -Ss <package name>
then select the one you need and install by:
> pacman -S <package name>

I am currently testing the Pi with the screen and 8 AA batteries to see how long it can work. I plan to build a case with a built-in 3.5" screen and a miniature keyboard so I can have a control unit for other RaspberryPi devices. My initial thoughts of using an old Gamboy case are now switched (with the encouragement of my electronic tec. guy David Shtibelman) to build a costume made case. Of course, David is the one to build it, since I have no skills to do that.

These are MY needed packages, feel free to follow my steps or ignore this.
1. Add rpmfusion to yum repo list:
> sudo yum localinstall --nogpgcheck http://download1.rpmfusion.org/free/fedora/rpmfusion-free-release-stable.noarch.rpm http://download1.rpmfusion.org/nonfree/fedora/rpmfusion-nonfree-release-stable.noarch.rpm
2.
sudo yum install -y yumex tcsh gmt ImageMagick ipython scipy scitools spyder gthumb nautilus-open-terminal bluecurve-icon-theme gnome-applets gnome-themes bluecurve-gnome-theme filezilla gimp swig lshw lshw-gui gnome-tweak-tool libreoffice gv gcc gcc-gfortran gcc-c++ dconf-editor tigervnc-server
install chrome from google
install virtualbox from oracle site.
copy panels:
on old machine: dconf dump /org/gnome/gnome-panel/ > panels
on new machine: dconf load /org/gnome/gnome-panel/ < panels
change with dconf-editor:
/org/gnome/desktop/interface/icon-theme to Bluecurve
/org/gnome/desktop/wm/preferences/theme Clearlooks
/org/gnome/desktop/background/show-desktop-icons True
change user shell to tcsh: sudo usermod -s /bin/tcsh USERNAME
install microsoft fonts:
wget http://www.my-guides.net/en/images/stories/fedora12/msttcore-fonts-2.0-3.noarch.rpm
sudo rpm -ivh msttcore-fonts-2.0-3.noarch.rpm

6.Install obspy as root:
sudo easy_install -N obspy.core obspy.mseed obspy.sac obspy.gse2 obspy.imaging obspy.signal obspy.arclink obspy.datamark obspy.db obspy.earthworm obspy.fissures obspy.iris obspy.neries obspy.realtime obspy.segy obspy.seg2 obspy.seisan obspy.seishub obspy.sh obspy.taup obspy.wav obspy.xseed

Finally, I finished my PhD. and now it time to move forward.
Taking my Professors advice, I'm making the first steps in the world of seismology and I'm going to spend the next few years in Tel-Aviv University, working with Dr. Alon Ziv and Dr. Hillel Wust-Bloch.
The main work is building and maintain a local mini seismic-array network. I am going to actually place the seismometers, connect the dataloggers and build the server to collect the data. We are going to use different devices - Lennartz LE-3Dlite1DV MkII seismometers, Trillium broadband seismometers, TSA-SMA accelerometers and Summit M Hydra dataloggers. On the server side, We are going to use Seiscomp3 as our main framework, but we probably research for new algorithms and Seiscomp3 modules to better suit it to our needs.

Me getting the PhD. diploma from Prof. Rivka Carmi, presedent of BGU

The JW24F14 device is accessible as a HID only as root. In order to access it as a user, I had to write a UDEV rule with the help of this post. more info about udev rules can be found here.
issue the command:
> sudo vi /etc/udev/rules.d/10-JoyWarrior24F14.rules
and when in editing mode (i) add:

# Rules for accessing JoyWarrior24F14 HID as user:
SUBSYSTEM=="usb", ATTRS{idVendor}=="07C0", ATTRS{idProduct}=="0140", MODE:="666", GROUP="users"
SUBSYSTEM=="usb", ATTRS{idVendor}=="07c0", ATTRS{idProduct}=="1116", MODE:="666", GROUP="users"

I wanted to have a file change listener such that if new data is added to a file I could run a code.
Following this post, I've written a modification so now the listener is working in multiprocessing.

#!/usr/bin/python
# A code to stream file as it grows
# by Ran Novitsky Nof Oct 25, 2012

import os,time
from multiprocessing import Process,Pipe

# change this file name
FileName = 'stream.txt'

# adapted from http://www.valuedlessons.com/2008/04/events-in-python.html

# the Event class take care of adding, removing and executing the functions
class Event:
def __init__(self):
self.handlers = set()

def handle(self, handler):
self.handlers.add(handler)
return self

def unhandle(self, handler):
try:
self.handlers.remove(handler)
except:
raise ValueError("Handler is not handling this event, so cannot unhandle it.")
return self

def fire(self, *args, **kargs):
for handler in self.handlers:
handler(*args, **kargs)

def getHandlerCount(self):
return len(self.handlers)

__iadd__ = handle
__isub__ = unhandle
__call__ = fire
__len__ = getHandlerCount

# the stream class takes care of the
# function to fire in case of an event
class Stream:
def __init__(self,fname):
self.fname = fname
self.fd = open(fname,'r')
self.get_stream()
def get_stream(self):
# this is the actual function to fire
print self.fd.read(),
def restart(self):
self.fd.seek(0)

# the MockFileWatcher creats the multiprocessing deamon to alert on file change
class MockFileWatcher:
def __init__(self,source_path,sleeptime=0.1):
self.fileChanged = Event()
self.source_path = source_path
self.running=False
self.deamon=False
self.fd = None
self.sleep = sleeptime
# the next two lines can be set outside the class
self.stream = Stream(source_path)
self.fileChanged += self.stream.get_stream # add event handler
def watchFiles(self):
# open the file in non-blocking read mode
self.fd = os.open(self.source_path,os.O_RDONLY | os.O_NONBLOCK)
# go to the end
os.lseek(self.fd,0,2)
self.running=True
self.deamon=False
while self.running:
# try to read
if os.read(self.fd,1)!='':
# New data was added
os.lseek(self.fd,0,2) # go to the end of file
self.fileChanged() # fire the event
# see if the parent process has a message
if self.chiled_conn.poll():
# execute the message
eval(self.chiled_conn.recv())
# sleep for a while
time.sleep(self.sleep)
def watchFilesDeamon(self):
# start the deamon
self.parent_conn,self.chiled_conn = Pipe()
p = Process(target=self.watchFiles)
p.start()
self.p = p
self.deamon=True
return p
def stop(self):
# stop deamon/process
if self.deamon:
self.parent_conn.send('self.stop()')
print "Hold while stopping child process..."
# hold depends on the sleeping time
if self.parent_conn.recv():
pass
#print "chiled is stopped"
self.p.terminate() # terminate multiprocessing
self.deamon=False
# stop process
self.running=False
self.fd=None
# if this is the process make sure parent know we are done
if not self.deamon: self.chiled_conn.send(True)
return True

if __name__=='__main__':
watcher = MockFileWatcher(FileName)
p = watcher.watchFilesDeamon()
print 'from this point the listner is running at the background'
print 'add more code here so you can see how it works.'

Try this by changing the file name at the top and run it. then add new lines to the file.

After installing Seiscomp3 and running it for a while, now it is time to play with the data.
As a first step I want to be able to see the events outside the SeisComp3 environment.
This will later serve me when exporting the data to a web page or reports.
The steps described here are: 1) get a list of events, 2) export the data to an xml file and 3) print the data as a bulletin.
1) Get a list of events:
scevtls -d mysql://sysop:sysop@localhost/seiscomp3
Will give a list of event codes like:

Wrong 'begin' format -> setting to None
Setting start to 1970-01-01 00:00:00
Setting end to 2012-10-24 08:17:13
gfz2012uizn
gfz2012atlf
...

where gfz will be the agency prefix and the last 4 letters will be random code.
You can also set begin and end time to get only a time window events. see scevtls documentation.

2) For this example I select one event (in csh environment):
set EVENT = gfz2012uizn
Now we can export the event data to an xml file:
scxmldump -d mysql://sysop:sysop@localhost/seiscomp3 -E $EVENT -P -M -A -o test.xml
see scxmldump for more information

3) Print the data to stdout:
scbulletin -i test.xml
Note that it is possible to skip the xmp part by using:
scbulletin -d mysql://sysop:sysop@localhost/seiscomp3 -E $EVENT
also adding -3 as a flag will output even more details.
see scbulletin for more info.

I'll update later with a python script to export the data to a Google Earth kml file.

After installing Fedora 17 I will now describe how to install MySQL server and SeisComp3.

first, install MySQL (following this):
> sudo yum install -y mysql mysql-server

enable and start MySQL service:
> sudo systemctl enable mysqld.service
> sudo systemctl start mysqld.service

Next, secure the server:
> sudo /usr/bin/mysql_secure_installation
(answer yes to all questions and change the MySQL root password)
If you need to change the database location from the default follow this (source):
Stop the service:
> systemctl stop mysqld.service
Copy from default to new location:
> cp -R -p /var/lib/mysql /path/to/new/datadir
Clean unneeded files:
> rm /path/to/new/datadir
(don't worry about error messages - the undeleted files are the databases)
Edit the configuration file:
> gedit /etc/mysql/my.cnf
(change keyword value of datadir to the new location)
Restart the service:
> systemctl restart mysqld.service

Now to the SeisComp3 install (Zurich version for fedora16 modified for fedora 17):
Add needed packages:
> sudo yum install -y boost-filesystem boost-iostreams boost-program-options boost-serialization boost-signals
boost-thread boost-wave boost-regex boost-devel ncurses python-numeric numpy alsa-utils festival
Link installed libraries (so.1.48) to the needed packages (so.1.47):

> sudo ln -s /usr/lib64/libboost_iostreams.so.1.48.0 /usr/lib64/libboost_iostreams.so.1.47.0

> sudo ln -s /usr/lib64/libboost_filesystem.so.1.48.0 /usr/lib64/libboost_filesystem.so.1.47.0
> sudo ln -s /usr/lib64/libboost_system.so.1.48.0 /usr/lib64/libboost_system.so.1.47.0
> sudo ln -s /usr/lib64/libboost_regex.so.1.48.0 /usr/lib64/libboost_regex.so.1.47.0
> sudo ln -s /usr/lib64/libboost_thread-mt.so.1.48.0 /usr/lib64/libboost_thread-mt.so.1.47.0
> sudo ln -s /usr/lib64/libboost_program_options.so.1.48.0 /usr/lib64/libboost_program_options.so.1.47.0
> sudo ln -s /usr/lib64/libboost_signals.so.1.48.0 /usr/lib64/libboost_signals.so.1.47.0
Put key files in ~/.seiscomp3/key/ directory
Install using the setup and seiscomp config procedures described in the SeisComp3 documentations
add a source to env.sh file in your .bashrc:
source ~/seiscomp3/lib/env.sh

For the seattle release use:
>sudo yum install -y boost boost-devel python-numeric

Don't forget to download and extract to the seiscomp directory, the correct fedora seiscomp package (e.g. seiscomp3-zurich-xxxx.xxx.xx-fedora16-x86_64.tar.gz), the maps (e.g. seiscomp3-maps.tar.gz) and the methadata of GE network (e.g. seiscomp3-metadata-GE-2010.tar.gz).

how to use python to write to a certain location (x,y) on the stdout (terminal):

copied form here

import sys
def print_there(x, y, text):
sys.stdout.write("\x1b7\x1b[%d;%df%s\x1b8" % (x, y, text))
sys.stdout.flush()

1. format the hard drives to raid using fdisk (gdisk in fedora 17)
or with parted (for 3TB disks):

> parted /dev/sdx
(parted) mklabel gpt
(parted) unit TB
(parted) mkpart primary 0 3001G
(parted) print
(parted) quit

2. Create the raid device (change level and devices as needed):
mdadm --create --level=1 --raid-devices=2 /dev/md0 /dev/sd[bc]1
3. format the file system of the new device:
mkfs.ext4 /dev/md0
4. create the mdadm.conf file:
mdadm --detail --scan > /etc/mdadm.conf
5. add mail alert:
echo "MAILADDR your_mail@somthing.com" >> /etc/mdadm.conf
echo "MAILFROM your_server" >> /etc/mdadm.conf
6. test mail alert:
mdadm --monitor --scan --test --oneshot
7. start mail alert deamon:
systemctl enable mdmonitor.service
systemctl start mdmonitor.service
8. edit the mail aliases file to forward root messages to your email:
vi /etc/aliases
add your mail under the line: # Person who should get root's mail

run newaliases command

see also at:
http://www.iceteks.com/articles.php/linuxmdadmraid/

* after name change, create the /etc/mdadm.conf file and run:
sudo update-initramfs -u
or
sudo dracut -f /boot/initramfs-[current].img

based on my-guides post-install post

1. Add to sudo:
> su --login -c 'visudo'
add under the line:
root ALL=(ALL) ALL
the line:
username ALL=(ALL) ALL
and remove the sign # from the line with %wheel
save and test using the line:
> sudo whoami
(answer should be root)
2. Disable SeLinux and firewall:
> sudo vi /etc/selinux/config
change enabled to disabled
> sudo systemctl disable iptables.service
> sudo systemctl stop iptables.service
3. change autologin for user:
sudo vi /etc/gdm/custom.conf
under [daemon], add:
TimedLoginEnable=true
TimedLogin=yourUserName
TimedLoginDelay=1
4. Switch to gnome3 fallback mode:
> gsettings set org.gnome.desktop.session session-name 'gnome-fallback'
5. enable sshd:
>sudo systemctl enable sshd.service
> sudo systemctl start sshd.service

6. auto rotate screen on login:
add to .profile on home:
xrandr --output VGA-1 --rotate left
7.add daily update report (as root):
following this:
> vi /etc/cron.daily/yum_update_note
#!/bin/sh

if [ ! -f /var/lock/subsys/yum ] ; then
msg=`/usr/bin/yum check-update 2>&1`
i=`echo "$msg" | wc -l`
if ((i>1)); then
echo "$msg" | /bin/mail -s "yum check-update output" root
fi
fi

and chnage mode for executable:

> sudo chmod +x /etc/cron.daily/yum_update_note

See a blog post about my recent paper:
Dead Sea rising: how the lowest place on Earth is getting (slightly) higher - GeoSpace - AGU Blogosphere

In order to run a file whenever a certain mail arrives to my inbox I created a filter using the Thunderbird add-on FiltaQuilla. The add-on expand the messages filtering options and add an option to "run a file". when a certain mail match the filter, Thunderbird will run the file specified. I use a python code - code.py so in the filter configuration under Linux (Fedora core 14) I use this line:

/usr/bin/python,/home/UserName/code.py,-a

The line is a comma separated 3 parts line. The first part is the python interperter, the second is the script name to run and the last part is an argument needed for the script to run. Omitting the first part and the filter wont work.

A second version of ASTER GDEM was released in Oct. 17 2011 by The Ministry of Economy, Trade and Industry of Japan (METI) and the National Aeronautics and Space Administration (NASA). The GDEM is in a 30m pixel resolution and 1x1 degree tiles. The data is distributed as zipped GeoTIFF files. When a larger than 1x1 deg. DEM is needed for Interferometry processing, you can use my GDEM.py script available from my script page. This python script will need a system with pylab and PIL python modules.

Simply:

Download the GDEM zipped files (see here on how to do that).

Unzip the zip files in the working directory.

The working directory should contain sub-directories containing the GeoTIFF files.

Run the script giving -h for help or the final DEM name as arguments.

The Final DEM is aimed for using as DEM in ROI-PAC interferometry processing.

If a -g is added as an argument, the produced DEM will be suitable for interferometric processing using Gamma software.

In order to create a contact-sheet or a collage of many files under linux, I prefer to use a single command-line. I create a directory with all the images and use ImageMagick montage function:

montage -font Arial-Narrow -pointsize 10 -label '%t' -resize 100x150 *.jpg -geometry +3+3 -tile 4x4 -shadow collage.jpg

The Problem:
I have a Linux machine (Fedora 14) with no screen attached. I view the graphics via VNC server. I had to adjust the system to use the Nvidia driver and tweak it for no screen.

The solution:
First we need to install the Nvidia driver. Then, we need to adjust the X configuration for no screen and our preferred resolution.

for more details see: http://forums.fedoraforum.org/showthread.php?t=204752

1. Update the kernel:

su
yum update kernel*
reboot

2. Add rpmfusion to your repository.

su
rpm -Uvh  http://download1.rpmfusion.org/free/fedora/rpmfusion-free-release-stable.noarch.rpm   http://download1.rpmfusion.org/nonfree/fedora/rpmfusion-nonfree-release-stable.noarch.rpm

3. Depending on your card add nvidia akmod to the kernel and the proper driver. for a Gforce 9500GT on my 64bit machine:

yum install akmod-nvidia-260.19.36-1.fc14.x86_64 xorg-x11-drv-nvidia-libs-260.19.36-1.fc14.x86_64

4. Edit the xorg.conf file. don't forget to adjust the Nvidia card for no screen attached by adding:
Option "ConnectedMonitor" "CRT"
to the device section and adjust the screen resolution under the Screen section:

Section "Device"
    Identifier     "Card0"
    Driver         "nvidia"
	Option "ConnectedMonitor" "CRT"
EndSection

Section "Screen"
    Identifier     "Screen0"
    Device         "Card0"
    Monitor        "Monitor0"
    SubSection     "Display"
        Viewport    0 0
        Modes      "1280x1024"
    EndSubSection
EndSection

Using JAXA satellite "DAICHI" (ALOS) for SAR Interferometry might raise the need for calculating perpendicular baselines. One can use AUIG site to locate the needed images for interferometry and use PALSAR baseline calculator CalcBP for calculating perpendicular baselines. more info is available here.

Well, after a long time using Fedora 11 I'm moving on to Fedora 14.
Here I'll detail the steps I used for the installation.
I generally follow my Fedora 10 installation guide with some changes of course.

LiveCD:
I use a liveCD (download page) to install.
after loading the system via cd,
I use the "Install to Hard Drive" icon.
next,next,next.
Giving the Hostname, setting location, giving a good root password.

File System:
Since I use several hard drives I need to configure the file system. I use one 1.5TB disk as my home, another one for data and processing and one for the system so all the installation will only deal with the last one. later I'll mount all the rest of the disks.
so the type of installation for me is "Replace Existing Linux System(s)", but with "Review and modify layout" switch on and I select the system disk as "Install Target Devices".
Next I change the partitions and mount points of all disks as needed.
Finally, The boot loader is set up and the system is installed on the hard drive.

Updates and softwares:
As a first step few softwares should be installed and update (this might take forever):

system-config-network
yumex
gconf-editor
gedit-plugins
gthumb
nautilus0open-terminal
tango-icon-theme
tango-icon-extras
fonts for Hebrew
all stix fonts
cm-lgc fonts
blender
gimp
gv
ImageMagick
inkscape
filezilla
thunderbird
emacs
open office
vim-enhanced
gcc
gmt
grace
HDF
HDF5
lapack
netcdf
numpy
python support for open office
compat-gcc
glade3
ipython
pango-devel
wget
tcsh
system-config-samba
system-config-users
system-config-nfs
dos2unix
selinux configuration gui
matplotlib
scipy
scitools
kdocker
numlockx
bluecurve icons
ssmtp

Customizing:
remove SElinux policy
disable filrewall
configure remote desktop
add user to sudo:

su --login -c 'visudo'
Now go below the line:
root ALL=(ALL) ALL
press a, and type the following
username ALL=(ALL) ALL
where username is the username you use to login. Next press Escape.
Now, if you want to be prompted for your root password each time you use the sudo command go to this line:
# %wheel ALL=(ALL) ALL
and with your cursor pointing on # press x
If you don't want to be prompted for your root password each time you use the sudo command go to this line:
# %wheel ALL=(ALL) NOPASSWD: ALL
and with your cursor pointing on # press x
Next, press :wq to save and exit. (type the : as well)
* Adding your user to the wheel group
Just type:
su -c 'gpasswd -a username wheel'
* Testing sudo
To test if you have done it correctly as a simple user type:
$ sudo whoami
If everything is working properly the command should return the word 'root'.

change autologin for user:
sudo vi /etc/gdm/custom.conf
under [daemon], add:
AutomaticLoginEnable=true
AutomaticLogin=yourUserName

edit /etc/fstab for mounting external HD (LAN):
//XXX.XX.XX.XX/DATA /mountpoint cifs auto,owner,rw,username=XXX,uid=XXX,gid=XXX 0 0
NOTE! - it will not mount on startup. need to do: "sudo mount mountpoint" to mount. one can add password=[user-passwd] after username=[username] or instead username and password : credentials=[path to a file containing user and password]
the file should contain the remote system properties:
username=[username]
password=[userpassword]

Enable Num Lock on GNOME startup:
su -c 'yum install numlockx'
su -c 'cp /etc/gdm/Init/Default /etc/gdm/Init/Default_backup'
su -c 'gedit /etc/gdm/Init/Default'
And add the following lines at the end of the file before the exit 0.
if [ -x /usr/bin/numlockx ]; then
/usr/bin/numlockx on
fi

set the samba server

set the mail server:
sudo vi /etc/ssmtp/ssmtp.conf
change the line:
mailhub=mail
to:
mailhub=[your outgoing mail server]

Set the screen resolution for vnc and no monitor:
create xorg.conf file
sudo Xorg :1 -configure
sudo cp /root/xorg.conf.new /etc/X11/xorg.conf
and add lines to the xorg.conf file:

Section "Monitor"
Identifier "Monitor0"
ModelName "LCD Panel 1280x1024"
HorizSync 31.5 - 64.0
VertRefresh 56.0 - 65.0
Option "dpms"
EndSection

Section "Device"
   ### Available Driver options are:-
   ### Values: <i>: integer, <f>: float, <bool>: "True"/"False",
   ### <string>: "String", <freq>: "<f> Hz/kHz/MHz",
   ### <percent>: "<f>%"
   ### [arg]: arg optional
   #Option "SWcursor" # [<bool>]
   #Option "HWcursor" # [<bool>]
   #Option "NoAccel" # [<bool>]
   #Option "ShadowFB" # [<bool>]
   #Option "UseFBDev" # [<bool>]
   #Option "Rotate" # [<str>]
   #Option "VideoKey" # <i>
   #Option "FlatPanel" # [<bool>]
   #Option "FPDither" # [<bool>]
   #Option "CrtcNumber" # <i>
   #Option "FPScale" # [<bool>]
   #Option "FPTweak" # <i>
   #Option "DualHead" # [<bool>]
Identifier "Card0"
Driver "vesa"
BusID "PCI:1:0:0"
EndSection

Section "Screen"
Identifier "Screen0"
Device "Card0"
Monitor "Monitor0"
SubSection "Display"
Viewport 0 0
Depth 24
Modes "1280x1024"
EndSubSection
EndSection

set the panels and desktop background as you like it.

Install RPMfusion repositories following: www.rpmfusion.org/Configuration

add true type fonts from windows

RoiView was originaly designed to replace DGX/MDX software which require bureaucracy to obtain. This software was used in order to view ROI-PAC interferometry processing results but producing an image from file was not so strait forward. In time support was added for ERDAS ER Mapper header files and for user provided parameters with or without header files.

Current RoiView version (0.75) is available for download at

See also RoiView's help files

Screenshot:

I needed to to get the text of a specific div elemnt in an html file. I tried to use python's standard library modules for markup processing such as htmllib etc. but I couldn't figure how to use them. I've created my own module just for getting the text from an html element:

#!/bin/env python
from htmlentitydefs import entitydefs as ent
import string
# This module enable you to extract text from a certian HTML element
#
# by Ran Novitsky Nof, 2010
# ran.nof@gmail.com
#
# example of use:
# say we want to get the text in an element of type tag (e.g. 'div','a','span' etc.)
# who has an attribute key (e.g. "id","class","href" etc.) with a value of val
# for example in order to extract the text of a div element with id of textdiv from a file htmlfile.html:
#
# <html>
# <head>
#   :
# </head>
# <body>
#   :
# <div "id"="textdiv">I will not buy this <a href="spam">record</a> it is scratched. </div>
#   :
# </body>
# </html>
#
# use:
# from htmlparser import Parser
# htmlfile='htmlfile.html'
# tag,key,val = ('div','id','textdiv')
# text=Parser(htmlfile).getText(tag,key,val)
# print(text)
#
class Element():
def __init__(self):
    self.startTag = -1
    self.endTag = -1
    self.attrib = {}
    self.keys = self.attrib.keys()
    self.innerHTML = ''
    self.tag = ''
    self.start = -1
    self.end = -1

class Parser():
def __init__(self,infile):
    self._root = open(infile).read()
    self._root = self._root[self._root.find('<body'):self._root.find('</body')]
    self._root = self._root[self._root.find('>')+1:]
    self.tags = set()
    i,j=0,0
    self.tagstarts = {}
    self.tagends = {}
    self.elements = []
    while i<len(self._root):
      i = self._root.find('<',i)
      j = self._root.find('>',i)
      if not j>i: break
      tag = self._root[i+1:j].split()[0]
      if tag.startswith('/'):
        tag = tag[1:]
        self.tagends[tag][-1-self.tagends[tag][::-1].index(None)]=((i,j))
      else:
        self.tags.add(tag)
        if tag in self.tagstarts:
          self.tagstarts[tag].append((i,j))
          self.tagends[tag].append(None)
        else:
          self.tagstarts[tag]=[(i,j)]
          self.tagends[tag]=[None]
      i=j+1
    self.getElements()
def getElements(self):
    for tag in self.tags:
      if not tag in ['img']:
        for i in range(len(self.tagstarts[tag])):
          element = Element()
          element.startTag = self.tagstarts[tag][i][0]
          element.endTag = self.tagstarts[tag][i][1]
          tagData = self._root[element.startTag+1:element.endTag].replace("\"","")
          element.tag = tag
          element.attrib=dict([a.split('=') for a in tagData.split() if '=' in a])
          element.start = element.startTag
          element.end = self.tagends[tag][i][1]
          element.innerHTML = (self.tagstarts[tag][i][1]+1,self.tagends[tag][i][0])
          self.elements.append(element)
def getText(self,tag,key,val):
    element = [element for element in self.elements if element.tag==tag and element.attrib[key]==val]
    if len(element):
      element = element[0]
      start,end = element.innerHTML
      text = self._root[start:end]
      i,j=0,0
      while i<len(text):
        i = text.find('<',i)
        j = text.find('>',i)
        if i<0 or j<0: break
        text = text[:i]+text[j+1:]
      i,j=-1,-1
      while i<len(text):
        i = text.find('&',i+1)
        j = text.find(';',i+1)
        if i<0 or j<0: break
        if text[i+1:j] in ent.keys(): text = text[:i]+ent[text[i+1:j]]+text[j+1:]
    else:
      text=None
    return text

Note the code does not check if the html code is correct. it also work only for the body part, and ignore img tags.

When you wish to login to a remote system through an ssh session without using any password you'll need to produce a key on the local computer and on the remote computer in order to allow password - free entrance.
first you'll need to produce the key in the local computer:
ssh-keygen -t dsa -b 1024 -P ""
press Enter on the questions followed by the command.
the output will be something like:
Your identification has been saved in /home/yourname/.ssh/id_dsa.
Your public key has been saved in /home/yourname/.ssh/id_dsa.pub.
The key fingerprint is:
4f:9e:4c:17:8f:19:37:32:bc:45:4f:2a:a2:8e:1c:15 yourname@your.loacl.machine
The key's randomart image is:
+--[ DSA 1024]----+
|           .o    |
|           .oo. |
|          . o+   |
|      o    *.+. .|
|       .E . + . |
|      o .+.. . .|
|       +o.o. + |
|       o+E       |
|       o=.       |
+-----------------+

This will produce two files in you .ssh directory: id_dsa and id_dsa.pub

cat id_dsa.pub into ~/.ssh/authorized_keys on the remote system you want to log into.
change mod to 600:
chmod 600 ~/.ssh/authorized_keys

(source of the tip)

The problem:
A table of xyz data points can be viewed spatially with different methods:
GMT and ArcView are just a couple of options. I often get tables of spatial data points in Microsoft Excel (xp ver.) format or CSV of some kind. An easy (and free) way to visualize the data points on Windows would be to export the data and view it on Google Earth.
The solution:
points2kml (115kb).
In order to ease the work-flow I wrote a VBA macro add-in for Excel that exports selected fields (longitude, latitude, name and data) to a Google Earth KML file. the point2kml add-in and can be found on my scripts page.

The Problem:
How to format a hard drive on Linux.

The Solution:
in order to format a HD on Linux you'll need super user privileges.
you can use sudo if you have sudo privilages or in the terminal
switch to super user:

>su
(give super user pwd)
and ignore the sudo in the following commands. if you're system can't find the commands while in su mode, try to add /sbin/ before the command.
follow the steps below: marked with •
run commnads marked with >
commands output looks like this.
Warning: formatting your HD will erase all data on it!

Connect the HD to the computer, make sure it's on.
get the hard drive's (device) name on your system:

>df -h

Filesystem            Size  Used Avail Use% Mounted on

/dev/mapper/vg_insarlablin-lv_root

                      908G   77G  786G   9% /

/dev/sda1             194M   31M  153M  17% /boot

/dev/sdb1             1.4T  1.1T  193G  86% /home/novitsky

/dev/sdc1             1.4T  1.1T  193G  86% /home/novitsky/data

Now lets say our HD is sdb1.

lets build it's partition:

>sudo fdisk /dev/sdb

The number of cylinders for this disk is set to 182401.

There is nothing wrong with that, but this is larger than 1024,

and could in certain setups cause problems with:

1) software that runs at boot time (e.g., old versions of LILO)

2) booting and partitioning software from other OSs

   (e.g., DOS FDISK, OS/2 FDISK)

Command (m for help):

pressing m will get this menu:

Command action

   a   toggle a bootable flag

   b   edit bsd disklabel

   c   toggle the dos compatibility flag

   d   delete a partition

   l   list known partition types

   m   print this menu

   n   add a new partition

   o   create a new empty DOS partition table

   p   print the partition table

   q   quit without saving changes

   s   create a new empty Sun disklabel

   t   change a partition's system id

   u   change display/entry units

   v   verify the partition table

   w   write table to disk and exit

   x   extra functionality (experts only)

pressing p will print the current partition:

Disk /dev/sdb: 1500.3 GB, 1500301910016 bytes

255 heads, 63 sectors/track, 182401 cylinders

Units = cylinders of 16065 * 512 = 8225280 bytes

Disk identifier: 0x000923e6

   Device Boot      Start         End      Blocks   Id  System

/dev/sdb1               1      182401  1465136001   83  Linux

pressing d will delete the current partition. depending on the partition you want to delete (this case 1) press the number of the partition to be deleted. if all the disk should be formatted and there are several partitions repeat pressing d and the number of the partitions until you delete them all.
pressing n will create a new partition.
pressing p will create a primary partition. select 1 as the number of the partition. use the defaults (just press Enter) for the start and end cylinders.
pressing t will set the file system type.
pressing l will list the available file systems.
select 83 if you want a Linux file system
pressing w will write and exit to shell.

Now we will create a new file system on the new partition:
> sudo mkfs -t ext2 /dev/hdb1
remember to change the number at the end if your HD device number is different.
if your using different fs than ext2 (newer Linux systems have ext3 and fc11 has ext4) replace the ext2 in the command with the correct fc name.
the mkfs may take some time to run (depending on the computer and HD size).

To add a label to the partition use:
e2label /dev/sdb1 NEWLABEL
where sdb should be changed to the correct device and NEWLABEL changed to the new label.

Now add the device to the fstab file in order to mount it:
> sudo vi /etc/fstab
and add a line like:
LABEL=NEWLABEL MOUNTPOINT ext3 defaults 1 2
where NEWLABEL is the disk labe and MOUNTPOINT is where you want to mount the disk.
mount it using:
> sudo mount -a
Note you might need to change permissions and owner of the mount point:
> sudo chown username:usergroop MOUNTPOINT
when finnished - reboot.

When using more memory than your RAM Linux will use swap which is placed on your hard drive and is much slower. Sometimes after using swap there might be a situation where the RAM is free but some of your data is on the swap. There is a way to relocate the swap data to the RAM.
more details and a nice script to do is are available at: Ubuntu Community Documentation Swap FAQ page.

sudo swapoff -a

sudo swapon -a

The commands simply shut off swap and turning it on back again.
You'll have to make sure there is enough space on the RAM to accommodate the Swap data, this can be done using:

free

the output will show you the used space of Swap and free space of RAM:

             total       used       free     shared    buffers     cached
Mem:       8047896     772268    7275628          0      19540     129372
-/+ buffers/cache:     623356    7424540
Swap:     10108920          0   10108920

The example above show that no data (0) is stored on Swap.

Matplotlib module enables hillshade method (shade) using a LightSource class (v 0.99).
The problem is it uses the data itself as intensity and data. It is very useful for viewing a DEM but sometimes you would like the DEM as intensity underlying some other data. Another problem is that the shade method is producing a very light colored image sometimes even white where intensity is high.
I used as an example a DEM derived from SRTM v4 data acquired at the International Centre for Tropical Agriculture (CIAT - http://srtm.csi.cgiar.org) the hillshade production was made using LightSource class with azimuth of 165 deg. and altitude of 45 deg.)


DEM - gist_earth color scheme		Hill-shade (azdeg-165,altdeg-45)

matplotlib shade method		My shade method

The difference in the shading colors derived from the method used to produce it. While the matplotlib method uses "hard light" method I use a "soft light" method. the matplotlib is converting the RGB colors to HSV and then calculate the new saturation and value according to the intensity. I use a formula based on the description of ImageMagick's pegtop_light.which is much faster as it is a single formula. Another advantage is the option to use a separate layer as the intensity and another as the data used for colors.
The modified functions are hillshade and set_shade as follows:

#!/bin/env python
from pylab import *
def set_shade(a,intensity=None,cmap=cm.jet,scale=10.0,azdeg=165.0,altdeg=45.0):
''' sets shading for data array based on intensity layer
or the data's value itself.
inputs:
a - a 2-d array or masked array
intensity - a 2-d array of same size as a (no chack on that)
                    representing the intensity layer. if none is given
                    the data itself is used after getting the hillshade values
                    see hillshade for more details.
cmap - a colormap (e.g matplotlib.colors.LinearSegmentedColormap
              instance)
scale,azdeg,altdeg - parameters for hilshade function see there for
              more details
output:
rgb - an rgb set of the Pegtop soft light composition of the data and
           intensity can be used as input for imshow()
based on ImageMagick's Pegtop_light:
http://www.imagemagick.org/Usage/compose/#pegtoplight'''
if intensity is None:
# hilshading the data
    intensity = hillshade(a,scale=10.0,azdeg=165.0,altdeg=45.0)
else:
# or normalize the intensity
    intensity = (intensity - intensity.min())/(intensity.max() - intensity.min())
# get rgb of normalized data based on cmap
rgb = cmap((a-a.min())/float(a.max()-a.min()))[:,:,:3]
# form an rgb eqvivalent of intensity
d = intensity.repeat(3).reshape(rgb.shape)
# simulate illumination based on pegtop algorithm.
rgb = 2*d*rgb+(rgb**2)*(1-2*d)
return rgb

def hillshade(data,scale=10.0,azdeg=165.0,altdeg=45.0):
''' convert data to hillshade based on matplotlib.colors.LightSource class.
    input:
         data - a 2-d array of data
         scale - scaling value of the data. higher number = lower gradient
         azdeg - where the light comes from: 0 south ; 90 east ; 180 north ;
                      270 west
         altdeg - where the light comes from: 0 horison ; 90 zenith
    output: a 2-d array of normalized hilshade
'''
# convert alt, az to radians
az = azdeg*pi/180.0
alt = altdeg*pi/180.0
# gradient in x and y directions
dx, dy = gradient(data/float(scale))
slope = 0.5*pi - arctan(hypot(dx, dy))
aspect = arctan2(dx, dy)
intensity = sin(alt)*sin(slope) + cos(alt)*cos(slope)*cos(-az - aspect - 0.5*pi)
intensity = (intensity - intensity.min())/(intensity.max() - intensity.min())
return intensity

Example of use:
One can save the code to a file named say: shading.py
now say we have a 4 byte float DEM data in a 560 lines 420 samples binary file. in a python code:

from pylab import *
from shading import set_shade
from shading import hillshade
dem = fromfile('DEM.dem',dtype=float32).reshape(560,420)
rgb = set_shade(dem,cmap=cm.gist_earth)
imshow(rgb)

will produce the "my shade method" image as above.

say we have a data to be plot using the DEM data as intensity:

replace the line before last with:
rgb = set_shade(data,intensity=hillshade(dem),cmap=cm.gist_earth)

After creating a DEM file using SRTM.py script, now say we need to convert it to a GMT .grd file.

1. convert the binary file to a grd file:
xyz2grd [binary file] -G[output_grd_file] -Zh -I[pixel_size] -R[west]/[south]/[east]/[north]r
- don't forget the r at the end.
- in the -Zh, h is for 2 byte integer. change if needed
2. Sample to a different pixel size and sub-region:
grdsample [input_grd_file] -G[output_grd_file] -I[new_pixel_size] -R[west]/[south]/[east]/[north]r
- if a known number of pixels is desired, say in east-west, use: east = west+[needed_pixel_number - 1]*[new_pixel_size]
3. put zeros instead of null values -9999:
grdmath [input_grd_file] -9999 GT [input_grd_file] MUL = [output_grd_file]
- note that if you sample some of the null pixels will be with different values so use say -500 instead of -9999 as a threshold.

When reading GTM file in GRD format you actually read a NetCDF format.
using python and pylab to read it enable matplotlib to plot it.
say we have a grd file named z.grd which contains 3 variables: x,y and z.
reading z values in python is done by the following three lines:

from pylab import *
from scipy.io import netcdf_file as netcdf
data = netcdf('z.grd','r').variables['z'][::-1]

Now data is a numpy array containing z values in the file dimensions.
for more info look at: http://gfesuite.noaa.gov/developer/netCDFPythonInterface.html

In order to replace a regexp in multiple files using the command line use:

sed -i -e 's/regexp/new_regexp/' file_name

where file_name can be also a regexp (e.g. *.txt)

July 22, 2009: One of my hard drives turned out to be dead.
When I disconnected it from the computer - all is well, the computer is working fine without the missing HD, only I had very important data on it and must resurrect it.
when I connect it - the whole system is freezing So I belive it is something to do with the logic board. I googled it up and came upon this.
Basically, I need another HD with the same properties (manufacturer, model,firmware) and just switch the logic board which should be easy and hopefully work.
I now need to find a similar HD for the replacement.

July 23, 2009: Well, no luck in finding a similar HD and sending the HD to data recovery is more than 150$. Before I spend so much money I would like to check the cheapest option:
Freezing the HD. lets see what happens...

July 23, 2009, supplemental: it worked!
I was able to recover the data after few hours in the freezer. I hardly believed it, but it actually worked. I'll now spend my money buying a larger HD.
wooohooo!

The Problem:
I can't get mdx/dgx in order to view my ROI PAC files.
The Solution:
a Python script called roiview (Get RoiView at

) wich is one of my other scripts.
This script uses pylab and can read, display and save images of ROI PAC quite similar to mdx/dgx.
it can handle *.unw, *.cor, *.slc, *.int and *.dem files.
just download, extract and use:
"roiview -h" for manual.
Update! Oct 13, 2010: A new version of RoiView (0.75) is now available at:

Thursday, August 29, 2013

Thursday, August 1, 2013

Sunday, July 7, 2013

Sunday, June 30, 2013

Saturday, June 1, 2013

Wednesday, January 16, 2013

Tuesday, January 1, 2013

Sunday, December 16, 2012

Thursday, October 25, 2012

Wednesday, October 24, 2012

Monday, September 24, 2012

Tuesday, August 28, 2012

Monday, August 27, 2012

Tuesday, May 29, 2012

Thursday, May 10, 2012

Monday, November 21, 2011

Saturday, July 30, 2011

Thursday, July 28, 2011

Thursday, July 14, 2011

Thursday, March 31, 2011

Wednesday, October 13, 2010

Tuesday, September 28, 2010

Sunday, July 18, 2010

Friday, June 25, 2010

Wednesday, June 2, 2010

Wednesday, April 28, 2010

Wednesday, April 21, 2010

Thursday, December 10, 2009

Friday, October 2, 2009

Tuesday, September 29, 2009

Wednesday, July 22, 2009

Monday, May 25, 2009