giovedì 26 novembre 2015

Primi passi con Octave

Durante il dottorato ho avuto il modo di utilizzare (a livello decisamente base) Matlab apprezzandone la sintassi compatta. Purtroppo non ho piu' disponibilita' della licenza universitaria e quindi ho deciso di provare Octave, una versione open source simile a Matlab (anche se chi lavora in modo serio con Matlab mi dice che Octave e' ancora ad anni luce di distanze per alcune funzionalita')

L'aspetto interessante e' che con la versione 3.8 e' stata implementata una GUI simile a Matlab (prima si lavorava praticamente a linea di comando da shell)

Il programma e' stato compilato da sorgenti su Centos 6
Le dipendenze sono
yum install gcc gcc-c++ kernel-devel make mercurial libtool libtool-ltdl-devel libtool-ltdl autoconf cmake lapack-devel \ lapack pcre-devel readline-devel readline fftw-devel glpk-devel suitesparse suitesparse-devel gnuplot libcurl-devel zlib-devel \ flex texlive gperf fltk-devel qhull-devel hdf5-devel gl2ps-devel qrupdate-devel arpack-devel qscintilla-devel llvm-devel qt-devel \ bison ghostscript-devel librsvg2-tools icoutils readline pcre
Dopo la compilazione (un po' lunghina per la verita' ...ma la mia macchina e' un po' vecchiotta) si lancia il programma con

octave --force-gui



Questo e' uno script base su Octave che dovrebbe funzionare senza troppi problemi su Octave
----------------------------------------------------------------------------
clear all
close all
echo off

function rad = radians(degree)
    rad = degree .* pi / 180;
end;

function [a,c,dlat,dlon]=haversine(lat1,lon1,lat2,lon2)
    dlat = radians(lat2-lat1);
    dlon = radians(lon2-lon1);
    lat1 = radians(lat1);
    lat2 = radians(lat2);
    a = (sin(dlat./2)).^2 + cos(lat1) .* cos(lat2) .* (sin(dlon./2)).^2;
    c = 2 .* asin(sqrt(a));
    #arrayfun(@(x) printf("distance: %.8f m\n",6372800 * x), c);
end;

[volt,lat,lon,quota]=textread("tetto_6_11_ridotto.csv", "%u;%u;%u;%u");
mlat = mean(lat);
mlon = mean(lon);
minlat = min(lat);
maxlat = max(lat);
minlon = min(lon);
maxlon = max(lon);
stdlat = std(lat);
stdlon = std(lon);
[r,c] = size(lat);
fprintf("Media Lat: %d\n",mean(lat));
fprintf("Media Lon: %d\n",mean(lon));
fprintf("Std Lat: %d\n",std(lat));
fprintf("Std Lon: %d\n",std(lon));
fig = figure;
plot(lat,lon);
axis([minlat maxlat minlon maxlon]);
title('Originali')
print(fig,'gps','-dpng');
centro = centroid(lat,lon);
fprintf("Centroide %d %d\n",centro);
#plot(volt)
#print(fig,'volt','-dpng')
disp ("-----------------")
mlat = double(mlat)/10000000.0
mlon = double(mlon)/10000000.0
#fprintf("MLat: %f\n",mlat);

ripuliti = zeros(r,2);
indice = 1;

for n = 1:r
  #fprintf("MLat: %f\n",mlat);
  #fprintf("MLon: %f\n",mlon);
  #fprintf("Lat: %f\n",double(lat(n))/10000000.0);
  #fprintf("Lon: %f\n",double(lon(n))/10000000.0);
  [a,c,dlat,dlon] = haversine(mlat,mlon,double(lat(n))/10000000.0,double(lon(n))/10000000.0);
  #fprintf("Distanza %i :%f\n",n,double(c)*6378160.0);
  if (c*6378160.0) < 5.0
        ripuliti(indice,1) = double(lat(n))/10000000.0;
        ripuliti(indice,2) = double(lon(n))/10000000.0;
        indice = indice + 1;
        end
end

fig2 = figure;
plot(ripuliti(1:indice-1,1),ripuliti(1:indice-1,2))
axis([double(minlat)/10000000.0 double(maxlat)/10000000.0 double(minlon)/10000000.0 double(maxlon)/10000000.0]);
title("Filtrati")
mflat = mean(ripuliti(1:indice-1,1));
mflon = mean(ripuliti(1:indice-1,2));
devflat = std(ripuliti(1:indice-1,1));
devflon = std(ripuliti(1:indice-1,2));

fprintf("Media Filtrata Lat: %.7f\n",mflat);
fprintf("Media Filtrata Lon: %.7f\n",mflon);
fprintf("DevSt Filtrata Lat: %.7f\n",devflat);
fprintf("DevSt Filtrata Lon: %.7f\n",devflon);
----------------------------------------------------------------------------


Octave ha un corrispettivo dei toolbox di Matlab, ovvero di pacchetti di estensioni dei comandi che possono essere scaricati, a seconda delle necessita' da Octave Forge
Se per esempio si vuole leggere/scrivere una porta seriale si deve scaricare il pacchetto instrument control in formato tar.gz e si installa con il comando

pkg install instrument-control-0.2.1.tar.gz

Di seguito un semplice esempio per leggere i dati seriali (da Arduino). Per far utilizzare la porta seriale ad un normale utente, quest'ultimo deve essere aggiunto al gruppo dialout (adduser luca dialout)
----------------------------------------------------------------------------
#! /usr/local/bin/octave -qf
pkg load instrument-control
s1 = serial("/dev/ttyACM0");
srl_flush(s1);
get(s1);
set(s1,"baudrate",9600);
set(s1,"bytesize",8);
set(s1,"parity",'N');
set(s1,"stopbits",1);
get(s1);
while true
   valore = srl_read(s1,6);
   char(valore)
   pause(1);
endwhile
srl_close(s1)



----------------------------------------------------------------------------

Per quanto riguarda la connettivita' con i database la versione attuale si interfaccia solo con PostgreSql mediante il pacchetto database 

pkg install database-2.3.2.tar.gz
(deve essere installato il pacchetto devel di postgresql, esempio ripreso da qui. Prima di arrivare ad una connessioni con Postgres c'e' un post da smanettare con la configurazione del server, in particolare con il file pg_hba.conf)

----------------------------------------------------------------------------
pkg load database
conn = pq_connect (setdbopts ("dbname", "test")); 
pq_exec_params (conn, "create table testtable (t text, i int2, b bytea);") 
pq_exec_params (conn, "insert into testtable values ($1, $2, $3);", {"name1", 1, uint8([2, 4, 5])}) 
 pq_exec_params (conn, "select * from testtable;")
 pq_exec_params (conn, "drop table testtable;")
 pq_close (conn);
----------------------------------------------------------------------------


I programmi in Octave possono essere lanciati anche all'esterno dell'ambiente di sviluppo rendendo eseguibile il file .m ed inserendo il testa al file .m la riga
#!/usr/local/bin/octave -qf