probe_image_extract.C

Go to the documentation of this file.

/* Gnome gxsm - Gnome X Scanning Microscopy
 * universal STM/AFM/SARLS/SPALEED/... controlling and
 * data analysis software
 * 
 * Gxsm Plugin Name: probe_image_extract.C
 * ========================================
 * 
 * Copyright (C) 1999 The Free Software Foundation
 *
 * Authors: Percy Zahl <zahl@fkp.uni-hannover.de>
 * additional features: Andreas Klust <klust@fkp.uni-hannover.de>
 *
 * This program is free software; you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation; either version 2 of the License, or
 * (at your option) any later version.
 *
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with this program; if not, write to the Free Software
 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307, USA.
 */

/* Please do not change the Begin/End lines of this comment section!
 * this is a LaTeX style section used for auto generation of the PlugIn Manual 
 * Chapter. Add a complete PlugIn documentation inbetween the Begin/End marks!
 * All "% PlugInXXX" commentary tags are mandatory
 * All "% OptPlugInXXX" tags are optional
 * --------------------------------------------------------------------------------
% BeginPlugInDocuSection
% PlugInDocuCaption: Probe Image Extract
% PlugInName: probe_image_extract
% PlugInAuthor: Bastian Weyers
% PlugInAuthorEmail: weyers@users.sf.net
% PlugInMenuPath: Math/Probe/Probe_Image_Extract

% PlugInDescription
Extract selected probe data/values and generates a new image from that.
Select a \GxsmEmph{Rectangle}.

% PlugInUsage
Call \GxsmMenu{Math/Probe/Probe\_Image\_Extract}

% OptPlugInSources
The active channel is used as data source.

% OptPlugInObjects
A rectangle object is needed. Probe events must be present in the active scan.

% OptPlugInDest
The computation result is placed into an existing math channel, else into a new created math channel.

% OptPlugInConfig
The PlugIn configurator allows to set a default index andan chanel. If this is non
zero the user is not prompted!

% EndPlugInDocuSection
 * -------------------------------------------------------------------------------- 
 */

#include <gtk/gtk.h>
#include <glib.h>
#include <math.h>
#include "config.h"
#include "gxsm/plugin.h"


static void probe_image_extract_about( void );
static void probe_image_extract_configuration( void );
static gboolean probe_image_extract_run( Scan *Src, Scan *Dest );


GxsmPlugin probe_image_extract_pi = {
        NULL,
        NULL,
        0,
        NULL,
        "Probe_Image_Extract-M1S-PX",
        NULL,
        NULL,
        "Bastian Weyers",
        "_Math/_Probe/",
        N_("Image Extract"),
        N_("Probe Image Extract of selection of scan"),
        "no more info",
        NULL,
        NULL,
        NULL,
        NULL,
        probe_image_extract_about,
        probe_image_extract_configuration,
        NULL,
        NULL
};

GxsmMathOneSrcPlugin probe_image_extract_m1s_pi = {
  probe_image_extract_run
};

static const char *about_text = N_("Gxsm (Raster-) Probe Image Extraction Math-Plugin\n\n"
                                   "Probe Event Data extraction:\n"
                                   "A new image is generated from specified probe data.");

double IndexDefault = 0.;
double IndexLast = 0.;

GxsmPlugin *get_gxsm_plugin_info ( void ){ 
        probe_image_extract_pi.description = g_strdup_printf(N_("Gxsm MathOneArg probe_image_extract plugin %s"), VERSION);
        return &probe_image_extract_pi; 
}

GxsmMathOneSrcPlugin *get_gxsm_math_one_src_plugin_info( void ) { 
        return &probe_image_extract_m1s_pi; 
}

static void probe_image_extract_about(void)
{
        const gchar *authors[] = { probe_image_extract_pi.authors, NULL};
        gtk_widget_show(gnome_about_new ( probe_image_extract_pi.name,
                                          VERSION,
                                          N_("(C) 2000 the Free Software Foundation"),
                                          about_text,
                                          authors,
                                          NULL, NULL, NULL
                                ));
}

static void probe_image_extract_configuration( void ){
        probe_image_extract_pi.app->ValueRequest("Enter Number", "Index", 
                                                 "Default set index for data extraction.",
                                                 probe_image_extract_pi.app->xsm->Unity, 
                                                 0., 20., ".0f", &IndexDefault);
}

gint find_probe (gpointer *se, gpointer *default_probe_event){
        return (((EventEntry*) ((ScanEvent*)se)->event_list->data)->description_id () == 'P' ? 0:-1);
}

gint compare_events_distance (ScanEvent *a, ScanEvent *b, double *xy){
        double da = a->distance2 (xy);
        double db = b->distance2 (xy);
        if (da < db) return -1;
        if (da > db) return 1;
        return 0;
}

gint check_probe (ProbeEntry *a, ProbeEntry *b, char *Source2){
        //if (! 
        for (int j=0; j<b->get_chunk_size (); j++){
                if (! strncmp (b->get_label (j), Source2, 4)) return j;
        }

        return -1;
}

void ChangeValue (GtkComboBox* Combo,  gchar **string){
        GtkTreeIter iter;
        GtkTreeModel *model;
        gtk_combo_box_get_active_iter (GTK_COMBO_BOX (Combo), &iter);
        model = gtk_combo_box_get_model (GTK_COMBO_BOX (Combo));
        gtk_tree_model_get (model, &iter, 0, &(*string), -1);
        //std::cout << "Channel: " << *string << std::endl;
}

void ChangeIndex (GtkSpinButton *spinner, double *index){
        *index = gtk_spin_button_get_value (spinner);

        //std::cout << "Change Value " << index << std::endl;
}

void Dialog (ProbeEntry* defaultprobe, gchar **Source1, gchar **Source2, double *index){

        GtkWidget *dialog;
        GtkWidget *label;
        GtkWidget *vbox;
        GtkWidget *hbox;
        GtkWidget *source1_combo;
        GtkWidget *source2_combo;
        GtkWidget *spinner;
        GtkAdjustment *spinner_adj;
        
        char* Source1_lookup[]  = { "Zmon", "Umon", "Time", "XS", "YS", 
                                     "ZS", "Bias", "Phase", NULL};

        char* Source2_lookup[]  = { "ADC0", "ADC1", "ADC2", "ADC3", "ADC4", "ADC5-I","ADC6","ADC7",
                                     "Zmon", "LockIn0", "LockIn1stA", "LockIn1stB", "LockIn2ndA",
                                     "LockIn2ndB", NULL};
        
        dialog = gtk_dialog_new_with_buttons (_("Settings"),
                                                   NULL,
                                                   (GtkDialogFlags)(GTK_DIALOG_MODAL | GTK_DIALOG_DESTROY_WITH_PARENT),
                                                   GTK_STOCK_OK,
                                                   GTK_RESPONSE_ACCEPT,
                                                   NULL);

        vbox = gtk_vbox_new (FALSE, 0);
        gtk_widget_show (vbox);
        gtk_box_pack_start(GTK_BOX(GNOME_DIALOG(dialog)->vbox),
                     vbox, FALSE, FALSE, GNOME_PAD);

        hbox = gtk_hbox_new (FALSE, 0);
        gtk_container_set_border_width (GTK_CONTAINER (hbox), 5);
        gtk_widget_show (hbox);
        gtk_box_pack_start (GTK_BOX (vbox), hbox, FALSE, FALSE, 0);
                        
        label = gtk_label_new (_("Source for layer: "));        
        gtk_widget_show (label);                        
        gtk_box_pack_start (GTK_BOX (hbox), label, FALSE, FALSE, 0); 

        source1_combo = gtk_combo_box_new_text ();
        gtk_combo_box_append_text (GTK_COMBO_BOX (source1_combo), "Index");
        for (int i=0; i<defaultprobe->get_chunk_size (); i++){
                int j=0;
                while (Source1_lookup[j]){
                        if (!strncmp (defaultprobe->get_label (i), Source1_lookup[j], 4)){
                                gtk_combo_box_append_text (GTK_COMBO_BOX (source1_combo), defaultprobe->get_label (i));
                                break;
                        }
                        j++;
                }
        }
        gtk_combo_box_set_active (GTK_COMBO_BOX (source1_combo), 0);
        ChangeValue(GTK_COMBO_BOX (source1_combo), &(*Source1));                //init Source1
        gtk_widget_show (source1_combo);                        
        gtk_signal_connect(GTK_OBJECT (source1_combo), "changed", GTK_SIGNAL_FUNC (ChangeValue), &(*Source1));
        gtk_box_pack_start (GTK_BOX (hbox), source1_combo, FALSE, FALSE, 0);
         

        hbox = gtk_hbox_new (FALSE, 0);
        gtk_container_set_border_width (GTK_CONTAINER (hbox), 5);
        gtk_widget_show (hbox);
        gtk_box_pack_start (GTK_BOX (vbox), hbox, FALSE, FALSE, 0);
                        
        label = gtk_label_new (_("Channel to map: "));  
        gtk_widget_show (label);                        
        gtk_box_pack_start (GTK_BOX (hbox), label, FALSE, FALSE, 0); 

        source2_combo = gtk_combo_box_new_text ();
        for (int i=0; i<defaultprobe->get_chunk_size (); i++){
                int j=0;
                while (Source2_lookup[j]){
                        if (!strncmp (defaultprobe->get_label (i), Source2_lookup[j], 4)){
                                gtk_combo_box_append_text (GTK_COMBO_BOX (source2_combo), defaultprobe->get_label (i));
                                break;
                        }
                        j++;
                }
        }
        gtk_combo_box_set_active (GTK_COMBO_BOX (source2_combo), 0);
        ChangeValue(GTK_COMBO_BOX (source2_combo), &(*Source2));                //init Source2
        gtk_widget_show (source2_combo);                        
        gtk_signal_connect(GTK_OBJECT (source2_combo), "changed", GTK_SIGNAL_FUNC (ChangeValue), &(*Source2));
        gtk_box_pack_start (GTK_BOX (hbox), source2_combo, FALSE, FALSE, 0);

        hbox = gtk_hbox_new (FALSE, 0);
        gtk_container_set_border_width (GTK_CONTAINER (hbox), 8);
        gtk_widget_show (hbox);
        gtk_box_pack_start (GTK_BOX (vbox), hbox, FALSE, FALSE, 0);
                        
        label = gtk_label_new (_("Image Filter\n0 (sharp) ... 10 (smooth)"));   
        gtk_widget_show (label);                        
        gtk_box_pack_start (GTK_BOX (hbox), label, FALSE, FALSE, 0); 

        spinner_adj = (GtkAdjustment *) gtk_adjustment_new (IndexLast, 0.0, 10.0, 0.1, 1.0, 1.0);
        spinner = gtk_spin_button_new (spinner_adj, 0.1, 1);
        gtk_widget_show (spinner);                      
        gtk_signal_connect(GTK_OBJECT (spinner), "changed", GTK_SIGNAL_FUNC (ChangeIndex), &(*index));
        gtk_box_pack_start (GTK_BOX (hbox), spinner, FALSE, FALSE, 0);

/*      
        gtk_widget_set_usize (VarName, 100, -1);
        gtk_widget_show (VarName);
        gtk_label_set_justify(GTK_LABEL(VarName), GTK_JUSTIFY_LEFT);
        gtk_misc_set_alignment (GTK_MISC (VarName), 0.0, 0.5);
        gtk_misc_set_padding (GTK_MISC (VarName), 5, 0);
        gtk_box_pack_start (GTK_BOX (hbox), VarName, TRUE, TRUE, 0);

*/      
        gtk_dialog_run(GTK_DIALOG(dialog));
        gtk_widget_destroy (dialog);

}


// PlugIn run function, called by menu via PlugIn mechanism.
static gboolean probe_image_extract_run(Scan *Src, Scan *Dest)
{
        PI_DEBUG (DBG_L2, "Probe_Image_Extract Scan");

        int numsets = Dest->mem2d->GetNv ();
        int avgvectors = 5;             // minimum number of vectors used for interpolation
        double radius;                  // the catchment area of each point, usually there should be #acgvectors in it.
        double xy[2];                   // real world position of a probe_event
        double index = IndexLast;
        gchar *Source1;                 // is the channel for layers V
        gchar *Source2;                 // is the channel, that will be mapped
        int iSource1 = -1;              // -1 means "Index" that is not available as a real channel
        int iSource2 = -1;              // -1 means no channel to map -> break
        int jSource2 = -1;              // -1 means no channel to map -> break
        int Nx = Src->mem2d->GetNx ();  // total Number of points in x
        int Ny = Src->mem2d->GetNy ();  // total Number of points in y
        int xmin, xmax, ymin, ymax;     // region that is taken into account for averaging


        GSList *ProbeEventList = g_slist_find_custom (Src->mem2d->scan_event_list, NULL, (GCompareFunc)find_probe);
        ScanEvent* default_probe_event = (ScanEvent*) ProbeEventList->data;
        if (!default_probe_event) return MATH_OK;  // in case of no probe event in the scan
        ProbeEntry* defaultprobe = (ProbeEntry*) default_probe_event->event_list->data; // define a default pattern all probes should match

        Dialog (defaultprobe, &Source1, &Source2, &index);  // Asks for channels to use
        IndexLast = index; 

        // if no valid channel is found for Source1, we will use a numerator #
        PI_DEBUG (DBG_L2, "Source1 = " << Source1 << std::endl);
        PI_DEBUG (DBG_L2, "Source2 = " << Source2 << std::endl);
        PI_DEBUG (DBG_L2, "Index: " << index << std::endl);
        
        
        for (int j=0; j<defaultprobe->get_chunk_size (); j++){
                // std::cout << "Label = " << defaultprobe->get_label (j) << std::endl;
                if (! strncmp (defaultprobe->get_label (j), Source1, 4)) iSource1 = j;
                if (! strncmp (defaultprobe->get_label (j), Source2, 4)) iSource2 = j;
        }
        if (iSource2 < 0) return MATH_OK;  // No valid channel found for Source2, so stop at this point

        numsets = defaultprobe->get_num_sets ();             // #points of a probe = #layers
        Dest->mem2d->Resize(Dest->mem2d->GetNx (), Dest->mem2d->GetNy (), numsets, ZD_FLOAT);
        Dest->mem2d->data->MkVLookup(0., (double)numsets);
        Dest->data.s.nvalues=numsets;
        Dest->data.s.ntimes=1;
        Dest->data.s.dz=1.;
        
        if (Source1 < 0){
                UnitObj *NUnit = new UnitObj( " "," ",".0f", "#");
                Dest->data.SetVUnit (NUnit);
                delete NUnit;
                for (int n=0; n<numsets; n++) Dest->mem2d->data->SetVLookup(n, n);
        } else {
                UnitObj *NUnit = new UnitObj( defaultprobe->get_unit_symbol(iSource1), defaultprobe->get_unit_symbol(iSource1),".3f", defaultprobe->get_label (iSource1));
                Dest->data.SetVUnit (NUnit);
                delete NUnit;
                for (int n=0; n<numsets; n++) Dest->mem2d->data->SetVLookup(n, defaultprobe->get (n, iSource1));
        }               

        UnitObj *VoltUnit = new UnitObj( defaultprobe->get_unit_symbol (iSource2), defaultprobe->get_unit_symbol (iSource2),".3f", defaultprobe->get_label (iSource2));
        Dest->data.SetZUnit (VoltUnit);
        delete VoltUnit;


        double smoothness = pow(2., index-1.);
        radius = sqrt( (double)(Dest->mem2d->GetNx () * Dest->mem2d->GetNy () * avgvectors)/(double)(3*g_slist_length (Src->mem2d->scan_event_list)))+smoothness;
        PI_DEBUG (DBG_L2, "radius = " << radius << std::endl);

        PI_DEBUG (DBG_L2, "Initialising lookup table for " << g_slist_length (Src->mem2d->scan_event_list) << " probe events." << std::endl);

        // a 'matrix' with information about the position of probe events.
        int* matrix = new int [Nx * Ny];
        for (int i=0; i < Nx * Ny; i++) matrix[i]=-1;   // set all elements to -1, after analysing all probe events -1 will mean no probe at this point

        for (int n=0; n < g_slist_length (Src->mem2d->scan_event_list); n++){
                ScanEvent* se = (ScanEvent*) g_slist_nth_data (Src->mem2d->scan_event_list, n); // get the nth scan_event
                if (!se) {
                        PI_DEBUG (DBG_L2, "Empty Scan_Event found!");
                }
                else if (((EventEntry*) (se->event_list->data))->description_id () == 'P'){     // let's look if it is a valid probe_event
                                int xn,yn = 0;
                                ProbeEntry* pe = (ProbeEntry*) se->event_list->data;
                                se->get_position(xy[0], xy[1]);
                                Dest->World2Pixel(xy[0],xy[1], xn , yn);
                                if ( xn >= 0 && xn < Nx && yn >= 0 && yn < Ny ){
                                        matrix[yn*Nx+xn] = n;
                                        //PI_DEBUG (DBG_L2, "PE" << n << " @ ( " << xn << " | " << yn << " ) " );
                                }
        }}
        


        double avgvalues[numsets], avgnorm;
        int dead_pixel = 0;     //dead pixel indicator. Did not found a vektor in the near? Shit we have a dead pixel with all values equal to zero :(

        gapp->progress_info_new ("Converting Probe Data...", 1);
        gapp->progress_info_set_bar_fraction (0., 1);
        gapp->progress_info_set_bar_text (" ", 1);
//      gapp->progress_info_set_bar_fraction (0., 2);
//      gapp->progress_info_set_bar_text ("0 %% ", 2);

        for (int yi=0; yi < Dest->mem2d->GetNy (); yi++){
        // std::cout << (double) yi/(Dest->mem2d->GetNy ())*100. << std::endl;
        gapp->progress_info_set_bar_fraction ( (double) (yi+1)/Dest->mem2d->GetNy () , 1);
        gapp->progress_info_set_bar_text (g_strdup_printf("Line %i of %i", yi+1, Dest->mem2d->GetNy ()), 1);
/*
        std::cout << "Total rows: " << Dest->mem2d->GetNy () << std::endl;
        std::cout << "Total columns: " << Dest->mem2d->GetNx () << std::endl;
        std::cout << "Total spectra: " << g_slist_length (Src->mem2d->scan_event_list) << std::endl;
*/
        int counter = 0;
        int pe_number = -1;             // position of the probe event in the event list

        for (int xi=0; xi < Dest->mem2d->GetNx (); xi++){
                //gapp->progress_info_set_bar_fraction ( (double) (xi+1)/(Dest->mem2d->GetNy ()) , 2);
                //gapp->progress_info_set_bar_text (g_strdup_printf("Point %i of %i", xi+1, Dest->mem2d->GetNx ()), 2);

                int found_a_vector = 0;         // number of vectors found in the given perimeter round (xi,yi)
                for (int i=0; i<numsets; i++) avgvalues[i] = avgnorm = 0.;      // initialises avg vectors and the kalibration value

                if ( (xmin = xi - (int)radius) < 0) xmin = 0;
                if ( (xmax = xi + (int)radius) > Nx) xmax = Nx;
                if ( (ymin = yi - (int)radius) < 0) ymin = 0;
                if ( (ymax = yi + (int)radius) > Ny) ymax = Ny;

                for (int xpos=xmin; xpos < xmax; xpos++){
                for (int ypos=ymin; ypos < ymax; ypos++){
                        if ( (pe_number = matrix[ypos*Nx+xpos]) >= 0){
                                counter++;
                                // get the probe data of the scan_event found with matrix
                                ProbeEntry* pe = (ProbeEntry*) ((ScanEvent*) g_slist_nth_data (Src->mem2d->scan_event_list, pe_number))->event_list->data;      
                                jSource2 = check_probe(defaultprobe, pe, Source2);
                                if ( jSource2 >= 0 ){
                                        double norm = exp(-sqrt((double)((xpos-xi)*(xpos-xi)+(ypos-yi)*(ypos-yi))/smoothness));
                                        avgnorm += norm;
                                        for (int zi=0; zi < numsets; zi++){
                                                avgvalues[zi] += norm * pe->get(zi, jSource2);
                                        }
                                        found_a_vector++;
                                }
                        }
                        
                }}
                
/*              for (int n=0; n < g_slist_length (Src->mem2d->scan_event_list); n++){   // have a look on all available vectors and check if they are matching
//                      std::cout << "Counter: " << counter << std::endl;
                        counter++;
                        ScanEvent* se = (ScanEvent*) g_slist_nth_data (Src->mem2d->scan_event_list, n); // get the nth scan_event
                        if (!se) {
                                PI_DEBUG (DBG_L2, "Empty Scan_Event found!");
                        }
                        else if (((EventEntry*) (se->event_list->data))->description_id () == 'P'){     // let's look if it is a probe_event
                                int xn,yn = 0;
                                ProbeEntry* pe = (ProbeEntry*) se->event_list->data;
                                se->get_position(xy[0], xy[1]);
                                Dest->World2Pixel(xy[0],xy[1], xn , yn);
                                jSource2 = check_probe(defaultprobe, pe, Source2);
                                if ( jSource2 >= 0 && (xn-xi)*(xn-xi)+(yn-yi)*(yn-yi) < (int) (radius*radius)){
                                        double norm = exp(-sqrt((double)((xn-xi)*(xn-xi)+(yn-yi)*(yn-yi))/smoothness));
                                        avgnorm += norm;
                                        for (int zi=0; zi < numsets; zi++){
                                                avgvalues[zi] += norm * pe->get(zi, jSource2);
                                        }
                                        found_a_vector++;
                                }
                                //std::cout << ((EventEntry*) (se->event_list->data))->description_id () << std::endl;
                                //std::cout << "x: " << xy[0] << "  y: " << xy[1] << std::endl;
                                //std::cout << "x: " << xn << "  y: " << yn << std::endl;
                        }
                }
*/
                //std::cout << "found_a_vector: " << found_a_vector << std::endl;
                for (int i=0; i<numsets; i++) {
                        if (found_a_vector == 0){
                                dead_pixel++;
                                PI_DEBUG (DBG_L2, "dead_pixel at X" << xi << " Y" << yi << std::endl);
                        }
                        else {
                                for (int zi=0; zi<numsets; zi++){
                                        Dest->mem2d->PutDataPkt ( avgvalues[zi]/avgnorm, xi, yi, zi);
                                }
                        }
                //std::cout << "avgvalue: " << avgvalues[0] << "avgnorm: " << avgnorm << std::endl;
                
                }
        }}

        delete [] matrix;

        gapp->progress_info_close ();
        PI_DEBUG (DBG_L2, "Total number of dead_pixel: " << dead_pixel << std::endl);
        return MATH_OK;
}

---