X-Git-Url: https://git.creatis.insa-lyon.fr/pubgit/?p=CreaPhase.git;a=blobdiff_plain;f=octave_packages%2Fsecs2d-0.0.8%2FUtilities%2FUise2pde.m;fp=octave_packages%2Fsecs2d-0.0.8%2FUtilities%2FUise2pde.m;h=2fc5389c80b36f3403d22a6f059de8bac1a9239f;hp=0000000000000000000000000000000000000000;hb=f5f7a74bd8a4900f0b797da6783be80e11a68d86;hpb=1705066eceaaea976f010f669ce8e972f3734b05

diff --git a/octave_packages/secs2d-0.0.8/Utilities/Uise2pde.m b/octave_packages/secs2d-0.0.8/Utilities/Uise2pde.m
new file mode 100644
index 0000000..2fc5389
--- /dev/null
+++ b/octave_packages/secs2d-0.0.8/Utilities/Uise2pde.m
@@ -0,0 +1,672 @@
+function [mesh,data_v]=Uise2pde(grid_file,pref,data_file_v,load_data,out_file)
+
+
+## [mesh,data]=ise2pde3(grid_file,pref,data_file,load_data,out_file)
+## ise2pde3
+## estrae dati dal formato DF-ISE di ISE a pdetool di Matlab
+## grid_file contiene il nome del file di griglia da estrarre
+## pref un prefisso che verra' dato ai files temporanei creati da grep
+## data_file e' un cell array delle file da estrarre
+## load_data e' un cell array che contiene i nomi delle grandezze da estrarre 
+## out_file e' il nome del file matlab opzionale per salvare i dati estratti
+##
+## 17-3-2004 ver 3.1 
+## Marco Bellini marco_bellini_1@yahoo.it
+## 14.02.2007 ver 3.2
+## Octave porting and bug fixes Carlo de Falco 
+
+##%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+##
+##  CHANGES:
+##   - riconosce i nomi dei contatti e delle regioni
+##   - tratta anche dispositivi con una sola regione
+##   - fornisce informazioni sulla conversione
+##
+##%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+
+##%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+##
+##  TO DO LOG:
+##
+##%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+
+
+##esempio per debug
+## prova a costruire una grid per pdetool
+##data_file="/home/marco/IseDB/jfet/jsdd2_mdr.dat";
+##grid_file="/home/marco/IseDB/jfet/jsdd2_mdr.grd";
+##pref="/home/marco/IseDB/jfet/jsdd2";
+
+## esempio per struttura onda viaggiante
+## data_file={"/home/marco/IseDB/prova/spice/on1_onda1_tr_des.dat"};
+## grid_file="/home/marco/IseDB/prova/spice/on_mdr.grd";
+## pref="/home/marco/IseDB/prova/spice/on";
+## load_data={"ElectrostaticPotential"};
+
+## dati da estrarre
+##load_data={"DopingConcentration","BoronActiveConcentration","PhosphorusActiveConcentration"};
+##load_data={"ElectrostaticPotential","eDensity","hDensity","eCurrentDensity","BeamGeneration"};
+
+##pdeplot(p1,e1,t1,"xydata",log10(abs(data(1,:)))',"zdata",log10(abs(data(1,:))));
+##[p1,e1,t1]=ise2pde(grid_file,pref);
+
+## leggo i vertici
+## i punti sono ordinati per regione
+## se la prima regione e' l'ossido ed ha 269 punti
+## i primi 269 valori in p1 si riferiscono ai punti dell'ossido
+## 
+## nei file di des.cmd di simulazione i punti vengono dati per ogni regione
+## quindi i punti della frontiera tra ossido e silicio vengono ripetuti due
+## volte (tra l'ossido e tra il silicio): e' necessario prenderli una volta
+## sola
+## siccome con piu' di due materiali diversi il discorso e' piu' complicato
+## si sceglie di considerare in questo programma solo un ossido ed il 
+## silicio
+
+
+## p1 conterra' le coordinate x e y dei vertici
+[s,w]=unix(sprintf("grep Vertices %s",grid_file));
+n_vert_str = regexp(w,'([0-9]+)',"Tokens");
+n_vert=str2num(n_vert_str{1}{1});
+unix(sprintf("grep Vertices -A %d %s | tail -n+2 > %s_vertex.txt",n_vert,grid_file,pref));
+p1=load(strcat(pref,"_vertex.txt"));
+unix(sprintf("rm %s_vertex.txt",pref));
+p1(:,2)=-p1(:,2);
+
+
+fprintf("Found %d vertex\n",n_vert);
+
+##leggo gli edge
+## el conterra' l'indice dei vertici degli edges degli elementi
+## cioe' l'indice dei vertici dei lati dei triangoli e l'indice dei vertici
+## dei segmenti
+[s,w]=unix(sprintf("grep Edges %s",grid_file));
+n_edges_str = regexp(w,'([0-9]+)',"Tokens");
+n_edges=str2num(n_edges_str{1}{1});
+unix(sprintf("grep Edges -A %d %s | tail -n+2  > %s_edges.txt",n_edges,grid_file,pref));
+el=load(strcat(pref,"_edges.txt"));
+unix(sprintf("rm %s_edges.txt",pref));
+fprintf("Found %d edges\n",n_edges);
+clear n_edges;
+el=el+1;
+
+
+##leggo gli elementi triangolari
+## el_tr contiene gli indici degli edge dei triangoli
+[s,w]=unix(sprintf("grep Elements %s",grid_file));
+n_els_str = regexp(w,'([0-9]+)',"Tokens");
+n_els=str2num(n_els_str{1}{1});
+## leggo solo gli elementi che iniziano per 2 e che corrispondono ai
+## triangoli, creando un file che contiene solo gli elementi triangolari
+unix(sprintf("grep Elements -A %d %s | head -n %d  | tail -n %d | awk '$1==2 {print $2,$3,$4}'  > %s_elements2.txt",n_els,grid_file,n_els+1,n_els,pref));
+el_tr=load(strcat(pref,"_elements2.txt"));
+unix(sprintf("rm %s_elements2.txt",pref));
+fprintf("Found %d triangular elements out of %d elements\n",length(el_tr),n_els);
+
+
+## creo un file che contiene gli elementi "segmenti"
+unix(sprintf("grep Elements -A %d %s | head -n %d  | tail -n %d | awk '$1==1 {print $2,$3}'  > %s_elements1.txt",n_els,grid_file,n_els+1,n_els,pref));
+##el_lin=load(strcat(pref,"_elements1.txt"));
+##unix(sprintf("rm %s_elements1.txt",pref));
+
+
+## creo un indice che dice se l'elemento e' un triangolo o una linea 
+## e lo salvo nel file pref_linee1.txt; ci sara' un 1 se l'elemento i-esimo
+## e' una linea o un 2 se e' un triangolo
+unix(sprintf("grep Elements -A %d %s | head -n %d  | tail -n %d  | awk ' {print $1} ' > %s_linee1.txt",n_els,grid_file,n_els+1,n_els,pref));
+clear n_els;
+
+## leggo le regions
+## cosi posso distinguere tra silicio, contatti e ossido
+## al max e' possibile leggere 60 regions
+[s,w]=unix(sprintf("grep regions -n %s  | grep -v nb_regions | tr -d []=",grid_file));
+w2=w(11:length(w));
+n_regions=0;
+l2=length(w2);cont=1;num=1;str="";begin=0;
+while (cont<l2)
+    c=w2(cont);
+    cont=cont+1;
+    if (c=="""")
+        if begin==0
+            begin=1;
+            str="";
+        else
+            begin=0;
+            regions{num}=str;
+            num=num+1;  n_regions= n_regions+1;
+        end; %begin
+    end; ## "
+    if ((c~="""") & (c~=" ") & (begin==1) )
+        str=strcat(str,c);
+    end; ## aggiungo carattere
+end ## while
+
+fprintf("Found %d regions/contacts:\n",n_regions);
+for cont=1:n_regions
+    fprintf("%d: %s\n",cont,char(regions(cont)));
+end
+
+## leggo le locations
+## trova le righe tra cui Ú compreso Locations
+[s,w]=unix(sprintf("grep Locations -n %s  | awk '/Locations/ {print $1}'",grid_file));
+loc_u=sscanf(w,"%d:Locations\n");
+[s,w]=unix(sprintf("grep Elements -n %s  | awk '/Elements/ {print $1}'",grid_file));
+tmp=sscanf(w,"%d:Elements\n");
+loc_d=tmp(1);
+## generazione file
+unix(sprintf(" head -n %d %s | tail -n+%d | tr [:cntrl:] ' ' | tr -d }| tr ief 012 > %s_loc.txt",...
+    loc_d-1,grid_file ,loc_u+1,pref) );
+loc=load(sprintf("%s_loc.txt",pref));
+## loc contiene
+## 0 per gli elementi interni
+## 1 per gli esterni
+## 2 per la frontiera
+unix(sprintf("rm %s_loc.txt",pref));
+clear loc_u;clear loc_d;
+
+
+
+
+##%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+##   CONVERSIONE DELLA GRIGLIA
+##%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+n_ed=max(size(el));
+n_el=max(size(el_tr));
+
+##%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+## calcolo di E1
+##%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+## trova gli elementi sul bordo (esterni)
+## del dispositivo (quelli rossi)
+ind=find(loc==1);
+e1=el(ind,:);
+
+##%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+## calcolo di T1
+##%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+el_si=sign(el_tr)+1;
+el_si=sign(el_si);
+## el_si  contiene le posizioni degli elementi positivi
+t1=zeros(n_el,4);
+
+## elementi tutti positivi: vertice di testa
+ind=find( el_si(:,:)==1);
+t1(ind)=el(el_tr(ind)+1,1);
+
+## elementi negativi: vertice di coda
+ind=find( el_si(:,:)==0);
+t1(ind)=el(-el_tr(ind),2);
+
+## inserisco il numero di subdominio
+
+## trovo il delimitatore inferiore di una regione
+[s,w]=unix(sprintf("grep material -n %s | grep -v materials  | awk '/material/ {print $1 }'",grid_file));
+mat_pos=sscanf(w,"%d:\n");
+[s,w]=unix(sprintf("grep material -n %s | grep -v materials  | awk '/material/ {print $2 }' |tr -d =",grid_file));
+mat_type_lim=findstr(w,"\n");
+if isempty(mat_type_lim)
+  mat_type=w;
+else
+  for ii=1:length(mat_type_lim)-1
+    mat_type{ii}=w(mat_type_lim(ii)+1:mat_type_lim(ii+1)-1);
+  end
+  mat_type{ii+1}=w(mat_type_lim(ii)+1:end);
+end
+
+[s,w]=unix(sprintf("grep material -n %s | grep -v materials  | awk '/material/ {print $4 }' |tr -d =",grid_file));
+mat_name_lim=findstr(w,"\n");
+if isempty(mat_name_lim)
+  mat_name=w;
+else
+  for ii=1:length(mat_name_lim)-1
+    mat_name{ii}=w(mat_name_lim(ii)+1:mat_name_lim(ii+1)-1);
+  end
+  mat_name{ii+1}=w(mat_name_lim(ii)+1:end);
+end
+
+
+## leggo il tipo di elemento
+## se e' un materiale i triangoli appartenenti a quel materiale otterranno
+## indice pari all'indice della regione. t1(4,:)' perche' non ho ancora
+## trasposto
+
+## se e' un contatto, i segmenti appartenenti a quel contatto saranno
+## indicati con un numero in e(5,:)' (non trasposto)
+
+el_type=load(strcat(pref,"_linee1.txt"));
+indice_tr=zeros(size(el_type));
+## indice_tr e' 1 nelle posizioni che corrispondono agli elementi triangolari nella lista
+## totale degli elementi nel file *.grd
+indice_tr(find(el_type==2))=1;
+## ora faccio la somma cumulativa in modo che nella posizione i-esima del
+## vettore degli elementi totali ci sia l'indice che identifica il triangolo
+## i-esimo in el_tr : ho costruito una tabella di conversione dall'indice di
+## ISE unico per triangoli e segmenti -> indice dei soli triangoli.
+indice_tr=cumsum(indice_tr);
+clear ind;
+
+## leggo gli elementi "segmento"
+try
+el_lin=load(strcat(pref,"_elements1.txt"));
+catch
+el_lin=-1;
+end_try_catch
+
+## aggiungo 1 perche' gli elementi partono da 0
+el_lin=el_lin+1;
+unix(sprintf("rm %s_elements1.txt",pref));
+## costruisco una mappa di conversione anche per gli elementi lineari (indicati da 1 in el_tr).
+indice_lin=zeros(size(el_type));
+indice_lin(find(el_type==1))=1;
+indice_lin=cumsum(indice_lin);
+##
+## creo un vettore che indica il tipo di elemento di frontiera
+## 1 per i lati del dispositivo (elementi esterni)
+## 2 per l'interfaccia ossido-silicio
+## da 3 in poi per i contatti
+
+e_fron=ones(max(size(e1)),1);
+
+## AGGIUNGO LA FRONTIERA
+##
+clear ind1;
+ind1=find(loc==2);
+name_contact{1}="External";
+name_contact{2}="I. Frontier";
+
+if (isempty(ind1)==0)
+## aggiungo gli elementi di frontiera trovati al vettore degli edges e1
+## mettendo indice 2
+    e1=[e1;el(ind1,:)];
+    e_fron=[e_fron;2*ones(length(ind1),1)];
+end
+
+clear el_type
+
+## numero del contatto
+n_contact=3;
+## name_material contiene il nome del materiale i-esimo
+## name_contact contiene il nome del contatto i-esimo
+
+
+for n=1:n_regions
+    
+    if strcmp(mat_type{n},"Contact")==0
+        
+## leggo gli elementi che costituiscono una regione
+        if (n~=n_regions)        
+            unix(sprintf(" head -n %d %s | tail -n+%d |  tr -d [:cntrl:]}  > %s_tmp2.txt", mat_pos(n+1)-2,grid_file,mat_pos(n)+2,pref));
+        else
+## trattare l'ultima regione che e' diversa
+            unix(sprintf(" tail -n+%d %s |  tr -d [:cntrl:]}  > %s_tmp2.txt",mat_pos(n)+2, grid_file,pref) );
+        end ## e' l'ultima regione?
+        tmp_el=load(strcat(pref,"_tmp2.txt"));
+        tmp_el=tmp_el+1;
+## trova i triangoli che fanno parte della regione n
+        t1(indice_tr(tmp_el),4)=n;
+        name_material{n}=mat_name{n};        
+    else
+## la regione e' un contatto
+        
+        if (n~=n_regions)        
+            unix(sprintf(" head -n %d %s | tail -n+%d |  tr -d [:cntrl:]}  > %s_tmp2.txt", mat_pos(n+1)-2,grid_file,mat_pos(n)+2,pref));
+        else
+## trattare l'ultima regione che e' diversa
+            unix(sprintf(" tail -n+%d %s |  tr -d [:cntrl:]}  > %s_tmp2.txt", mat_pos(n)+2,grid_file,pref) );
+        end ## e' l'ultima regione?
+        tmp_el=load(strcat(pref,"_tmp2.txt"));
+## tmp_el contiene gli indici dei vertici che appartengono al
+## contatto: aggiungo 1 perche' ise parte da 0.
+        tmp_el=tmp_el+1;
+        
+## predo i vertici la cui posizione e' in tmp_el e la converto da
+## dessis a pde
+        e1=[e1;el_lin(indice_lin(tmp_el),:)];
+        
+## aggiungo il numero del contatto
+        tmp_el=tmp_el';
+        e_fron=[e_fron;n_contact*ones(size(tmp_el))];
+        
+        name_contact{n_contact}=regions(n);        
+        n_contact=n_contact+1;
+        
+        
+    end ## non e' un contatto
+end ## fine scansione regioni
+
+unix(sprintf("rm %s_linee1.txt",pref));
+
+
+## trasposizione
+t1=t1';
+e1=e1';
+p1=p1';
+
+
+
+## individuo a che regione appartengono dei vertici
+## fondamentale per leggere i set di dati prodotti dalle simulazioni
+
+tmpt1=[t1(1,:),t1(2,:),t1(3,:);t1(4,:),t1(4,:),t1(4,:)];
+stmpt1=sortrows(tmpt1',1)';
+dtmpt1=diff(stmpt1(1,:));
+##clear stmpt1;
+ind1=find(dtmpt1==1);
+##%%% questa prende i valori dei vertici in ordine solo per TEST 
+##%%% reg_vert=[stmpt1(1,ind1),stmpt1(1,length(dtmpt1)+1)];
+
+## reg_vert contiene la regione a cui il punto appartiene
+reg_vert=[stmpt1(2,ind1),stmpt1(2,length(dtmpt1)+1)];
+clear tmpt1 stmpt1 dtmpt1 ind1;
+
+
+## individuo i vertici appartenenti alla frontiera
+##  le righe precedenti hanno assegnato questi punti
+##  a una o all'altra regione. Nei file delle simulazioni questi
+##  punti sono assegnati casualmente ad una delle regioni secondo
+##  un criterio arbitrario. E' necessario individuare i punti alla
+##  frontiera per tenerne conto.
+
+## costruisco e1 completo
+e1=[e1;zeros(size(e1))];
+e1=[e1;e_fron'];
+
+## trovo i punti di frontiera
+ind1=find(loc==2);
+## nel caso la frontiera non esista non fa niente
+
+
+
+if (isempty(ind1)~=1)
+## controlla se c'e' una frontiera tra diversi elementi
+    
+    el=el';
+    tmpe1=[el(1,ind1),el(2,ind1)];
+    stmpe1=sortrows(tmpe1',1)';
+    dtmpe1=diff(stmpe1(1,:));
+    clear ind1;
+## non e' detto che i vertici di frontiera siano ad incrementi unitari
+    ind1=find(dtmpe1>0);
+    
+## in teoria dovrei controllare che non ci sia un solo elemento di 
+## frontiera
+##    if (length(ind1)>1)
+    front_vert=[stmpe1(1,ind1),stmpe1(1,length(dtmpe1)+1)];
+##   else
+## c'e' solo 1 punto di frontiera
+##front_vert=[stmpe1(1,ind1)];
+##end
+    
+    
+    
+## inserisco anche nel reg_vert l'informazione se un punto e' di frontiera
+## in questo modo il valore assoluto indica la regione alla quale e' 
+## attribuito il valore.
+## i punti con valore negativo sono punti di frontiera
+    
+    reg_vert(front_vert)=-reg_vert(front_vert);
+    
+    clear tmpe1 stmpe1 dtmpe1 ind1;
+end % esiste la frontiera
+
+## % % test
+##   ind1=find(reg_vert==2);
+##   plot(p1(1,ind1),p1(2,ind1),'b.')
+##   whos ind1
+##   hold on;
+##   ind1=find(reg_vert==1);
+##   whos ind1
+##   plot(p1(1,ind1),p1(2,ind1),'g.')
+##   plot(p1(1,front_vert),p1(2,front_vert),'ro')
+##   whos front_vert
+
+
+
+## es. testato: ok
+## ind=find(reg_vert==2);
+## plot(p1(1,:),p1(2,:),'b.')
+## hold on;
+## plot(p1(1,ind1),p1(2,ind1),'g.')
+
+
+##%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+##         DATASETS             %
+##%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+if (exist("data_file_v")~=0)
+## guardo quanti file devo leggere
+    n_data_file=max(size(data_file_v));
+    
+    
+    for nf=1:n_data_file
+        
+        data_file=char(data_file_v(nf));
+        
+## leggo i nomi dei datasets
+        [s,w]=unix(sprintf("grep datasets %s",data_file));
+        w2=w(18:length(w)-1);
+        
+        ndatasets=0;
+        datasets=cell(30,1);
+        l2=length(w2);cont=1;num=1;str="";begin=0;
+        while (cont<l2)
+            c=w2(cont);
+            cont=cont+1;
+            if (c=="""")
+                if begin==0
+                    begin=1;
+                    str="";
+                else
+                    begin=0;
+## controllo che non ci siano due regioni in cui e'
+## presente lo stesso dato
+##%  if (num==1) || (strcmp(str,datasets(num-1))==0)
+##% tolgo cosi' ho le label con la molteplicita' giusta
+                    datasets(num)={str};
+                    num=num+1;  ndatasets= ndatasets+1;
+##%  end % non e' lo stesso dataset
+                end; %begin
+            end; ## "
+            if ((c~="""") & (c~=" ") & (begin==1) )
+                str=strcat(str,c);
+            end; ## aggiungo carattere
+        end ## while
+        
+        
+## estrae i valori di tutte le variabili in un unico file
+        unix(sprintf("grep Values  -A %d %s > %s_tmp.txt",n_vert,data_file,pref));
+        
+## estrae le righe che delimitano i dati da estrarre 
+        
+## del contiene le righe sotto le quali si trovano i valori sui vertici
+## dei triangoli (Values)
+        [s,w]=unix(sprintf("grep Values -n %s_tmp.txt  | awk '/Values/{print $1}' | tr -d : > %s_del.txt",pref,pref));
+        del=load(sprintf("%s_del.txt",pref));
+        unix(sprintf("rm %s_del.txt",pref));
+        
+## uel contiene le righe sopra le quali si trovano i valori sui vertici
+## dei triangoli (})
+        [s,w]=unix(sprintf("grep } -n %s_tmp.txt  | awk '/}/ {print $1}' | tr -d :  > %s_del.txt",pref,pref));
+        uel=load(sprintf("%s_del.txt",pref));
+        unix(sprintf("rm %s_del.txt",pref));
+        
+## nval contiene il numero valori sui vertici contenuto in questa
+## sezione del file
+        [s,w]=unix(sprintf("grep Values %s_tmp.txt |  tr -d  [:alpha:]{\\(\\) > %s_nval.txt",pref,pref));
+        nval=load(sprintf("%s_nval.txt",pref));
+        unix(sprintf("rm %s_nval.txt",pref));
+        
+        nload=max(size(load_data));
+        data=zeros(nload,n_vert);
+        
+##calcolo gli elementi ripetuti alla frontiera
+        n_rip=sum(loc==2)+1;
+        
+        nl=1;nd=1;
+        while(nl<=nload)
+            while(nd<=ndatasets)
+                if strcmp(load_data(nl),datasets(nd))==1
+                    
+## controllo se il numero dati sui vertici e' < del
+## numero totale dei vertici: accade quando c'e' dell'
+## ossido
+                    
+                    if (nval(nd)<n_vert)
+## caso in presenza di ossido
+                        
+## ho trovato il set
+## controllo se e' un caso particolare (ho un set di dati per l'ossido e uno per il silicio)
+                        if strcmp(load_data(nl),"ElectrostaticPotential")==1 | ...
+                                strcmp(load_data(nl),"ElectricField")==1  ...
+                                | strcmp(load_data(nl),"LatticeTemperature")==1
+                            
+                            
+## leggo la prima regione 
+                            dup=del(nd)+1;
+                            c=1; ## trova il delimitatore successivo 
+                            while uel(c)<=dup
+                                c=c+1;
+                            end;
+                            ddo=uel(c);
+                            [s,w]=unix(sprintf(" head -n %d %s_tmp.txt | tail -n+%d |  tr -d [:cntrl:]}  > %s_tmp2.txt", ddo, pref,dup,pref));
+                            tmp=load(strcat(pref,"_tmp2.txt")); ## tmp2 e' il vettore riga dei dati cercati
+## inserisco i primi elementi: ossido, cioe'
+## regione 2: gli elementi nel file .dat (var: tmp) sono
+## ossido interno + frontiera. 
+## ind1 e' appunto dato da ossido + frontiera
+                            
+                            ind1=find(abs(reg_vert)==2 | reg_vert<0);
+                            
+                            data(nl,ind1)=tmp;
+                            
+## leggo la seconda regione
+                            nd=nd+1;
+                            dup=del(nd)+1;
+                            c=1; ## trova il delimitatore successivo 
+                            while uel(c)<=dup
+                                c=c+1;
+                            end;
+                            ddo=uel(c);
+                            [s,w]=unix(sprintf(" head -n %d %s_tmp.txt | tail -n+%d |  tr -d [:cntrl:]}  > %s_tmp2.txt", ddo, pref,dup,pref));
+                            tmp=load(strcat(pref,"_tmp2.txt")); ## tmp2 e' il vettore riga dei dati cercati
+## inserisco i secondi elementi
+## silicio + frontiera
+                            clear ind;
+                            ind=find(abs(reg_vert)==1 | reg_vert<0 );
+                            data(nl,ind)=tmp;
+                            
+                        else
+## caso non particolare
+## metto a zero i dati della regione mancante
+                            ind1=find(abs(reg_vert)==2 | reg_vert<0);
+                            data(nl,ind1)=0;
+                            
+## leggo la regione 
+                            dup=del(nd)+1;
+                            c=1; ## trova il delimitatore successivo 
+                            while uel(c)<=dup
+                                c=c+1;
+                            end;
+                            ddo=uel(c);
+                            [s,w]=unix(sprintf(" head -n %d %s_tmp.txt | tail -n+%d |  tr -d [:cntrl:]}  > %s_tmp2.txt", ddo, pref,dup,pref));
+                            tmp=load(strcat(pref,"_tmp2.txt")); ## tmp2 il vettore riga dei dati cercati
+                            
+## inserisco i primi elementi
+                            ind=find(abs(reg_vert)==1 | reg_vert<0 );
+                            data(nl,ind)=tmp;
+                            
+##n_prev_region=length(tmp);
+##data(nl,1:n_prev_region)=tmp;
+                            
+                        end ## trattamento casi particolari
+                        
+                    else
+## il numero dei dati sui vertici e' = al numero
+## dei vertici: uso lo stesso codice di ise2pde 
+                        dup=del(nd)+1;
+                        c=1; ## trova il delimitatore successivo 
+                        while uel(c)<=dup
+                            c=c+1;
+                        end;
+                        ddo=uel(c);
+                        [s,w]=unix(sprintf(" head -n %d %s_tmp.txt | tail -n+%d |  tr -d [:cntrl:]}  > %s_tmp2.txt", ddo, pref,dup,pref));
+                        data(nl,:)=load(strcat(pref,"_tmp2.txt")); ## tmp2 e' il vettore riga dei dati cercati
+                    end ## numero dei vertici =
+                end %if trovato il set
+                nd=nd+1;   
+            end %nd
+            nd=1;
+            nl=nl+1;
+        end %nl
+        
+        clear nd nl;
+        
+        data_v(:,:,nf)=data(:,:);
+    end ## fine lettura n-esimo data_file
+    
+end ## if c'e' un data_file
+clear n_vert w2 l2 cont num str c;
+
+mesh.p=p1;
+mesh.e=e1;
+mesh.t=t1;
+mesh.materials=name_material;
+mesh.contacts=name_contact;
+
+for cont=1:length(name_material)
+    fprintf("Extracted region %d : %s\n",cont,char(name_material{cont}));
+end;
+for cont=3:length(name_contact)
+    fprintf("Extracted contact %d : %s\n",cont,char(name_contact{cont}));
+end;
+
+
+if exist("out_file")
+    if (exist("data_file_v")~=0)
+        save (out_file,"mesh","data");
+        fprintf("mesh and data saved");
+    else
+        save (out_file,"mesh");
+        fprintf("mesh saved");
+    end;
+end
+
+unix(sprintf("rm %s_tmp.txt",pref));
+unix(sprintf("rm %s_tmp2.txt",pref));
+
+##%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+## Cdf: Fix edge matrix and build output structure
+##%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+return
+mesh.e(7,end)=0;
+mesh = Umeshproperties(mesh);
+
+alledges = [mesh.t([1,2,4],:),mesh.t([2,3,4],:),mesh.t([3,1,4],:)];
+
+
+for iedge = 1:size(mesh.e,2)
+    whatedgeL = find((mesh.e(1,iedge)==alledges(1,:)& mesh.e(2,iedge)==alledges(2,:)));
+    whatedgeR = find((mesh.e(1,iedge)==alledges(2,:)& mesh.e(2,iedge)==alledges(1,:)));
+    
+    if (length(whatedgeL)==1)
+        mesh.e(6,iedge)=alledges(3,whatedgeL);
+    end
+    
+    if (length(whatedgeR)==1)
+        mesh.e(7,iedge)=alledges(3,whatedgeR);
+    end
+    
+end
+
+maxedge  = max(mesh.e(5,:));
+intedges = find((mesh.e(6,:)~=0)&((mesh.e(7,:)~=0)&((mesh.e(7,:)~=(mesh.e(6,:)) ) )));
+mesh.e (5,intedges) = maxedge+1;
+mesh.contacts{end+1} = "Interface";
+
+return;
\ No newline at end of file