Initial commit

GUITITAN_CHECK.pro

+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;; GUITITAN_CHECK ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;
+; This routine is called to manage automatically the dependecies
+; of different TITAN keywords and avoid bad/incompatible 
+; combinations. Is called each time a keyword is selected or
+; deselected
+
+PRO GUITITAN_CHECK, event ; event structure, passed automatically to the procedure when invoked  - INPUT
+
+WIDGET_CONTROL, event.top, get_uvalue = GuiTitan
+ntag = 11
+c_ids = BYTARR(ntag)
+s_ids = BYTARR(ntag)
+
+; let's ask cw_bgroup which is the current combination and store the results in 
+; c_ids (selected ones) and s_ids (sensitive ones)
+FOR i = 0,ntag-1 do BEGIN
+  
+  c_ids[i] = WIDGET_INFO(GuiTitan.button_ids[i], /button_set)
+  s_ids[i] = WIDGET_INFO(GuiTitan.button_ids[i], /sensitive)
+
+ENDFOR
+
+; that we adjust sensitivities properly
+CASE event.id OF
+
+  GuiTitan.button_merge: BEGIN
+
+    IF event.select THEN BEGIN
+
+      c_ids[1:*] = 0 & s_ids[1:*] = 0
+
+    ENDIF ELSE BEGIN
+
+      c_ids[0] = 0
+      s_ids[0:5] = 1 & s_ids[9:10] = 1
+
+    ENDELSE
+
+  END
+
+  GuiTitan.button_firstrun: BEGIN
+
+    IF event.select THEN BEGIN
+
+      c_ids[*] = 0 & s_ids[*] = 0
+      c_ids[2] = 1 & s_ids[2] = 1
+
+    ENDIF ELSE BEGIN
+
+      c_ids[2] = 0
+      s_ids[0:5] = 1 & s_ids[9:10] = 1
+
+    ENDELSE
+
+  END
+
+  GuiTitan.button_analysis: BEGIN
+
+    IF event.select THEN BEGIN
+
+      c_ids[0] = 0 & s_ids[0] = 0
+      c_ids[2] = 0 & s_ids[2] = 0
+
+    ENDIF ELSE BEGIN
+
+      s_ids[0] = 1
+      s_ids[2] = 1
+      c_ids[3] = 0
+
+    ENDELSE
+
+  END
+
+  GuiTitan.button_coincidenceTi: BEGIN
+
+    IF event.select THEN BEGIN
+
+      c_ids[0] = 0 & s_ids[0] = 0
+      c_ids[2] = 0 & s_ids[2] = 0
+      s_ids[6:8] = 1
+
+    ENDIF ELSE BEGIN
+
+      IF c_ids[5] EQ 0 THEN BEGIN
+
+        s_ids[0] = 1
+        s_ids[2] = 1
+        s_ids[6:8] = 0
+        c_ids[6:8] = 0
+
+      ENDIF
+
+    ENDELSE
+
+  END
+
+  GuiTitan.button_coincidenceAl: BEGIN
+
+    IF event.select THEN BEGIN
+
+      c_ids[0] = 0 & s_ids[0] = 0
+      c_ids[2] = 0 & s_ids[2] = 0
+      s_ids[6:8] = 1
+
+    ENDIF ELSE BEGIN
+
+      IF c_ids[4] EQ 0 THEN BEGIN
+
+        s_ids[0] = 1
+        s_ids[2] = 1
+        s_ids[6:8] = 0
+        c_ids[6:8] = 0
+
+      ENDIF
+
+    ENDELSE
+
+  END
+
+  ELSE:BEGIN
+
+  END
+
+ENDCASE
+
+FOR i = 0, ntag-1 DO BEGIN
+
+  WIDGET_CONTROL, GuiTitan.button_ids[i], sensitive =s_ids[i]
+  WIDGET_CONTROL, GuiTitan.button_ids[i], set_button=c_ids[i]
+
+ENDFOR
+
+END
\ No newline at end of file

GUITITAN_CLEANUP.pro

+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;; GUITITAN_CHECK ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;
+; This routine is called when the user press the X button of the
+; GUITITAN window
+
+PRO GUITITAN_CLEANUP, widgetid ; widget iden
+
+WIDGET_CONTROL, widgetid, /destroy
+
+END
\ No newline at end of file

GUITITAN_EVENT.pro

+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;; GUITITAN_EVENT ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;
+; Standard routine automatically invoked when an event occurred
+; in tlb main base of GUITITAN, if there is not an event_pro
+; specified in the definition of the objet that raised the event
+
+PRO GUITITAN_EVENT, event ; event structure, passed automatically to the procedure when invoked  - INPUT
+
+; get the uvalue of the event parent
+WIDGET_CONTROL, event.top, get_uvalue = GuiTitan
+
+; case of different ids
+CASE event.id OF
+  
+  ; the user clicked on the change dir button
+  GuiTitan.button_chgdir: BEGIN
+
+    nomedir = DIALOG_PICKFILE(/directory)
+
+    IF nomedir NE '' THEN BEGIN
+
+      WIDGET_CONTROL, GuiTitan.label.label_dir, set_value = nomedir
+
+    ENDIF
+
+  END
+  
+  ; the user clicked on the change detector buttons
+  GuiTitan.detector_button: BEGIN
+
+    CASE event.value OF
+      
+      ; case gem90
+      'gem90': BEGIN
+
+        workdir = 'n6724data\'
+        datadir = 'n6724data_original\'
+        WIDGET_CONTROL, GuiTitan.label.label_coinc_time0, set_value = '0'
+        WIDGET_CONTROL, GuiTitan.label.label_coinc_time1, set_value = '14000'
+        WIDGET_CONTROL, GuiTitan.label.label_coinc_time0, /sensitive
+        WIDGET_CONTROL, GuiTitan.label.label_coinc_time1, /sensitive
+        WIDGET_CONTROL, GuiTitan.label.label_Ge_domain1, set_value = '32000'
+        WIDGET_CONTROL, GuiTitan.label.label_association_limit, set_value = '3'
+
+      END
+       
+      ; case gem150
+      'gem150': BEGIN
+
+        workdir = 'camacdata\'
+        datadir = 'camacdata_original\'
+        WIDGET_CONTROL, GuiTitan.label.label_coinc_time0, set_value = ''
+        WIDGET_CONTROL, GuiTitan.label.label_coinc_time1, set_value = ''
+        WIDGET_CONTROL, GuiTitan.label.label_coinc_time0, sensitive = 0
+        WIDGET_CONTROL, GuiTitan.label.label_coinc_time1, sensitive = 0
+        WIDGET_CONTROL, GuiTitan.label.label_Ge_domain1, set_value = '16000'
+        WIDGET_CONTROL, GuiTitan.label.label_association_limit, set_value = '2'
+
+      END
+
+    ENDCASE
+
+    WIDGET_CONTROL, GuiTitan.label.label_workdir, set_value = workdir
+    WIDGET_CONTROL, GuiTitan.label.label_datadir, set_value = datadir
+
+  END
+  
+  ; the user clicked on EXECUTE TITAN
+  GuiTitan.button_execute: BEGIN
+    
+    ; let's read all the parameters in tlb
+    nomi    = TAG_NAMES(GuiTitan.label)
+    nomivar = TAG_NAMES(GuiTitan.var)
+    nvar    = N_ELEMENTS(nomivar)
+
+    FOR i = 0, nvar-1 DO BEGIN
+
+      nlen = STRLEN(nomivar[i])
+      ii   = WHERE(STRMID(nomi, 6, nlen) EQ nomivar[i])
+
+      IF ii[0] ne -1 THEN BEGIN
+
+        s = size(ii)
+        
+        ; if label[i] is a scalar
+        IF s[1] EQ 1 THEN BEGIN
+
+          res = EXECUTE('WIDGET_CONTROL, GuiTitan.label.' + nomi[ii[0]] +', get_value=valore')
+          res = EXECUTE('GuiTitan.var.' + nomivar[i] + '=valore')
+
+        ENDIF ELSE BEGIN
+          
+          ; if label[i] is a vector
+          FOR j = 0, N_ELEMENTS(ii)-1 DO BEGIN
+
+            res   = EXECUTE('WIDGET_CONTROL, GuiTitan.label.' + nomi[ii[j]] + ', get_value=valore')
+            index = STRMID(nomi(ii(j)), STRLEN(nomi[ii[j]])-1, 1)
+            res   = EXECUTE('GuiTitan.var.' + nomivar[i] + '[' + index + ']=valore')
+
+          ENDFOR
+
+        ENDELSE
+
+      ENDIF
+
+    ENDFOR
+
+    WIDGET_CONTROL, GuiTitan.detector_button, get_value = a_det
+    IF a_det EQ 0 THEN GuiTitan.var.det = 'gem90'
+    IF a_det EQ 1 THEN GuiTitan.var.det = 'gem150'
+
+    ; let's read the keywords in option panel
+    strkeyword = ''
+    keywords   = ['/MERGE',           $
+                  '/ALLOW_ANTICOINC', $
+                  '/FIRSTRUN',        $
+                  '/ANALYSIS',        $
+                  '/COINCIDENCETI',   $
+                  '/COINCIDENCEAL',   $
+                  '/C511',            $
+                  '/C1022',           $
+                  '/CTOT',            $
+                  '/FINAL',           $
+                  '/STOP'             $
+                  ]
+
+    FOR i = 0, N_ELEMENTS(keywords)-1 DO BEGIN
+
+      isc = WIDGET_INFO(GuiTitan.button_ids[i], /button_set)
+      IF isc THEN strkeyword = strkeyword + ', ' + keywords[i]
+
+    ENDFOR
+
+    IF GuiTitan.var.nomemet EQ '' THEN BEGIN
+
+      error =  DIALOG_MESSAGE("Null analysis ID, can't execute!")
+      RETURN
+
+    ENDIF
+
+    stringa = ['DESCRIPTION:', '', '   TITAN PROCESSING']
+    WIDGET_CONTROL, GuiTitan.label_description, set_value = stringa
+    WIDGET_CONTROL, /HOURGLASS
+    res     = EXECUTE('titan' + strkeyword + ',forwidget=guititan.var')
+    stringa = ['DESCRIPTION:', '', '']
+    WIDGET_CONTROL, GuiTitan.label_description, set_value = stringa
+
+  END
+  
+  ; the user clicked on the save par. menu
+  GuiTitan.menu_save_par: BEGIN
+
+    nomi   = TAG_NAMES(GuiTitan.label)
+    nlabel = N_ELEMENTS(nomi)
+    str    = '{'
+
+    FOR i = 0, nlabel-2 DO BEGIN
+
+      str += nomi[i]+ ":" + "''" + ","
+
+    ENDFOR
+
+    str += nomi[nlabel-1] + ":" + "''" + ",det:''" + "}"
+
+    res = EXECUTE("labelsave = " + str)
+
+    FOR i = 0, nlabel-1 DO BEGIN
+
+      res = EXECUTE('WIDGET_CONTROL, guititan.label.' + nomi[i] + ', get_value=valore')
+      res = EXECUTE('labelsave.' + nomi[i] + ' = valore')
+
+    ENDFOR
+
+    WIDGET_CONTROL, GuiTitan.detector_button, get_value = a_det
+    IF a_det EQ 0 THEN labelsave.det = 'gem90'
+    IF a_det EQ 1 THEN labelsave.det = 'gem150'
+
+    nomedir = DIALOG_PICKFILE(/directory, title = 'Please select the destination of TITAN par. sav file:')
+    CD, nomedir, current=oldir
+    SAVE, labelsave, filename = 'guititan_param.sav'
+    CD, oldir
+
+  END
+  
+  ; the user clicked on the load par. menu
+  GuiTitan.menu_load_par:BEGIN
+
+    nomesav = DIALOG_PICKFILE(title = 'Please select your .sav file to be imported:')
+
+    IF nomesav NE '' THEN BEGIN
+
+      RESTORE, nomesav
+      nomilabelsave = TAG_NAMES(labelsave)
+      nomi          = TAG_NAMES(GuiTitan.label)
+      nlabel        = N_ELEMENTS(nomi)
+
+      FOR i = 0, nlabel-1 DO BEGIN
+
+        ii  = WHERE(nomilabelsave EQ nomi[i])
+        res = EXECUTE('valore = labelsave.' + nomilabelsave[ii[0]])
+        res = EXECUTE('WIDGET_CONTROL, guititan.label.' + nomi(i) + ', set_value=valore')
+
+      ENDFOR
+
+      IF labelsave.det EQ 'gem90' THEN BEGIN
+
+        WIDGET_CONTROL, GuiTitan.detector_button, set_value = 0
+
+      ENDIF
+
+      IF labelsave.det EQ 'gem150' THEN BEGIN
+
+        WIDGET_CONTROL, GuiTitan.detector_button, set_value = 1
+
+      ENDIF
+
+    ENDIF
+
+  END
+  
+  ; the user clicked on the WTITAN menu
+  GuiTitan.menu_wtitan: BEGIN
+
+    WIDGET_CONTROL, GuiTitan.label.label_dir, get_value=dir
+    WIDGET_CONTROL, GuiTitan.label.label_workdir, get_value=workdir
+    WIDGET_CONTROL, GuiTitan.label.label_nomemet, get_value=nomemet
+    WIDGET_CONTROL, GuiTitan.label.label_info, get_value=info
+
+    IF STRLEN(nomemet) NE 0 THEN BEGIN
+           
+      answer = DIALOG_MESSAGE('Do you want to upload WTITAN for ' + nomemet + info + ' ?', /question)
+  
+      IF answer EQ 'Yes' THEN WTITAN, for_gui=dir + workdir + nomemet + info + '\'
+      IF answer EQ 'No' THEN WTITAN
+
+    ENDIF ELSE BEGIN
+    
+      WTITAN
+    
+    ENDELSE
+
+  END
+  
+  ; the user clicked on the WGESPECTRUM menu
+  GuiTitan.menu_wgespectrum:BEGIN
+    
+    WIDGET_CONTROL, GuiTitan.label.label_dir, get_value = dir
+    WIDGET_CONTROL, GuiTitan.label.label_workdir, get_value = workdir
+    WIDGET_CONTROL, GuiTitan.label.label_nomemet, get_value = nomemet
+    WIDGET_CONTROL, GuiTitan.label.label_info, get_value = info
+    
+    IF STRLEN(nomemet) NE 0 THEN BEGIN
+      
+      ff = FILE_SEARCH(dir + workdir + nomemet + info, '*_spectra.sav')
+      
+      IF ff[0] NE '' THEN BEGIN
+        
+        junk   = STRSPLIT(ff[0], '\', /extract)
+        sav    = junk[N_ELEMENTS(junk)-1]      
+        answer = DIALOG_MESSAGE('Do you want to upload WGESPECTRUM for ' + nomemet + info + '?', /question)
+        
+        IF answer EQ 'Yes' THEN WGESPECTRUM, for_gui=ff[0]
+        IF answer EQ 'No' THEN WGESPECTRUM
+      
+      ENDIF ELSE BEGIN
+        
+        WGESPECTRUM
+        
+      ENDELSE
+      
+    ENDIF ELSE BEGIN
+      
+      WGESPECTRUM
+      
+    ENDELSE
+
+  END
+
+  ELSE: BEGIN ; all the other events -> description visualization in the 
+              ; specified box
+
+    ii = WHERE(TAG_NAMES(event) EQ 'ENTER')
+
+    IF ii[0] NE -1 THEN BEGIN
+
+      IF event.enter THEN BEGIN
+
+        WIDGET_CONTROL, event.id, get_uvalue = stringa23
+        WIDGET_CONTROL, event.id, get_value = nome
+        stringa = ['DESCRIPTION: {' + nome + '}', stringa23]
+
+      ENDIF ELSE BEGIN
+
+        stringa = ['DESCRIPTION:', '', '']
+
+      ENDELSE
+
+    ENDIF ELSE BEGIN
+
+      stringa = ['DESCRIPTION:', '', '']
+
+    ENDELSE
+
+    WIDGET_CONTROL, GuiTitan.label_description, set_value = stringa
+
+  END
+
+ENDCASE
+
+; set the uvalue of the event parent
+WIDGET_CONTROL, event.top, set_uvalue = GuiTitan
+
+END

TITAN/ANALYSIS_INFO.pro

+FUNCTION ANALYSIS_INFO
+
+;;;;;;;;;;;;;;;;;;;;
+;; TITAN KEYWORDS ;;
+;;;;;;;;;;;;;;;;;;;;
+; merge               -> merge .dat files to idl binary files and calculate realtime/livetime
+; firstrun            -> first execution keyword
+; analysis            -> analysis routine keyword
+; allow_anticoinc     -> allow anticoincidence counts on H2D
+; coincidenceTI       -> coincidence routine on Ti range
+; coincidenceAL       -> coincidence routine on Al range
+; c511                -> 511 keV coincidence
+; c1022               -> 1022 keV coincidence
+; ctot                -> arbitrary coincidence of NaI channel range
+; save                -> save coincidence data
+; final               -> execute final normalization of Ti counts to c/n ratio
+; stop                -> stop at the end of the program
+                                                                              
+; DIRECTORIES AND .DAT FILE LIST
+dir     = 'E:\Lavoro\StageOATO\IDLWorkspace84\IDL_Spectrum_Analysis\'  ; general working directory
+nomemet = 'Ejby'                                                       ; name of the meteorite
+info    = ''                                                           ; other info for this analysis
+
+det = 'gem150'                                                         ; detector of this measurement (gem150 or gem90)
+
+IF det EQ 'gem150' THEN BEGIN
+
+  workdir       = 'camacdata\'                                         ; working directory for gem150 analysis of this meteorite
+  datadir       = 'camacdata_original\'                                ; .dat directory for gem150 analysis of this meteorite
+  coinc_time    = 0.                                                   ; not used for GEM150
+
+ENDIF
+
+IF det EQ 'gem90' THEN BEGIN
+
+  workdir       = 'n6724data\'                                          ; working directory for gem90 analysis of this meteorite
+  datadir       = 'n6724data_original\'                                 ; .dat directory for gem90 analysis of this meteorite
+  coinc_time    = [0,14000.]                                            ; min and max time (in 10ns units) for gem90 coincidence
+
+ENDIF
+
+; GENERAL PARAMETERS FOR ANALYSIS
+winwidth = 30                                                          ; half width of window used for fitting (total width= 2*winwidth+1)
+sigmamin = 5                                                           ; minimum sigma for line detection
+width = 100                                                            ; window width for noise (used for line detection)
+orderch2kev = 2                                                        ; order of ch to kev calibration
+orderch2sigma = 1                                                      ; order of ch to sigma calibration
+lim = 100                                                              ; minimum counts limits for a peak to be included in calibrations
+association_limit = 2                                                  ; limit, in channels, of distance for line association
+coinc_peaks=[1808.65, 1120.29]                                         ; peaks for 511-1022 keV peak identification on NaI spectrum
+Ge_domain = [0,16000]                                                  ; channel limits for the definition of H2D [Ge, NaI] histogram
+NaI_domain = [0,8000] 
+
+; COINCIDENCE PARAMETERS                                                
+fixbkg_TiBi = [ROUND(30),ROUND(30)]
+fixbkg_Al = [ROUND(60),ROUND(60)]
+fitmodel_TiBi = 1                                                      ; fit model fot Ti-Bi lines
+fitmodel_Al = 0.                                                       ; fit model fot Al line
+
+; PARAMETERS FOR ESTIMATION OD COINCIDENCE WINDOWS PARAMETERS
+dw511 = 1500                                                           ; lower ch limit for gaussian fit of 511 peak
+up511 = 2500                                                           ; upper ch limit for gaussian fit of 511 peak
+dw1022 = 3000                                                          ; lower ch limit for gaussian fit of 1022 peak
+up1022 = 4800                                                          ; upper ch limit for gaussian fit of 1022 peak
+estimates511 = [10, 2000, 150, 0., 0., 0.]                              ; initial estimates parameters for 511 keV peak
+estimates1022 = [3, 3800, 200, 0., 0., 0.]                             ; initial estimates parameters for 1022 keV peak
+
+Bi214 = 1155.210                                                       ; keV of Bi214 line
+Ti44 = 1157.020                                                        ; keV of Ti44 line
+Al26 = 1808.65                                                         ; keV of Al26 line
+
+; 511 KeV COINCIDENCE PARAMETERS
+nsxsigma_511 = 6                                                       ; left range explored of NaI channel for TiBi coincidence optimisation (in sigma units) for 511 keV
+ndxsigma_511 = 6                                                       ; right range explored of NaI channel for TiBi coincidence optimisation (in sigma units) for 511 keV
+steppeak_511 = 10                                                      ; step variation for the center of the coincidence window for 511 keV (channel units)
+stepdelta_511 = 10                                                     ; step variation for the half-width of the coincidence window for 511 keV (channel units)
+
+; 1022 KeV COINCIDENCE PARAMETERS
+nsxsigma_1022 = 6                                                      ; left range explored of NaI channel for TiBi coincidence optimisation (in sigma units) for 1022 keV
+ndxsigma_1022 = 6                                                      ; right range explored of NaI channel for TiBi coincidence optimisation (in sigma units) for 1022 keV
+steppeak_1022 = 10                                                     ; step variation for the center of the coincidence window for 1022 keV (channel units)
+stepdelta_1022 = 10                                                    ; step variation for the half-width of the coincidence window for 1022 keV (channel units)
+
+; TOT KeV COINCIDENCE PARAMETERS
+steppeak_tot = 10                                                      ; step variation for the center of the coincidence window for general coincidence (channel units) for TiBi
+stepdelta_tot = 10                                                     ; step variation for the half-width of the coincidence window for general coincidence (channel units) for TiBi
+
+retv = {                                                $
+        dir:dir,                                        $
+        nomemet:nomemet,                                $
+        info:info,                                      $
+        det:det,                                        $
+        workdir:workdir,                                $
+        datadir:datadir,                                $
+        coinc_time:coinc_time,                          $
+        winwidth:winwidth,                              $
+        sigmamin:sigmamin,                              $
+        width:width,                                    $
+        orderch2kev:orderch2kev,                        $
+        orderch2sigma:orderch2sigma,                    $
+        lim:lim,                                        $
+        association_limit:association_limit,            $
+        coinc_peaks:coinc_peaks,                        $
+        Ge_domain:Ge_domain,                            $
+        NaI_domain:NaI_domain,                          $
+        fixbkg_TiBi:fixbkg_TiBi,                        $
+        fixbkg_Al:fixbkg_Al,                            $
+        fitmodel_TiBi:fitmodel_TiBi,                    $
+        fitmodel_Al:fitmodel_Al,                        $
+        dw511:dw511,                                    $
+        up511:up511,                                    $
+        dw1022:dw1022,                                  $
+        up1022:up1022,                                  $
+        estimates511:estimates511,                      $
+        estimates1022:estimates1022,                    $
+        Bi214:Bi214,                                    $
+        Ti44:Ti44,                                      $
+        Al26:Al26,                                      $
+        nsxsigma_511:nsxsigma_511,                      $
+        ndxsigma_511:ndxsigma_511,                      $
+        steppeak_511:steppeak_511,                      $
+        stepdelta_511:stepdelta_511,                    $
+        nsxsigma_1022:nsxsigma_1022,                    $
+        ndxsigma_1022:ndxsigma_1022,                    $
+        steppeak_1022:steppeak_1022,                    $
+        stepdelta_1022:stepdelta_1022,                  $
+        steppeak_tot:steppeak_tot,                      $
+        stepdelta_tot:stepdelta_tot                     $
+        }
+        
+RETURN, retv
+
+END
+
+
+
+
+
+
+
+
+
+
+

TITAN/GAMMA_LIBRARY.pro

+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+; GAMMA_LIBRARY
+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+
+; creates gamma library
+
+FUNCTION GAMMA_LIBRARY, br_min = br_min
+
+; RADIONUCLIDES ARE LISTED FOLLOWING MASS-NUMBER ORDER. 
+; NUCLIDES WITH SAME MASS-NUMBER ARE IN ALPHABETICAL ORDER.
+
+nuc       = [           $
+               'Be7'    $
+             , 'Na22'   $
+;             , 'Na24'   $
+             , 'Al26'   $
+             , 'K40'    $
+             , 'K42'   $
+             , 'Sc44'   $
+             , 'Ti44'   $
+             , 'Sc46'   $
+             , 'Ca47'   $
+             , 'Sc47'   $
+;             , 'Cr48'   $
+             , 'V48'    $
+             , 'Cr51'   $
+             , 'Mn52'   $
+             , 'Mn54'   $
+             , 'Co57'   $
+;             , 'Ni57'   $
+             , 'Co58'   $
+             , 'Fe59'   $
+             , 'Co60'   $
+             , 'Cs137'  $
+             , 'Tl208'  $
+             , 'Bi212'  $
+             , 'Pb212'  $
+             , 'Bi214'  $
+             , 'Pb214'  $
+             , 'Ra224'  $ 
+             , 'Ra226'  $
+             , 'Ac228'  $
+             , 'Th228'  $
+             , 'Th234'  $
+             , 'Pa234m' $
+             , 'U235'   $
+             ]
+ 
+; FOR EACH NUCLIDE IN 'NUCLIDES' VECTOR WE DEFINE:
+;  
+; - keV_ list (gamma energy)
+; - int_ list (gamma perc. intensity)
+; - tau_      (half time)
+; - ftau_     (format of half time... hour, day, year)
+; - len_      (number of selected gammas from same nuclide)
+
+; FOR OUR PURPOSES, WE ADD ONLY:
+; 
+; - NUCLIDES WITH MEAN HALF TIME GREATER THAN DAYS (exception for nuclides in U-Th chains)
+; - GAMMA LINES WITH INTENSITY GREATER THAN 0.1% (exception for sum peak)    
+       
+; Be-7
+
+keV_Be7 = [477.6035] 
+int_Be7 = [10.44]
+tau_Be7 = 53.22
+ftau_Be7 = 'd'
+
+len_Be7 = N_ELEMENTS(keV_Be7)
+
+; Na-22
+
+keV_Na22 = [1274.537] 
+int_Na22 = [99.940]
+tau_Na22 = 2.6018 
+ftau_Na22 = 'y'
+
+len_Na22 = N_ELEMENTS(keV_Na22)
+
+;; Na-24
+;
+;keV_Na24 = [1368.626, 2754.007] 
+;int_Na24 = [99.9936, 99.855]
+;tau_Na24 = 14.997
+;ftau_Na24 = 'h'
+;
+;len_Na24 = N_ELEMENTS(keV_Na24)
+
+; Al-26
+
+keV_Al26 = [1129.67, 1808.65, 2938.]
+int_Al26 = [2.50, 99.76, 0.24]
+tau_Al26 = 7.17E+5
+ftau_Al26 = 'y'
+
+len_Al26 = N_ELEMENTS(keV_Al26)
+
+; K-40
+
+keV_K40 = [1460.820]
+int_K40 = [10.66]
+tau_K40 = 1.248E+9
+ftau_K40 = 'y'
+
+len_K40 = N_ELEMENTS(keV_K40)
+
+; K-42
+
+keV_K42 = [312.60, 1524.6]
+int_K42 = [0.336, 18.08]
+tau_K42 = 12.355
+ftau_K42 = 'h'
+
+len_K42 = N_ELEMENTS(keV_K42)
+
+; Sc-44
+
+keV_Sc44 = [1157.020, 1499.46, 2656.48]
+int_Sc44 = [99.9, 0.908, 0.112]
+tau_Sc44 = 3.97
+ftau_Sc44 = 'h'
+
+len_Sc44 = N_ELEMENTS(keV_Sc44)
+
+; Ti-44
+
+keV_Ti44 = [67.8679, 78.3234, 146.22]
+int_Ti44 = [93.0, 96.4, 0.092]
+tau_Ti44 = 59.2
+ftau_Ti44 = 'y'
+
+len_Ti44 = N_ELEMENTS(keV_Ti44)
+
+; Sc-46
+
+keV_Sc46  = [889.277, 1120.545]
+int_Sc46  = [99.984, 99.987]
+tau_Sc46  = 83.79
+ftau_Sc46 = 'd'
+
+len_Sc46 = N_ELEMENTS(keV_Sc46)
+
+; Ca-47
+
+keV_Ca47  = [489.23, 767.1, 807.86, 1297.09]
+int_Ca47  = [5.9, 0.18, 5.9, 67.]
+tau_Ca47  = 4.536
+ftau_Ca47 = 'd'
+
+len_Ca47 = N_ELEMENTS(keV_Ca47)
+
+; Sc-47
+
+keV_Sc47  = [159.381]
+int_Sc47  = [68.3]
+tau_Sc47  = 3.3492
+ftau_Sc47 = 'd'
+
+len_Sc47 = N_ELEMENTS(keV_Sc47)
+
+;; Cr-48
+;
+;keV_Cr48 = [112.31, 308.24]
+;int_Cr48 = [96.0, 100.]
+;tau_Cr48 = 21.56
+;ftau_Cr48 = 'h'
+;
+;len_Cr48 = N_ELEMENTS(keV_Cr48)
+
+; V-48
+
+keV_V48 = [802.90, 928.327, 983.525, 1312.106, 1437.52, 2240.396]
+int_V48 = [0.136, 0.783, 99.98, 98.2, 0.120, 2.333]
+tau_V48 = 15.9735
+ftau_V48 = 'd'
+
+len_V48 = N_ELEMENTS(keV_V48)
+
+; Cr-51
+
+keV_Cr51 = [320.0824]
+int_Cr51 = [9.910]
+tau_Cr51 = 27.704
+ftau_Cr51 = 'd'
+
+len_Cr51 = N_ELEMENTS(keV_Cr51)
+
+; Mn-52
+
+keV_Mn52 = [744.233, 935.544, 1434.092, 1727.53]
+int_Mn52 = [90.0, 94.5, 100.0, 0.216]
+tau_Mn52 = 5.591
+ftau_Mn52 = 'd'
+
+len_Mn52 = N_ELEMENTS(keV_Mn52)
+
+; Mn-54
+
+keV_Mn54 = [834.848]
+int_Mn54 = [99.9760]
+tau_Mn54 = 312.20
+ftau_Mn54 = 'd'
+
+len_Mn54 = N_ELEMENTS(keV_Mn54)
+
+; Co-57
+
+keV_Co57 = [14.4129, 122.06065, 136.47356, 692.41]
+int_Co57 = [9.16, 85.60, 10.68, 0.149]
+tau_Co57 = 271.74
+ftau_Co57 = 'd'
+
+len_Co57 = N_ELEMENTS(keV_Co57)
+
+;; Ni-57
+;
+;keV_Ni57 = [127.164, 1046.68, 1377.63, 1757.55, 1919.52]
+;int_Ni57 = [16.7, 0.134, 81.7, 5.75, 12.3]
+;tau_Ni57 = 35.60
+;ftau_Ni57 = 'h'
+;
+;len_Ni57 = N_ELEMENTS(keV_Ni57)
+
+; Co-58
+
+keV_Co58 = [810.7593, 863.951, 1674.725]
+int_Co58 = [99.450, 0.686, 0.517]
+tau_Co58 = 70.86
+ftau_Co58 = 'd'
+
+len_Co58 = N_ELEMENTS(keV_Co58)
+
+; Fe-59
+
+keV_Fe59 = [142.651, 192.343, 334.8, 1099.245, 1291.590]
+int_Fe59 = [1.02, 3.08, 0.270, 56.5, 43.2]
+tau_Fe59 = 44.495
+ftau_Fe59 = 'd'
+
+len_Fe59 = N_ELEMENTS(keV_Fe59)
+
+; Co-60
+ 
+keV_Co60 = [1173.228, 1332.492, 2505.692]
+int_Co60 = [99.85, 99.9826, 2.0E-6]
+tau_Co60 = 1925.28
+ftau_Co60 = 'd'
+
+len_Co60 = N_ELEMENTS(keV_Co60)
+
+; Cs-137
+
+keV_Cs137 = [31.817, 32.194, 36.304, 36.378, 37.255, 661.657]
+int_Cs137 = [1.99, 3.64, 0.348, 0.672, 0.213, 85.10]
+tau_Cs137 = 30.08
+ftau_Cs137 = 'y'
+
+len_Cs137 = N_ELEMENTS(keV_Cs137)
+
+; Tl-208
+ 
+keV_Tl208 = [72.805, 74.969, 84.45, 84.938, 87.3, 211.40, 233.36, 252.61, 277.371, 510.77, 583.187, 722.04, 763.13, 860.557, 927.6, 982.7, 1093.9, 2614.511]
+int_Tl208 = [2.01, 3.35, 0.404, 0.776, 0.283, 0.180, 0.310, 0.780, 6.6, 22.60, 85.0, 0.24, 1.79, 12.50, 0.125, 0.205, 0.430, 99.754]
+tau_Tl208 = 3.053
+ftau_Tl208 = 'm'
+
+len_Tl208 = N_ELEMENTS(keV_Tl208)
+
+; Bi-212
+
+keV_Bi212 = [39.857, 288.20, 328.03, 452.98, 727.330, 785.37, 893.408, 952.120, 1078.62, 1512.7, 1620.50]
+int_Bi212 = [1.06, 0.337, 0.125, 0.363, 6.67, 1.102, 0.378, 0.17, 0.564, 0.29, 1.47]
+tau_Bi212 = 60.55
+ftau_Bi212 = 'm'
+
+len_Bi212 = N_ELEMENTS(keV_Bi212)
+
+; Pb-212
+
+keV_Pb212 = [74.815, 77.107, 86.83, 87.349, 89.784, 115.183, 238.632, 300.087]
+int_Pb212 = [10.28, 17.1, 2.07, 3.97, 1.46, 0.596, 43.6, 3.30]
+tau_Pb212 = 10.64
+ftau_Pb212 = 'h'
+
+len_Pb212 = N_ELEMENTS(keV_Pb212)
+
+; Bi-214
+
+keV_Bi214 = [76.863, 79.29, 89.256, 89.807, 273.80, 348.92, 386.78, 388.89, 405.72, 454.79, 469.77, 609.320, 665.447, 703.11, 719.87, 752.85, 768.360, 786.35, 806.180, 821.18, 826.45, 934.056, 964.08, 1051.96, 1069.96, 1120.294, 1133.66, 1155.210, 1207.68, 1238.122, 1280.976, 1303.75, 1377.669, 1385.310, 1401.515, 1407.988, 1509.210, 1538.53, 1543.34, 1583.20, 1594.75, 1599.37, 1661.274, $
+             1684.012, 1729.595, 1764.491, 1838.36, 1847.429, 1873.16, 1896.05, 2118.514, 2204.059, 2293.38, 2447.70]
+int_Bi214 = [0.545, 0.907, 0.110, 0.210, 0.128, 0.104, 0.295, 0.402, 0.169, 0.292, 0.132, 45.49, 1.531, 0.472, 0.392, 0.128, 4.894, 0.32, 1.264, 0.161, 0.117, 3.107, 0.365, 0.313, 0.272, 14.92, 0.2512, 1.633, 0.451, 5.834, 1.434, 0.107, 3.988, 0.793, 1.330, 2.394, 2.130, 0.398, 0.303, 0.705, 0.267, 0.324, 1.047, 0.214, 2.878, 15.30, 0.350, 2.025, 0.214, 0.149, 1.160, 4.924, 0.307, 1.548]
+
+tau_Bi214 = 19.9
+ftau_Bi214 = 'm'
+
+len_Bi214 = N_ELEMENTS(keV_Bi214)
+
+; Pb-214
+
+keV_Pb214 = [53.2284, 74.815, 77.107, 86.83, 87.349, 89.784, 241.9950, 258.86, 274.80, 295.2228, 351.9321, 462.02, 480.432, 487.14, 533.66, 580.14, 785.96, 839.07]
+int_Pb214 = [1.075, 5.8, 9.7, 1.17, 2.24, 0.82, 7.251, 0.531, 0.355, 18.42, 35.60, 0.212, 0.337, 0.432, 0.181, 0.370, 1.06, 0.583]
+tau_Pb214 = 26.8
+ftau_Pb214 = 'm'
+
+len_Pb214 = N_ELEMENTS(keV_Pb214)
+
+; Ra-224
+
+keV_Ra224 = [81.069, 83.787, 240.986]
+int_Ra224 = [0.127, 0.209, 4.10]
+tau_Ra224 = 3.66
+ftau_Ra224 = 'd'
+
+
+len_Ra224 = N_ELEMENTS(keV_Ra224)
+
+; Ra-226
+
+keV_Ra226 = [81.069, 83.787, 186.211]
+int_Ra226 = [0.196, 0.323, 3.64]
+tau_Ra226 = 1600
+ftau_Ra226 = 'y'
+
+len_Ra226 = N_ELEMENTS(keV_Ra226)
+
+; Ac-228
+
+keV_Ac228 = [57.766, 89.957, 93.35, 99.509, 104.819, 105.604, 108.582, 129.065, 145.849, 153.977, 191.353, 199.407, 204.026, 209.253, 214.85, 270.245, 278.95, 321.646, 327.44, 328.000, 332.370, 338.320, 340.96, 409.462, 440.44, 463.004, 478.33, 503.823, 508.959, 523.131, 546.47, 562.500, 570.91, 572.14, 583.41, 674.750, 701.747, 707.41, 726.863, 755.315, 772.291, 782.142, 794.947, 830.486, 835.710, 840.377, 904.20, 911.204, $
+             947.982, 958.61, 964.766, 968.971, 1033.248, 1065.18, 1095.679, 1110.610, 1153.52, 1247.08, 1459.138, 1495.910, 1501.57, 1557.11, 1580.53, 1588.20, 1625.06, 1630.6270, 1638.281, 1666.523]
+int_Ac228 = [0.47, 1.9, 3.1, 1.26, 0.39, 0.74, 0.28, 2.42, 0.158, 0.722, 0.123, 0.315, 0.112, 3.89, 0.76, 3.46, 0.160, 0.226, 0.12, 2.95, 0.40, 11.27, 0.369, 1.92, 0.121, 4.40, 0.209, 0.182, 0.45, 0.103, 0.201, 0.87, 0.182, 0.150, 0.111, 2.1, 0.173, 0.155, 0.62, 1.00, 1.49, 0.485, 4.25, 0.540, 1.61, 0.91, 0.77, 25.8, 0.106, 0.28, 4.99, 15.80, 0.201, 0.132, 0.129, 0.285, 0.139, 0.50, 0.83, 0.86, 0.46, 0.178, 0.60, 3.22, 0.255, 1.51, 0.47, 0.178]
+tau_Ac228 = 6.15
+ftau_Ac228 = 'h'
+
+len_Ac228 = N_ELEMENTS(keV_Ac228)
+
+; Th-228
+
+keV_Th228 = [83.373, 131.613, 166.410, 215.983]
+int_Th228 = [1.19, 0.127, 0.101, 0.247]
+tau_Th228 = 1.9125
+ftau_Th228 = 'y'
+
+len_Th228 = N_ELEMENTS(keV_Th228)
+
+; Th-234
+
+keV_Th234 = [63.29, 92.38, 92.80, 112.81]
+int_Th234 = [3.7, 2.13, 2.10, 0.210]
+tau_Th234 = 24.10
+ftau_Th234 = 'd'
+
+len_Th234 = N_ELEMENTS(keV_Th234)
+
+; Pa-234m
+
+keV_Pa234m  = [766.38, 1001.03]
+int_Pa234m  = [0.294, 0.837]
+tau_Pa234m  = 24.10
+ftau_Pa234m = 'd'
+
+len_Pa234m = N_ELEMENTS(keV_Pa234m)
+
+; U-235
+
+keV_U235 = [72.7, 89.957, 93.35, 104.81, 105.604, 108.582, 109.19, 140.76, 143.76, 163.356, 182.62, 185.715, 194.94, 202.12, 205.316, 221.386]
+int_U235 = [0.120, 3.43, 5.54, 0.685, 1.31, 0.500, 1.66, 0.200, 10.96, 5.08, 0.39, 57.0, 0.630, 1.080, 5.02, 0.118]
+tau_U235 = 703.8E+6
+ftau_U235 = 'y'
+
+len_U235 = N_ELEMENTS(keV_U235)
+
+;---------------------------------------------------------------
+;---------------------------------------------------------------
+
+n = N_ELEMENTS(nuc)
+   
+str  = ''
+str1 = ''
+str2 = ''
+str3 = ''
+tot  = 0
+
+FOR i = 0, n-1 DO BEGIN
+  
+  str  = str  + ' ' + nuc[i] + ':keV_'  + nuc[i] + ','
+  str1 = str1 + ' ' + nuc[i] + ':int_'  + nuc[i] + ','
+  str2 = str2 + ' ' + nuc[i] + ':tau_'  + nuc[i] + ','
+  str3 = str3 + ' ' + nuc[i] + ':ftau_' + nuc[i] + ','
+  
+  res = EXECUTE('tot = tot + len_' + nuc[i])
+  
+ENDFOR
+
+str  = STRMID(str, 0, STRLEN(str)-1)
+str1 = STRMID(str1, 0, STRLEN(str1)-1)
+str2 = STRMID(str2, 0, STRLEN(str2)-1)
+str3 = STRMID(str3, 0, STRLEN(str3)-1)
+
+res = EXECUTE('keV = {'+str+'}')
+res = EXECUTE('int = {'+str1+'}')
+res = EXECUTE('tau = {'+str2+'}')
+res = EXECUTE('ftau = {'+str3+'}')
+
+listkeV = FLTARR(tot)
+listnuc = STRARR(tot)
+listint = FLTARR(tot)
+
+res = EXECUTE('listkeV[0:len_' + nuc[0] + '-1] = keV_' + nuc[0])
+res = EXECUTE('listnuc[0:len_' + nuc[0] + '-1] = nuc[0]')
+res = EXECUTE('listint[0:len_' + nuc[0] + '-1] = int_' + nuc[0])
+res = EXECUTE('inizio = len_' + nuc[0])
+
+FOR i = 1, n-1 DO BEGIN
+  
+  res = EXECUTE('listkeV[inizio:len_' + nuc[i] + '+inizio-1] = keV_'+nuc[i])
+  res = EXECUTE('listnuc[inizio:len_' + nuc[i] + '+inizio-1] = nuc[i]')
+  res = EXECUTE('listint[inizio:len_' + nuc[i] + '+inizio-1] = int_'+nuc[i])
+  res = EXECUTE('inizio = inizio + len_' + nuc[i])
+  
+ENDFOR
+
+ii = SORT(listkeV)
+listnuc = listnuc[ii]
+listkeV = listkeV[ii]
+listint = listint[ii]
+
+ret = {                 $
+       NUC:nuc,         $
+       KEV:keV,         $
+       INT:int,         $
+       TAU:tau,         $
+       FTAU:ftau,       $
+       LISTNUC:listnuc, $
+       LISTKEV:listkeV, $
+       LISTINT:listint  $
+      }
+       
+
+RETURN, ret
+
+END
\ No newline at end of file

TITAN/MPLINE_FIT.pro

+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;; MPLINE_FIT ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;
+; This procedure implement the LM algorithm provided by
+; MINPACK-1 (Markwardt IDL Library) to fit a peak on the spectrum
+; and is invoked by his wrapper TITAN_MPPEAKFIT
+
+FUNCTION MPLINE_FIT, xIn,                            $   ; X indipendent variable (e.g. channels)                       - INPUT
+                     yIn,                            $   ; Y dipendent variable (e.g. counts)                           - INPUT
+                     a,                              $   ; vector of fitted parameters                                  - OUTPUT
+                     CHISQ = chisq,                  $   ; chi square                                                   - OUTPUT
+                     ESTIMATES = est,                $   ; parameters estimates                                         - INPUT
+                     MEASURE_ERRORS = measureErrors, $   ; assumed errors on Y measure (to be gaussian -> w = 1/s^2)    - INPUT
+                     SIGMA = sigma,                  $   ; vector of uncertainties on fitted parameters                 - OUTPUT
+                     STATUS = status,                $   ; status of the fit (1 - OK, 0 - BAD)                          - OUTPUT
+                     FITA = fita,                    $   ; mask of parameters to be fitted (1 - FIT , 0 - NO FIT/FIXED) - INPUT
+                     COVAR = covar                       ; covariance matrix for parameters                             - OUTPUT
+
+COMPILE_OPT idl2
+ON_ERROR,2                   
+
+; number of points to fit
+n = N_ELEMENTS(yIn) 
+nMeas = N_ELEMENTS(measureErrors)
+
+; avoid floating point errors
+x = DOUBLE(xIn)
+y = DOUBLE(yIn)
+
+; treat bad measureErrors
+IF (nMeas GT 0) THEN BEGIN
+  
+  ii = WHERE(measureErrors EQ 0)
+  IF ii[0] NE -1 THEN measureErrors[ii] = 1.
+  
+ENDIF
+
+a  = DOUBLE(est)
+na = N_ELEMENTS(A)
+
+; compute number of lines
+nlines = na/6
+
+; parinfo structure definition
+parinfo = REPLICATE({fixed:0, limited:[0,0], limits:[0.,0.]}, N_ELEMENTS(a))
+parinfo.fixed = ~fita
+
+; set som limits for parameters
+FOR j = 0, nlines-1 DO BEGIN
+  
+  parinfo[0+6*j].limited = [1,0]
+  parinfo[0+6*j].limits  = [0.,0.]
+  parinfo[2+6*j].limited = [1,0]
+  parinfo[2+6*j].limits  = [0.,0.]
+  parinfo[4+6*j].limited = [1,1]
+  parinfo[4+6*j].limits  = [0.,1.]
+  parinfo[5+6*j].limited = [1,1]
+  parinfo[5+6*j].limits  = [0.,10.]
+  
+ENDFOR
+
+functargs = {X:x, Y:y, ERR:measureErrors}
+
+; call MPFIT
+param = MPFIT('MPLINE_FUNCT', a, FUNCTARGS = functargs, PERROR = sigma, BESTNORM = bestnorm, STATUS = status, $
+               COVAR = covar, AUTODERIVATIVE = 1, DOF = dof, PARINFO = parinfo, /QUIET)
+               
+IF dof LE 0 THEN dof=1
+
+; reducing status information to 0 (ok) or 1 (not ok)
+stat1 = ABS([-18,-16,-15,-14,-13,-12,-11,-10,-9,-8,-7,-6,-5,-4,-3,-2,-1,0]-status)
+stat0 = ABS([1,2,3,4,5,6,7,8,9]-status)
+IF MIN(stat1) EQ 0 THEN status = 1
+IF MIN(stat0) EQ 0 THEN status = 0
+
+; if ok
+IF status EQ 0 THEN BEGIN
+  
+  chisq = bestnorm/dof
+  sigma = FLOAT(sigma)*sqrt(chisq)
+  a     = FLOAT(param)
+  covar = FLOAT(covar)*chisq
+  
+ENDIF ELSE BEGIN ; else, allocate zeros to return
+  
+  sigma    = FLTARR(na)
+  covar    = FLTARR(na,na)
+  chisq    = 0. 
+  
+ENDELSE
+
+; computing the fitted function
+yfit = MPLINE_YFIT(x,a)
+
+RETURN, yfit
+
+END

TITAN/MPLINE_FUNCT.pro

+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;; MPLINE_FUNCT ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;
+; This is the must user-supplied function that define the z 
+; residuals to be called by MPFIT fitting procedure in 
+; order to define the chi2 to be minimized in the
+; LM algorithm
+
+FUNCTION MPLINE_FUNCT, A,        $  ; vector of parameters                             - INPUT
+                       X = x,    $  ; vector of indipendent variables (e.g. channels)  - INPUT
+                       Y = y,    $  ; vector of dipendent variables (e.g. counts)      - INPUT
+                       ERR = err    ; error for dipendent variable measurement         - INPUT
+
+COMPILE_OPT idl2, hidden
+ON_ERROR,2
+
+na = N_ELEMENTS(A)
+nlines = na/6
+nx = N_ELEMENTS(x)
+
+F = FLTARR(nx)
+
+for i=0,nlines-1 do begin
+
+  ; get Z
+  Z = (X-A[1+6*i])/(SQRT(2.)*A[2+6*i])
+
+  ; gaussian function
+  EZ = EXP(-Z^2)
+
+  ; step function
+  SF = A[3+6*i]*ERFC(Z)
+
+  ; tail function
+  T1 = (X-A[1+6*i])/A[5+6*i] + (A[2+6*i]/(SQRT(2.)*A[5+6*i]))
+  TF = (SQRT(!pi)*A[0+6*i]*A[4+6*i]/(SQRT(2.)*EXP(0.25)))*EXP(T1)*ERFC(T1)
+
+  ; sum of parts
+  F  =  F + A[0+6*i]*EZ + SF + TF
+
+ENDFOR
+
+; add background
+ibg = 1 + FINDGEN(nlines)*6
+mubg = MEAN(A[ibg])
+sibg = MEAN(A[ibg+1])
+ZBG = (X-mubg)/(SQRT(2.)*sibg)
+SFBG = A[na-1]*ERFC(ZBG)
+F = F + A[na-4] + A[na-3]*X + A[na-2]*X^2 + SFBG
+
+resid = (Y - F)/ERR
+
+RETURN, resid
+
+END
\ No newline at end of file

TITAN/MPLINE_YFIT.pro

+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;; MPLINE_YFIT ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;
+; This function is called by MPLINEFIT to evaluate yfit to return
+; at the end of the fitting procedure
+
+FUNCTION MPLINE_YFIT, x, $  ; vector of indipendent variable (e.g. channels)  - INPUT
+                      a     ; vector of parameters                            - INPUT
+ 
+COMPILE_OPT idl2, hidden
+ON_ERROR,2
+
+na = N_ELEMENTS(A)
+nlines = na/6 
+nx = N_ELEMENTS(x)
+
+F = FLTARR(nx)
+
+for i=0,nlines-1 do begin
+  
+  ; get Z
+  Z = (X-A[1+6*i])/(SQRT(2.)*A[2+6*i]) 
+ 
+  ; gaussian function
+  EZ = EXP(-Z^2)
+  
+  ; step function
+  SF = A[3+6*i]*ERFC(Z)
+  
+  ; tail function
+  T1 = (X-A[1+6*i])/A[5+6*i] + (A[2+6*i]/(SQRT(2.)*A[5+6*i]))
+  TF = (SQRT(!pi)*A[0+6*i]*A[4+6*i]/(SQRT(2.)*EXP(0.25)))*EXP(T1)*ERFC(T1)
+  
+  ; sum of parts
+  F  =  F + A[0+6*i]*EZ + SF + TF
+
+ENDFOR
+
+; add background
+ibg = 1 + FINDGEN(nlines)*6
+mubg = MEAN(A[ibg])
+sibg = MEAN(A[ibg+1])
+ZBG = (X-mubg)/(SQRT(2.)*sibg)
+SFBG = A[na-1]*ERFC(ZBG)
+F = F + A[na-4] + A[na-3]*X + A[na-2]*X^2 + SFBG
+
+RETURN, F
+
+END
\ No newline at end of file

TITAN/TITAN.pro

+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;; TITAN ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;
+; Titan is an ensemble of routines for reduction and analysis of 
+; data from HPGe gamma spectrometers at the Monte dei Cappuccini
+; laboratory in Turin (Italy).
+; This is the main routine that invokes all the subroutines
+; associated with each keyword. Titan is called from GUITITAN
+; when the 'execute' button is pushed.
+
+PRO TITAN, MERGE = MERGE,                       $    ; merge .dat files to idl binary files amd calculate realtime/livetime
+           FIRSTRUN = FIRSTRUN,                 $    ; first execution keyword
+           ANALYSIS = ANALYSIS,                 $    ; analysis routine keyword
+           ALLOW_ANTICOINC = ALLOW_ANTICOINC,   $    ; allow anticoincidence counts on H2D
+           COINCIDENCETI = COINCIDENCETI,       $    ; coincidence routine on Ti range 
+           COINCIDENCEAL = COINCIDENCEAL,       $    ; coincidence routine on Al range
+           C511 = C511,                         $    ; 511 keV coincidence
+           C1022 = C1022,                       $    ; 1022 keV coincidence
+           CTOT = CTOT,                         $    ; arbitrary coincidence of NaI channel range                
+           FINAL = FINAL,                       $    ; execute final normalization of Ti counts to c/n ratio
+           STOP = STOP,                         $    ; stop at the end of the program
+           FORWIDGET = FORWIDGET                     ; keyword thas specifies if the procedure was called from GUITITAN widget
+    
+!EXCEPT=0
+
+if KEYWORD_SET(forwidget) then begin ; if Titan is invoked from GuiTitan widget
+  
+  ; retrieve the variables structure from GuiTitan
+  var = forwidget
+  
+endif else begin
+  
+  ; load the variables structure from Analysis_Info (OBSOLETE) 
+  var = ANALYSIS_INFO()
+  
+endelse
+
+; execute a routine that exchange background colors for plot procedures
+SWAPPA_BG
+
+lavdir  = var.dir + var.workdir + var.nomemet + var.info + '\'              ; working directory for gem150 analysis of this meteorite
+dirdati = var.dir + var.datadir + var.nomemet + '\'                         ; .dat directory for gem150 analysis of this meteorite
+
+; define name of list file
+IF var.det EQ 'gem150' THEN nomefilelist = 'camacfiles_' + var.nomemet + var.info + '.txt' 
+IF var.det EQ 'gem90'  THEN nomefilelist = 'n6724files_' + var.nomemet + var.info + '.txt'
+
+var = create_struct(var, 'lavdir', lavdir, 'dirdati', dirdati, 'nomefilelist', nomefilelist, 'nomefile', var.nomemet + var.info + '.sav')
+
+; go to lavdir
+CD, var.lavdir
+
+; if the user set /MERGE keyword
+IF KEYWORD_SET(MERGE) THEN BEGIN
+  
+  TITAN_MERGE, var
+  RETURN
+ 
+ENDIF 
+
+; if the user set /FIRSTRUN or /ANALYSIS keyword, lets begin
+IF KEYWORD_SET(FIRSTRUN) OR KEYWORD_SET(ANALYSIS) THEN BEGIN
+
+  ; load binary file (.sav extension)
+  filesav = FILE_SEARCH(STRMID(var.nomefile, 0, STRLEN(var.nomefile)-4) + '.sav')
+  RESTORE, filesav
+  
+  ch  = Germ_t
+  NaI = NaI_t
+  
+  ; constructs various histograms
+  fGe     = HISTOGRAM(ch, min = var.Ge_domain[0], max = var.Ge_domain[1], binsize = 1)
+  fGe[0]  = 0  
+  fNaI    = HISTOGRAM(NaI, min = var.NaI_domain[0], max = var.NaI_domain[1], binsize=1)
+  fNaI[0] = 0 
+  xGe     = FINDGEN(N_ELEMENTS(fGe)) + var.Ge_domain[0]
+  xNaI    = FINDGEN(N_ELEMENTS(fNaI)) + var.NaI_domain[0]
+  
+  ; if the user allowed anticoincidence (/ALLOW_ANTICOINC)
+  IF KEYWORD_SET(ALLOW_ANTICOINC) THEN BEGIN
+    
+    jj = WHERE(ch GE var.Ge_domain[0] AND ch LE var.Ge_domain[1] AND nai GE var.NaI_domain[0] AND nai LE var.NaI_domain[1])
+    
+  ENDIF ELSE BEGIN
+    
+    ; standard behaviour for constructing H2D
+    jj = WHERE((ch*nai NE 0) AND ch GE var.Ge_domain[0] AND ch LE var.Ge_domain[1] AND nai GE var.NaI_domain[0] AND nai LE var.NaI_domain[1])
+    
+  ENDELSE
+  
+  ; constuct H2D
+  h2d = HIST_2D(ch[jj], NaI[jj], min1 = var.Ge_domain[0], max1 = var.Ge_domain[1], min2 = var.NaI_domain[0], max2 = var.NaI_domain[1])
+  
+  dimh2d = h2d.dim
+  dimGe = dimh2d[0]
+  dimNaI = dimh2d[1]
+  
+  ; create fit and line structure 
+  fit  = TITAN_FIT_STRUCTURE()
+  line = TITAN_LINE_STRUCTURE(fit)
+  
+  ; create calibration lines structure       
+  callines = TITAN_CREA_CALLINES()
+  
+  ; create spectra structure
+  spectra = {                            $
+             VAR:var,                    $
+             H2D:h2d,                    $
+             XGE:xGe,                    $
+             FGE:fGe,                    $
+             XNAI:xNaI,                  $
+             FNAI:fNaI,                  $
+             LINE:line,                  $
+             CALLINES:callines,          $
+             TIME:time                   $
+            }
+
+  ; if this is the first run of the meteorite
+  IF KEYWORD_SET(FIRSTRUN) THEN BEGIN
+    
+    TITAN_FIRSTRUN, spectra
+    
+  ENDIF
+  
+  ; if the user has selected /ANALYSIS keyword
+  IF KEYWORD_SET(ANALYSIS) THEN BEGIN
+    
+    TITAN_ANALYSIS, spectra
+    
+  ENDIF
+  
+ENDIF
+
+; if the user has selected a /COINCIDENCE keyword
+IF KEYWORD_SET(COINCIDENCETI) OR KEYWORD_SET(COINCIDENCEAL) THEN BEGIN
+  
+  ff = FILE_SEARCH(var.nomemet + var.info + '_spectra.sav')
+  
+  IF ff[0] NE '' THEN BEGIN
+    
+    RESTORE, ff[0]
+    
+    SPECTRA.var = var
+    
+  ENDIF ELSE BEGIN
+    
+    error = DIALOG_MESSAGE('Mandatory *_spectra.sav file missing, run ANALYSIS before.', /error)
+    return
+    
+  ENDELSE
+  
+  ; construct coincidence parameters structure
+  coinc_par = TITAN_COINC_PAR(spectra)
+  
+  if coinc_par.isNaIps EQ 0 THEN RETURN
+     
+  ; coincidence routines
+  IF KEYWORD_SET(COINCIDENCETI) THEN BEGIN 
+
+    IF KEYWORD_SET(C511) THEN var_TiBi_511 = TITAN_COINC(spectra, coinc_par, 'TiBi', '511')
+     
+    IF KEYWORD_SET(C1022) THEN var_TiBi_1022 = TITAN_COINC(spectra, coinc_par, 'TiBi', '1022')
+     
+    IF KEYWORD_SET(CTOT) THEN var_TiBi_tot = TITAN_COINC(spectra, coinc_par, 'TiBi', 'TOT')
+         
+  ENDIF
+   
+  IF KEYWORD_SET(COINCIDENCEAL) THEN BEGIN
+     
+    IF KEYWORD_SET(C511) THEN var_Al_511 = TITAN_COINC(spectra, coinc_par, 'Al', '511')
+     
+    IF KEYWORD_SET(C1022) THEN var_Al_1022 = TITAN_COINC(spectra, coinc_par, 'Al', '1022')
+     
+    IF KEYWORD_SET(CTOT) THEN var_Al_tot = TITAN_COINC(spectra, coinc_par, 'Al', 'TOT')
+    
+  ENDIF
+
+ENDIF
+
+; if the user has selected the /FINAL keyword
+IF KEYWORD_SET(FINAL) THEN TITAN_FINAL
+
+; go back to the original directory
+CD, var.dir
+
+; if the user selected the /STOP keyword
+IF KEYWORD_SET(stop) THEN STOP
+
+END

TITAN/TITAN_ANALYSIS.pro

+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;; TITAN_ANALYSIS ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;
+; This is the main routine for the analysis of data, if the 
+; /FIRSTRUN keyword was already selected
+
+PRO TITAN_ANALYSIS, spectra  ; SPECTRA structure passed from main TITAN, contains all the data needed - INPUT
+
+; lets call TITAN_MODILINES procedure:
+; the user can modify found lines with an interactive window  
+foundlines    = TITAN_MODILINES(spectra.fGe, spectra.var.nomefile)
+lines         = REPLICATE(spectra.line, N_ELEMENTS(foundlines))
+nlines        = N_ELEMENTS(lines)
+lines.id      = foundlines.id
+lines.peak_ch = foundlines.peak_ch
+
+; perform fit over every found lines by calling TITAN_FITLINES function
+fittedlines = TITAN_FITLINES(spectra.xGe, spectra.fGe, lines.peak_ch, spectra.var.winwidth)
+
+; results structure assignment (structure assignment in IDL is kind of 
+; tricky thing, doing it tag by tag we avoid bad assignment and possible 
+; errors that will cause the code to break)
+lines.multi      = fittedlines.multi
+lines.domainsize = fittedlines.domainsize
+lines.peak_ch    = fittedlines.peak_ch
+lines.fit.delta  = fittedlines.fit.delta
+lines.fit.bkg    = fittedlines.fit.bkg
+lines.fit.tail   = fittedlines.fit.tail
+lines.fit.par    = fittedlines.fit.par
+lines.fit.spar   = fittedlines.fit.spar
+lines.fit.covar  = fittedlines.fit.covar
+lines.fit.xGe    = fittedlines.fit.xGe
+lines.fit.hGe    = fittedlines.fit.hGe
+lines.fit.yfit   = fittedlines.fit.yfit
+lines.fit.E      = fittedlines.fit.E
+lines.fit.sE     = fittedlines.fit.sE
+lines.fit.chisq  = fittedlines.fit.chisq
+lines.fit.status = fittedlines.fit.status
+
+; READ existing *titan.txt file
+nomefiletxt = STRMID(spectra.var.nomefile, 0, STRLEN(spectra.var.nomefile)-4)+'_titan.txt'
+READCOL, nomefiletxt, id, multi, status, t, peak_ch, A, sA, x0, sx0, sigma1, ssigma, step, sstep, $
+                      tai, stai, bet, sbet, chisq, Energy, sigEnergy, CPM, sigCPM, Kev, sigKev, $
+                      lineName, format = '(I,I,I,I,I,F,F,F,F,F,F,F,F,F,F,F,F,F,F,F,F,F,F,F,A)', /silent
+
+; iterative way to find other peaks based on previous calibration and associations
+spectra.callines.Ch = INTERPOL(x0, Kev, spectra.callines.Kev)
+
+; lets exclude bad fits
+ee = WHERE(lines.fit.E LT 0)
+IF ee[0] NE -1 THEN lines[ee].fit.E = 0.
+ee = WHERE(lines.fit.sE LT 0)
+IF ee[0] NE -1 THEN lines[ee].fit.sE = 1.
+ee = WHERE(lines.fit.status EQ 1)
+
+IF ee[0] NE -1 THEN BEGIN
+
+  lines  = lines(WHERE(lines.fit.status EQ 0))
+  nlines = nlines - N_ELEMENTS(ee)
+
+ENDIF
+
+; define new variable for calibrations and
+; save some backups for other analysis and plots
+lines2    = lines[WHERE(lines.fit.E GE spectra.var.lim)]
+callines1 = spectra.callines
+
+; perform line association (maybe this can be done in a more efficient way 
+; implementing branching ratio analysis)
+FOR i = 0, N_ELEMENTS(spectra.callines)-1 DO BEGIN
+
+  ii = WHERE(ABS(lines.peak_ch-spectra.callines[i].ch) LE spectra.var.association_limit)
+  jj = WHERE(ABS(lines2.peak_ch-callines1[i].ch) LE spectra.var.association_limit)
+  spectra.callines[i].lineid = ii[0]
+
+  IF jj[0] NE -1 THEN BEGIN
+
+    callines1[i].lineid = ii[0]
+
+  ENDIF ELSE BEGIN
+
+    callines1[i].lineid = -1
+
+  ENDELSE
+
+ENDFOR
+
+; erasing multiple associations with calibration lines library
+FOR i=0, MAX(spectra.callines.lineid) DO BEGIN
+  
+  ii = WHERE(spectra.callines.lineid EQ i)
+  
+  IF N_ELEMENTS(ii) EQ 1 THEN CONTINUE
+  
+  spectra.callines[ii[1:*]].lineid = -1
+  
+ENDFOR
+
+FOR i=0, MAX(callines1.lineid) DO BEGIN
+
+  ii = WHERE(callines1.lineid EQ i)
+
+  IF N_ELEMENTS(ii) EQ 1 THEN CONTINUE
+
+  callines1[ii[1:*]].lineid = -1
+
+ENDFOR
+
+; save other backups for plots
+ncallines     = N_ELEMENTS(spectra.callines)
+lines.KeV     = 0.
+lines.sKev    = 0.
+callines_new  = callines1[WHERE(callines1.lineid ne -1)]
+callines_plot = spectra.callines[WHERE(spectra.callines.lineid ne -1)]
+lines_new     = lines[callines_new.lineid]
+lines_plot    = lines[callines_plot.lineid]
+
+IF N_ELEMENTS(lines_new.fit) LE 1 THEN BEGIN
+  
+  error = DIALOG_MESSAGE('Not enough associated lines to perform calibrations. Try to modify user-defined limits for calibration' + $
+                         ' or input more lines in _titan.txt and more precisely.', /error)
+  RETURN
+  
+ENDIF
+
+; perform various calibration
+ch2kevcal   = TITAN_CH2KEVCAL(lines_new.fit.par[1], callines_new.keV, spectra.var.orderch2kev)
+kev2chcal   = TITAN_KEV2CHCAL(callines_new.keV, lines_new.fit.par[1], spectra.var.orderch2kev)
+ch2sigmacal = TITAN_CH2SIGMACAL(lines_new.fit.par[1], lines_new.fit.par[2], spectra.var.orderch2sigma)
+
+; defines calibrated vectors
+FOR i = 0, spectra.var.orderch2kev DO BEGIN
+
+  lines.KeV      = lines.KeV + ch2kevcal.calpar[i]*(lines.fit.par[1])^i
+  lines2.KeV     = lines2.KeV + ch2kevcal.calpar[i]*(lines2.fit.par[1])^i
+  lines_new.Kev  = lines_new.KeV + ch2kevcal.calpar[i]*(lines_new.fit.par[1])^i
+  lines_plot.Kev = lines_plot.KeV + ch2kevcal.calpar[i]*(lines_plot.fit.par[1])^i
+
+  IF i GT 0 THEN BEGIN
+
+    lines.sKev      = lines.sKev + ch2kevcal.calpar[i]*i*(lines.fit.par[1])^(i-1) * lines.fit.spar[1]
+    lines2.sKev     = lines2.sKev + ch2kevcal.calpar[i]*i*(lines2.fit.par[1])^(i-1) * lines2.fit.spar[1]
+    lines_new.sKev  = lines_new.sKev + ch2kevcal.calpar[i]*i*(lines_new.fit.par[1])^(i-1) * lines_new.fit.spar[1]
+    lines_plot.sKev = lines_plot.sKev + ch2kevcal.calpar[i]*i*(lines_plot.fit.par[1])^(i-1) * lines_plot.fit.spar[1]
+
+  ENDIF
+
+ENDFOR
+
+; lets print various output files
+
+; _titan.ps
+; ps file with plots of: full spectrum, found lines and calibration curves
+SET_PLOT, 'ps'
+nomefileps = STRMID(spectra.var.nomefile, 0, STRLEN(spectra.var.nomefile)-4) + '_titan.ps'
+DEVICE,file = nomefileps, /port, /col, BITS=8, ys=26, yo=1, xs=20, xo=0.5
+icallines = WHERE(spectra.callines.lineid NE -1)
+WL = ROUND(FINDGEN(5)*(spectra.var.Ge_domain[1]-spectra.var.Ge_domain[0])/5) + spectra.var.Ge_domain[0]
+yrange = [1, MAX(spectra.fGe)]
+!p.multi = [0, 1, 2]
+
+FOR WN = 0, N_ELEMENTS(WL)-2 DO BEGIN
+
+  PLOT, spectra.xGe[0+WL[WN]:WL[WN+1]-1], spectra.fGe[0+WL[WN]:WL[WN+1]-1], title=spectra.var.nomemet + spectra.var.info, $
+        xtitle = 'Channel', ytitle = 'Counts', yrange = yrange, xstyle=1, /ylog, xtickformat='(I)'
+
+  FOR i = 0, nlines-1 DO BEGIN
+
+    IF lines[i].peak_ch GT WL[WN] AND lines[i].peak_ch LE WL[WN+1] THEN BEGIN
+
+      totalyrange = (!y.crange(1)-!y.crange(0))
+      ycrange     = 10.^([totalyrange*1./10, totalyrange*9./10])
+      xyoutsrange = 10.^([totalyrange*0.7/10, totalyrange*9.3/10])
+      
+      IF (i mod 2) THEN OPLOT, [lines[i].peak_ch, lines[i].peak_ch], [spectra.fGe[lines[i].peak_ch], ycrange(i MOD 2)], color = CGCOLOR("red") $
+      ELSE OPLOT, [lines[i].peak_ch, lines[i].peak_ch], [ycrange(i MOD 2), spectra.fGe[lines[i].peak_ch]], color=CGCOLOR("red")
+      ii = WHERE(i EQ spectra.callines.lineid)
+      IF ii[0] ne -1 THEN charthick=3.0 ELSE charthick=1.0
+      XYOUTS, lines[i].peak_ch, xyoutsrange[i mod 2], STRTRIM(i, 2), /data, color=CGCOLOR("red"), charsize=0.8, charthick=charthick
+
+    ENDIF
+
+  ENDFOR
+
+ENDFOR
+
+!p.multi=[0,3,6]
+
+FOR i = 0, nlines-1 DO BEGIN
+  
+  maxv = max(lines[i].fit.hge[0:lines[i].domainsize-1])*1.15
+  titolo = 'L.'+STRTRIM(i, 2) + ' - Ch.' + STRTRIM(lines[i].fit.par[1], 2) + ' - Kev ' + STRTRIM(lines[i].kev, 2)
+  ii = WHERE(i EQ spectra.callines.lineid)
+  IF ii[0] NE -1 THEN titolo = titolo + ' - '+spectra.callines[ii[0]].Name
+  PLOT, lines[i].fit.xge[0:lines[i].domainsize-1], lines[i].fit.hge[0:lines[i].domainsize-1], $
+        title=titolo, xtitle='Channel', ytitle='Counts', xticks=5, xtickformat='(I)', xstyle=1, ystyle=1, yrange=[0,maxv], psym=10
+  OPLOT, [lines[i].peak_ch,lines[i].peak_ch], [!y.crange(0),!y.crange(1)],color=cgcolor("red")
+  OPLOT, lines[i].fit.xge(0:lines[i].domainsize-1), lines[i].fit.yfit(0:lines[i].domainsize-1), color=cgcolor("green")
+
+ENDFOR
+
+!p.multi = [0,1,2]
+xch   = FINDGEN(ncallines)
+chan  =  FINDGEN(MAX(lines_plot.peak_ch)-MIN(lines_plot.peak_ch)+1)+MIN(lines_plot.peak_ch)
+kev   = FLTARR(MAX(lines_plot.peak_ch)-MIN(lines_plot.peak_ch)+1)
+sigma = FLTARR(MAX(lines_plot.peak_ch)-MIN(lines_plot.peak_ch)+1)
+
+FOR i = 0, spectra.var.orderch2kev DO BEGIN
+
+  kev = kev + ch2kevcal.calpar[i]*chan^i
+
+ENDFOR
+
+FOR i = 0, spectra.var.orderch2sigma DO BEGIN
+
+  sigma = sigma + ch2sigmacal.calpar[i]*chan^i
+
+ENDFOR
+
+IF icallines[0] NE -1 THEN BEGIN
+
+  PLOT, lines_plot.fit.par[1], callines_plot.kev,psym=4, title = 'Channel/KeV calibration', $
+        xtitle = 'channel', ytitle = 'KeV', color = CGCOLOR("black")
+  OPLOT, lines_new.fit.par[1], callines_new.kev, psym = 4, color = CGCOLOR("red")
+  OPLOT, chan, kev, color=CGCOLOR("red")
+
+  PLOT, lines_plot.fit.par[1], (lines_plot.KeV - callines_plot.KeV), $
+        psym = 4, title= 'Channel/KeV calibration', xtitle = 'channel', ytitle = 'Residuals (KeV)', color = CGCOLOR("black"), $
+        yrange = [MIN(lines_plot.KeV - callines_plot.KeV) - MAX(lines_plot.sKev), MAX(lines_plot.KeV - callines_plot.KeV) + MAX(lines_plot.sKev)]
+  ERRPLOT, lines_plot.fit.par[1], lines_plot.kev - callines_plot.kev, lines_plot.sKev, color = CGCOLOR("black")
+  OPLOT, lines_new.fit.par[1], (lines_new.KeV - callines_new.KeV), psym=4, color = CGCOLOR("red")
+  ERRPLOT,lines_new.fit.par[1], lines_new.kev-callines_new.kev, lines_new.sKev, color = CGCOLOR("red")
+
+  PLOT, lines_plot.fit.par[1], lines_plot.fit.par[2], title = 'Channel/Sigma calibration', $
+        xtitle = 'channel', ytitle = 'sigma (Ch)', psym = 4, yrange = [0, MAX(lines_plot.fit.par[2])+3], color = CGCOLOR("black")
+  ERRPLOT, lines_plot.fit.par[1], lines_plot.fit.par[2], lines_plot.fit.spar[2], color = CGCOLOR("black")
+  OPLOT, lines_new.fit.par[1], lines_new.fit.par[2], psym = 4, color = CGCOLOR("red")
+  ERRPLOT, lines_new.fit.par[1], lines_new.fit.par[2], lines_new.fit.spar[2], color = CGCOLOR("red")
+  OPLOT, chan, ch2sigmacal.calpar[0] + ch2sigmacal.calpar[1]*chan, color = CGCOLOR("red")
+
+ENDIF
+
+DEVICE, /CLOSE
+SET_PLOT, 'win'
+!p.multi = 0
+
+; _titan.txt
+; poorly readble file fot the user, but easily read by Titan during the next iteration
+nomefiletxt = STRMID(spectra.var.nomefile, 0, STRLEN(spectra.var.nomefile)-4) + '_titan.txt'
+OPENW, lun, nomefiletxt, /get_lun
+PRINTF, lun, spectra.var.nomefile
+PRINTF, lun, 'Titan parameters'
+PRINTF, lun, 'Line detection:'
+PRINTF, lun, '  Min sigma for line detection   = ' + STRTRIM(FIX(spectra.var.sigmamin), 2)
+PRINTF, lun, '  Window width for noise determ. = ' + STRTRIM(FIX(spectra.var.width), 2)
+PRINTF, lun, 'Line fitting:'
+PRINTF, lun, '  Window width for fit           = ' + STRTRIM(FIX(spectra.var.winwidth), 2)
+PRINTF, lun, 'Line calibration:'
+PRINTF, lun, '  Counts Limit                   = ' + STRTRIM(spectra.var.lim, 2)
+PRINTF, lun, '  Channel Limit                  = ' + STRTRIM(spectra.var.association_limit, 2)
+
+FOR i = 0, spectra.var.orderch2kev DO PRINTF, lun, '  ' + STRTRIM(STRING(i), 2)+' order                        = ', ch2kevcal.calpar[i]
+
+PRINTF, lun, 'Times:'
+PRINTF, lun, '  RealTime                       = ' + STRTRIM(spectra.time.realtime_m, 2) + ' min.'
+PRINTF, lun, '  LiveTime                       = ' + STRTRIM(spectra.time.livetime_m, 2) + ' min.'
+PRINTF, lun, '========================================================='
+PRINTF, lun, 'FOUND LINES'
+PRINTF, lun, '     M FS T   Ch          A       s_A        x0      s_x0     sigma   s_sigma        step      s_step        tail      s_tail        beta      s_beta      chisq    E[cnts]        s_E         cpm       s_cpm       KeV     s_KeV'
+
+FOR i = 0, nlines-1 DO BEGIN
+
+  IF lines[i].fit.par[9] EQ 0 THEN BEGIN
+
+    stringa=STRING(lines[i].id, lines[i].multi, lines[i].fit.status, lines[i].fit.tail, lines[i].peak_ch, $
+                   lines[i].fit.par[0], lines[i].fit.spar[0],                                             $
+                   lines[i].fit.par[1], lines[i].fit.spar[1],                                             $
+                   lines[i].fit.par[2], lines[i].fit.spar[2],                                             $
+                   lines[i].fit.par[3], lines[i].fit.spar[3],                                             $
+                   lines[i].fit.par[4], lines[i].fit.spar[4],                                             $
+                   lines[i].fit.par[5], lines[i].fit.spar[5],                                             $
+                   lines[i].fit.chisq,                                                                    $
+                   lines[i].fit.E, lines[i].fit.sE,                                                       $
+                   lines[i].fit.E/spectra.time.livetime_m, lines[i].fit.sE/spectra.time.livetime_m,       $
+                   lines[i].KeV, lines[i].sKev,                                                           $
+                   format='(2I3,I2,I3,I6,6F10.2,6E12.2,3F11.2,2E12.2,2F10.2)')
+
+    ii = WHERE(i EQ spectra.callines.lineid)
+    IF ii[0] NE -1 THEN stringa = stringa + '  '+spectra.callines[ii[0]].Name + ' ('+STRTRIM(spectra.callines[ii[0]].Kev, 2) + ')'
+    PRINTF, lun, stringa
+
+  ENDIF ELSE BEGIN
+
+    stringa=STRING(lines[i].id, lines[i].multi, lines[i].fit.status, lines[i].fit.tail, lines[i].peak_ch, $
+                   lines[i].fit.par[0], lines[i].fit.spar[0],                                             $
+                   lines[i].fit.par[1], lines[i].fit.spar[1],                                             $
+                   lines[i].fit.par[2], lines[i].fit.spar[2],                                             $ 
+                   lines[i].fit.par[9], lines[i].fit.spar[9],                                             $
+                   lines[i].fit.par[4], lines[i].fit.spar[4],                                             $
+                   lines[i].fit.par[5], lines[i].fit.spar[5],                                             $
+                   lines[i].fit.chisq,                                                                    $
+                   lines[i].fit.E, lines[i].fit.sE,                                                       $
+                   lines[i].fit.E/spectra.time.livetime_m, lines[i].fit.sE/spectra.time.livetime_m,       $
+                   lines[i].KeV, lines[i].sKev,                                                           $
+                   format='(2I3,I2,I3,I6,6F10.2,6E12.2,3F11.2,2E12.2,2F10.2)')
+
+    ii = WHERE(i EQ spectra.callines.lineid)
+    IF ii[0] NE -1 THEN stringa = stringa + '  ' + spectra.callines[ii[0]].Name + ' (' + STRTRIM(spectra.callines[ii[0]].Kev, 2) + ')'
+    PRINTF, lun, stringa
+
+  ENDELSE
+
+ENDFOR
+
+CLOSE, lun & FREE_LUN, lun
+
+; _report.txt
+; easily readble file for user to extrapolate results, but with less informations
+nomefiletxt = STRMID(spectra.var.nomefile,0,STRLEN(spectra.var.nomefile)-4) + '_report.txt'
+OPENW, lun, nomefiletxt,/get_lun
+PRINTF, lun, ''
+PRINTF, lun, spectra.var.nomemet + spectra.var.info + ' report'
+PRINTF, lun, ''
+PRINTF, lun, SYSTIME()
+PRINTF, lun, ''
+PRINTF, lun, '================================================================================'
+PRINTF, lun, ''
+PRINTF, lun, 'ch -> keV calib.'
+FOR i = 0, spectra.var.orderch2kev DO PRINTF, lun, '  ' + STRTRIM(STRING(i), 2) + ' order             = ', STRTRIM(STRING(ch2kevcal.calpar[i]), 2) + ' ± ' + STRTRIM(STRING(ch2kevcal.sigma[i]), 2)
+PRINTF, lun, ''
+PRINTF, lun, 'keV -> ch calib.'
+FOR i = 0, spectra.var.orderch2kev DO PRINTF, lun, '  ' + STRTRIM(STRING(i), 2) + ' order             = ', STRTRIM(STRING(kev2chcal.calpar[i]), 2) + ' ± ' + STRTRIM(STRING(kev2chcal.sigma[i]), 2)
+PRINTF, lun, ''
+PRINTF, lun, 'Times:'
+PRINTF, lun, '  RealTime               = ' + STRTRIM(STRING(spectra.time.realtime_d), 2) + ' d. = ' + STRTRIM(STRING(spectra.time.realtime_h), 2) + ' h. ='
+PRINTF, lun, '                         = ' + STRTRIM(STRING(spectra.time.realtime_m), 2) + ' min. = ' + STRTRIM(STRING(spectra.time.realtime_s), 2) + ' sec.'
+PRINTF, lun, ''
+PRINTF, lun, '  LiveTime               = ' + STRTRIM(STRING(spectra.time.livetime_d), 2) + ' d. = ' + STRTRIM(STRING(spectra.time.livetime_h), 2) + ' h. = '
+PRINTF, lun, '                         = ' + STRTRIM(STRING(spectra.time.livetime_m), 2) + ' min. = ' + STRTRIM(STRING(spectra.time.livetime_s), 2) + ' sec.'
+PRINTF, lun, ''
+PRINTF, lun, '================================================================================'
+PRINTF, lun, ''
+
+FOR i = 0, nlines-1 DO BEGIN
+
+  ii = WHERE(i EQ spectra.callines.lineid)
+
+  IF ii[0] ne -1 THEN BEGIN
+
+    stringa = STRING(lines[i].id, format='(I3)') + ' - ' + spectra.callines[ii[0]].Name + ' (' + STRTRIM(spectra.callines[ii[0]].Kev, 2)+' keV, ' + STRTRIM(spectra.callines[ii[0]].intensity, 2) + ' %)'
+
+  ENDIF ELSE BEGIN
+
+    stringa = STRING(lines[i].id, format='(I3)') +' - Unidentified line'
+
+  ENDELSE
+
+  stringa1 = 'centr.   = ' + STRTRIM(STRING(lines[i].fit.par[1]), 2) + ' ± ' + STRTRIM(STRING(lines[i].fit.spar[1]), 2) + ' ch'
+  stringa2 = 'sigma    = ' + STRTRIM(STRING(lines[i].fit.par[2]), 2) + ' ± ' + STRTRIM(STRING(lines[i].fit.spar[2]), 2) + ' ch'
+  stringa3 = 'Counts   = ' + STRTRIM(STRING(lines[i].fit.E), 2) + ' ± ' + STRTRIM(STRING(lines[i].fit.sE),2) + ' cnts'
+  stringa4 = 'CPM      = ' + STRTRIM(STRING(lines[i].fit.E/spectra.time.livetime_m),2) + ' ± ' + STRTRIM(STRING(lines[i].fit.sE/spectra.time.livetime_m), 2) + ' cnts/min'
+  stringa5 = 'Cal. En. = ' + STRTRIM(STRING(lines[i].keV), 2) + ' ± ' + STRTRIM(STRING(lines[i].skeV), 2) + ' keV'
+  stringa6 = 'Chisq.   = ' + STRTRIM(STRING(lines[i].fit.chisq), 2)
+
+  IF ((i-4) MOD 6) EQ 0 THEN BEGIN
+
+    PRINTF, lun, ''
+    PRINTF, lun, ''
+
+  ENDIF
+
+  PRINTF, lun, stringa
+  PRINTF, lun, ''
+  PRINTF, lun, stringa1
+  PRINTF, lun, stringa2
+  PRINTF, lun, stringa3
+  PRINTF, lun, stringa4
+  PRINTF, lun, stringa5
+  PRINTF, lun, stringa6
+  PRINTF, lun, ''
+  PRINTF, lun, ''
+
+  IF ((i-3) MOD 6) EQ 0 AND i NE 3 THEN BEGIN
+
+    PRINTF, lun, ''
+    PRINTF, lun, ''
+
+  ENDIF
+
+ENDFOR
+
+CLOSE, lun & FREE_LUN, lun
+
+kevGe = FLTARR(N_ELEMENTS(spectra.xGe))
+
+FOR i = 0, spectra.var.orderch2kev DO BEGIN
+
+  kevGe = kevGe + ch2kevcal.calpar[i]*spectra.xGe^i
+
+ENDFOR
+
+; lets save the spectra structure
+spectra = {                           $
+           VAR:spectra.var,           $
+           H2D:spectra.h2d,           $
+           XGE:spectra.xGe,           $
+           KEVGE:keVGe,               $
+           FGE:spectra.fGe,           $
+           XNAI:spectra.xNaI,         $
+           FNAI:spectra.fNaI,         $
+           LINES:lines,               $
+           TIME:spectra.time,         $
+           CH2KEVCAL:ch2kevcal,       $
+           KEV2CHCAL:kev2chcal,       $
+           CH2SIGMACAL:ch2sigmacal    $
+          }
+
+nomefilesav = STRMID(spectra.var.nomefile, 0, STRLEN(spectra.var.nomefile)-4) + '_spectra.sav'
+SAVE, spectra, filename = nomefilesav
+
+information = DIALOG_MESSAGE('Successful analysis. Please check results in report.txt and titan.ps files.', /information)
+
+END
\ No newline at end of file

TITAN/TITAN_CH2KEVCAL.pro

+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;; TITAN_CH2KEVCAL ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;
+; Simple function that performs ch->keV calibration for Titan
+  
+FUNCTION TITAN_CH2KEVCAL, linepos, $  ; position of lines to be used in calibration (centroid of the peak - channels)  - INPUT
+                          kev,     $  ; energy in keV of peak corresponding to linepos                                 - INPUT
+                          order       ; polynomial order of the calibration                                            - INPUT
+
+calpar    = POLY_FIT(linepos, kev, order, sigma = sigma, yfit = ycalfit)
+
+ch2kevcal = {                 $
+             ORDER:order,     $
+             CALPAR:calpar,   $
+             YCALFIT:ycalfit, $
+             SIGMA:sigma      $
+             }
+
+RETURN, ch2kevcal
+  
+END
\ No newline at end of file

TITAN/TITAN_CH2SIGMACAL.pro

+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;; TITAN_CH2SIGMACAL ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;
+; Simple function that performs ch->sigma calibration for Titan.
+
+FUNCTION TITAN_CH2SIGMACAL, linespos,   $  ; position of lines to be used in calibration (centroid of the peak - channels)  - INPUT
+                            linessigma, $  ; sigma in channels of peak corresponding to linepos                             - INPUT
+                            order          ; polynomial order of the calibration
+
+calpar      = POLY_FIT(linespos, linessigma, order, sigma = sigma, yfit = ycalfit)
+
+ch2sigmacal = {                 $
+               ORDER:order,     $
+               CALPAR:calpar,   $
+               YCALFIT:ycalfit, $
+               SIGMA:sigma      $
+              }
+
+RETURN, ch2sigmacal
+
+END
\ No newline at end of file

TITAN/TITAN_COINC.pro

+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;; TITAN_COINC ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+
+; this function performs different coincidence cycle by invoking
+; TITAN_COINCCYCLE_TIBI or _AL
+
+FUNCTION TITAN_COINC, spectra,    $  ; spectra structure, passed from Titan main, that contain evertything is needed - INPUT
+                      coinc_par,  $  ; coincidence parameters structure                                              - INPUT
+                      TiBiAl,     $  ; 'TiBi' or 'Al'                                                                - INPUT
+                      whichcoinc     ; '511', '1022' or 'tot'                                                        - INPUT
+
+; the user selected the TiBi coincidence keyword
+IF TiBIAl EQ 'TiBi' THEN BEGIN
+  
+  ; the user selected the 511 keV TiBI coincidence
+  IF whichcoinc EQ '511' THEN BEGIN
+    
+    dim_m_511    = FLOOR((spectra.var.nsxsigma_511+spectra.var.ndxsigma_511)* coinc_par.nomsigma_511/spectra.var.steppeak_511)+1
+    dim_w_511    = FLOOR((spectra.var.nsxsigma_511+spectra.var.ndxsigma_511)* coinc_par.nomsigma_511/(2*spectra.var.stepdelta_511))+2
+    mv_511       = FINDGEN(dim_m_511)
+    wv_511       = FINDGEN(dim_w_511)
+    mean_511     = coinc_par.nompeak_511 - spectra.var.nsxsigma_511*coinc_par.nomsigma_511 + spectra.var.steppeak_511*mv_511
+    width_511    = spectra.var.stepdelta_511*wv_511 + spectra.var.stepdelta_511
+    limitinf_511 = coinc_par.nompeak_511 - spectra.var.nsxsigma_511*coinc_par.nomsigma_511
+    limitsup_511 = coinc_par.nompeak_511 + spectra.var.ndxsigma_511*coinc_par.nomsigma_511
+
+    coinc_info = {                                         $
+                  COINC:'511',                             $
+                  DIM_M:dim_m_511,                         $
+                  DIM_W:dim_w_511,                         $
+                  MEAN:mean_511,                           $
+                  WIDTH:width_511,                         $
+                  LIMITINF:limitinf_511,                   $
+                  LIMITSUP:limitsup_511,                   $
+                  MINI_TIBI:coinc_par.mini_TiBi,           $
+                  MAXI_TIBI:coinc_par.maxi_TiBi,           $
+                  FIXPAR_TIBI:coinc_par.fixpar_TiBi,       $
+                  FIXBKG_TIBI:spectra.var.fixbkg_TiBi,     $
+                  FITMODEL_TIBI:spectra.var.fitmodel_TiBi, $
+                  NOMPEAK_511:coinc_par.nompeak_511,       $
+                  NOMPEAK_1022:coinc_par.nompeak_1022,     $
+                  NOMSIGMA_511:coinc_par.nomsigma_511,     $
+                  NOMSIGMA_1022:coinc_par.nomsigma_1022,   $
+                  CH_TI44:coinc_par.ch_Ti44,               $
+                  CH_BI214:coinc_par.ch_Bi214,             $
+                  SIGMA_TI44:coinc_par.sigma_Ti44,         $
+                  SIGMA_BI214:coinc_par.sigma_Bi214        $
+                 }
+
+    var_TiBi_511 = TITAN_COINCCYCLE_TIBI(spectra, coinc_info)
+
+    addinfo = STRTRIM(STRING(FIX(spectra.var.nsxsigma_511)), 1)+ '-' + STRTRIM(STRING(FIX(spectra.var.ndxsigma_511)), 1) + '-' + STRTRIM(STRING(FIX(spectra.var.steppeak_511)), 1) + '-' +  $
+      STRTRIM(STRING(FIX(spectra.var.stepdelta_511)), 1) + '_' + 'bkg' + STRTRIM(STRING(FIX(spectra.var.fixbkg_TiBi[0])), 1) + '-' + STRTRIM(STRING(FIX(spectra.var.fixbkg_TiBi[1])), 1)
+    nomefilesav = STRMID(spectra.var.nomefile, 0, STRLEN(spectra.var.nomefile)-4)+ '_TiBi_' + addinfo + '_' + coinc_info.coinc + '.sav'
+
+    SAVE, var_TiBi_511, filename = nomefilesav
+    
+    RETURN, var_TiBi_511
+  
+  ENDIF
+  
+  ; the user selected the 1022 keV TiBi coincidence
+  IF whichcoinc EQ '1022' THEN BEGIN
+    
+    dim_m_1022    = FLOOR((spectra.var.nsxsigma_1022+spectra.var.ndxsigma_1022)*coinc_par.nomsigma_1022/spectra.var.steppeak_1022)+1
+    dim_w_1022    = FLOOR((spectra.var.nsxsigma_1022+spectra.var.ndxsigma_1022)*coinc_par.nomsigma_1022/(2*spectra.var.stepdelta_1022))+2
+    mv_1022       = FINDGEN(dim_m_1022)
+    wv_1022       = FINDGEN(dim_w_1022)
+    mean_1022     = coinc_par.nompeak_1022 - spectra.var.nsxsigma_1022*coinc_par.nomsigma_1022 + spectra.var.steppeak_1022*mv_1022
+    width_1022    = spectra.var.stepdelta_1022*wv_1022 + spectra.var.stepdelta_1022
+    limitinf_1022 = coinc_par.nompeak_1022 - spectra.var.nsxsigma_1022*coinc_par.nomsigma_1022
+    limitsup_1022 = coinc_par.nompeak_1022 + spectra.var.ndxsigma_1022*coinc_par.nomsigma_1022
+
+    coinc_info = {                                         $
+                  COINC:'1022',                            $
+                  DIM_M:dim_m_1022,                        $
+                  DIM_W:dim_w_1022,                        $
+                  MEAN:mean_1022,                          $
+                  WIDTH:width_1022,                        $
+                  LIMITINF:limitinf_1022,                  $
+                  LIMITSUP:limitsup_1022,                  $
+                  MINI_TIBI:coinc_par.mini_TiBi,           $
+                  MAXI_TIBI:coinc_par.maxi_TiBi,           $
+                  FIXPAR_TIBI:coinc_par.fixpar_TiBi,       $
+                  FIXBKG_TIBI:spectra.var.fixbkg_TiBi,     $
+                  FITMODEL_TIBI:spectra.var.fitmodel_TiBi, $
+                  NOMPEAK_511:coinc_par.nompeak_511,       $
+                  NOMPEAK_1022:coinc_par.nompeak_1022,     $
+                  NOMSIGMA_511:coinc_par.nomsigma_511,     $
+                  NOMSIGMA_1022:coinc_par.nomsigma_1022,   $
+                  CH_TI44:coinc_par.ch_Ti44,               $
+                  CH_BI214:coinc_par.ch_Bi214,             $
+                  SIGMA_TI44:coinc_par.sigma_Ti44,         $
+                  SIGMA_BI214:coinc_par.sigma_Bi214        $
+                 }
+ 
+    var_TiBi_1022 = TITAN_COINCCYCLE_TIBI(spectra, coinc_info)
+
+    addinfo = STRTRIM(STRING(FIX(spectra.var.nsxsigma_1022)), 1)+ '-' + STRTRIM(STRING(FIX(spectra.var.ndxsigma_1022)), 1) + '-' + STRTRIM(STRING(FIX(spectra.var.steppeak_1022)), 1) + '-' +  $
+              STRTRIM(STRING(FIX(spectra.var.stepdelta_1022)), 1) + '_' + 'bkg' + STRTRIM(STRING(FIX(spectra.var.fixbkg_TiBi[0])), 1) + '-' + STRTRIM(STRING(FIX(spectra.var.fixbkg_TiBi[1])), 1)
+    nomefilesav = STRMID(spectra.var.nomefile, 0, STRLEN(spectra.var.nomefile)-4)+ '_TiBi_' + addinfo + '_' + coinc_info.coinc + '.sav'
+
+    SAVE, var_TiBi_1022, filename = nomefilesav
+    
+    RETURN, var_TiBi_1022
+    
+  ENDIF
+  
+  ; the user selected the TOT TiBi coincidence
+  IF whichcoinc eq 'TOT' THEN BEGIN
+    
+    dim_m_tot    = FLOOR((spectra.var.NaI_domain[1]-spectra.var.NaI_domain[0])/spectra.var.steppeak_tot)+1
+    dim_w_tot    = FLOOR((spectra.var.NaI_domain[1]-spectra.var.NaI_domain[0])/(2*spectra.var.stepdelta_tot))+2
+    mv_tot       = FINDGEN(dim_m_tot)
+    wv_tot       = FINDGEN(dim_w_tot)
+    mean_tot     = spectra.var.steppeak_tot*mv_tot
+    width_tot    = spectra.var.stepdelta_tot*wv_tot + spectra.var.stepdelta_tot
+    limitinf_tot = 0
+    limitsup_tot = spectra.var.NaI_domain[1]
+
+    coinc_info = {                             $
+                  COINC:'tot',                             $
+                  DIM_M:dim_m_tot,                         $
+                  DIM_W:dim_w_tot,                         $
+                  MEAN:mean_tot,                           $
+                  WIDTH:width_tot,                         $
+                  LIMITINF:limitinf_tot,                   $
+                  LIMITSUP:limitsup_tot,                   $
+                  MINI_TIBI:coinc_par.mini_TiBi,           $
+                  MAXI_TIBI:coinc_par.maxi_TiBi,           $
+                  FIXPAR_TIBI:coinc_par.fixpar_TiBi,       $
+                  FIXBKG_TIBI:spectra.var.fixbkg_TiBi,     $
+                  FITMODEL_TIBI:spectra.var.fitmodel_TiBi, $
+                  NOMPEAK_511:coinc_par.nompeak_511,       $
+                  NOMPEAK_1022:coinc_par.nompeak_1022,     $
+                  NOMSIGMA_511:coinc_par.nomsigma_511,     $
+                  NOMSIGMA_1022:coinc_par.nomsigma_1022,   $                        
+                  CH_TI44:coinc_par.ch_Ti44,               $
+                  CH_BI214:coinc_par.ch_Bi214,             $
+                  SIGMA_TI44:coinc_par.sigma_Ti44,         $
+                  SIGMA_BI214:coinc_par.sigma_Bi214        $
+                 }
+
+    var_TiBi_tot = TITAN_COINCCYCLE_TIBI(spectra, coinc_info)
+
+    addinfo = STRTRIM(STRING(FIX(spectra.var.steppeak_tot)), 1) + '-' + STRTRIM(STRING(FIX(spectra.var.stepdelta_tot)), 1) + '_' + 'bkg' + STRTRIM(STRING(FIX(spectra.var.fixbkg_TiBi[0])), 1) + '-' + $
+              STRTRIM(STRING(FIX(spectra.var.fixbkg_TiBi[1])), 1)
+    nomefilesav = STRMID(spectra.var.nomefile, 0, STRLEN(spectra.var.nomefile)-4)+ '_TiBi_' + addinfo + '_' + coinc_info.coinc + '.sav'
+
+    SAVE, var_TiBi_tot, filename = nomefilesav
+    
+    RETURN, var_TiBi_tot
+  
+  ENDIF
+  
+  
+ENDIF
+
+; the user selected the Al coincidence keyword
+IF TiBIAl EQ 'Al' THEN BEGIN
+  
+  ; the user selected the 511 keV Al coincidence
+  IF whichcoinc EQ '511' THEN BEGIN
+
+    dim_m_511    = FLOOR((spectra.var.nsxsigma_511+spectra.var.ndxsigma_511)*coinc_par.nomsigma_511/spectra.var.steppeak_511)+1
+    dim_w_511    = FLOOR((spectra.var.nsxsigma_511+spectra.var.ndxsigma_511)*coinc_par.nomsigma_511/(2*spectra.var.stepdelta_511))+2
+    mv_511       = FINDGEN(dim_m_511)
+    wv_511       = FINDGEN(dim_w_511)
+    mean_511     = coinc_par.nompeak_511 - spectra.var.nsxsigma_511*coinc_par.nomsigma_511 + spectra.var.steppeak_511*mv_511
+    width_511    = spectra.var.stepdelta_511*wv_511 + spectra.var.stepdelta_511
+    limitinf_511 = coinc_par.nompeak_511 - spectra.var.nsxsigma_511*coinc_par.nomsigma_511
+    limitsup_511 = coinc_par.nompeak_511 + spectra.var.ndxsigma_511*coinc_par.nomsigma_511
+
+    coinc_info = {                                         $
+                  COINC:'511',                             $
+                  DIM_M:dim_m_511,                         $
+                  DIM_W:dim_w_511,                         $
+                  MEAN:mean_511,                           $
+                  WIDTH:width_511,                         $
+                  LIMITINF:limitinf_511,                   $
+                  LIMITSUP:limitsup_511,                   $
+                  MINI_AL:coinc_par.mini_Al,               $
+                  MAXI_AL:coinc_par.maxi_Al,               $
+                  FIXPAR_AL:coinc_par.fixpar_Al,           $
+                  FIXBKG_AL:spectra.var.fixbkg_Al,         $
+                  FITMODEL_AL:spectra.var.fitmodel_Al,     $
+                  NOMPEAK_511:coinc_par.nompeak_511,       $
+                  NOMPEAK_1022:coinc_par.nompeak_1022,     $
+                  NOMSIGMA_511:coinc_par.nomsigma_511,     $
+                  NOMSIGMA_1022:coinc_par.nomsigma_1022,   $
+                  CH_AL26:coinc_par.ch_Al26,               $
+                  SIGMA_AL26:coinc_par.sigma_Al26,         $
+                  E_Al26_N:coinc_par.E_Al26_N,             $
+                  sE_Al26_N:coinc_par.sE_Al26_N            $
+                 }
+
+    var_Al_511 = TITAN_COINCCYCLE_AL(spectra, coinc_info)
+
+    addinfo = STRTRIM(STRING(FIX(spectra.var.nsxsigma_511)), 1)+ '-' + STRTRIM(STRING(FIX(spectra.var.ndxsigma_511)), 1) + '-' + STRTRIM(STRING(FIX(spectra.var.steppeak_511)), 1) + '-' +  $
+              STRTRIM(STRING(FIX(spectra.var.stepdelta_511)), 1) + '_' + 'bkg' + STRTRIM(STRING(FIX(spectra.var.fixbkg_Al[0])), 1) + '-' + STRTRIM(STRING(FIX(spectra.var.fixbkg_Al[1])), 1)
+    nomefilesav = STRMID(spectra.var.nomefile, 0, STRLEN(spectra.var.nomefile)-4)+ '_Al_' + addinfo + '_' + coinc_info.coinc + '.sav'
+
+    SAVE, var_Al_511, filename = nomefilesav
+    
+    RETURN, var_Al_511
+
+  ENDIF
+  
+  ; the user selected the 1022 keV Al coincidence
+  IF whichcoinc eq '1022' THEN BEGIN
+    
+    dim_m_1022    = FLOOR((spectra.var.nsxsigma_1022+spectra.var.ndxsigma_1022)*coinc_par.nomsigma_1022/spectra.var.steppeak_1022)+1
+    dim_w_1022    = FLOOR((spectra.var.nsxsigma_1022+spectra.var.ndxsigma_1022)*coinc_par.nomsigma_1022/(2*spectra.var.stepdelta_1022))+2
+    mv_1022       = FINDGEN(dim_m_1022)
+    wv_1022       = FINDGEN(dim_w_1022)
+    mean_1022     = coinc_par.nompeak_1022 - spectra.var.nsxsigma_1022*coinc_par.nomsigma_1022 + spectra.var.steppeak_1022*mv_1022
+    width_1022    = spectra.var.stepdelta_1022*wv_1022 + spectra.var.stepdelta_1022
+    limitinf_1022 = coinc_par.nompeak_1022 - spectra.var.nsxsigma_1022*coinc_par.nomsigma_1022
+    limitsup_1022 = coinc_par.nompeak_1022 + spectra.var.ndxsigma_1022*coinc_par.nomsigma_1022
+
+    coinc_info = {                                         $
+                  COINC:'1022',                            $
+                  DIM_M:dim_m_1022,                        $
+                  DIM_W:dim_w_1022,                        $
+                  MEAN:mean_1022,                          $
+                  WIDTH:width_1022,                        $
+                  LIMITINF:limitinf_1022,                  $
+                  LIMITSUP:limitsup_1022,                  $
+                  MINI_AL:coinc_par.mini_Al,               $
+                  MAXI_AL:coinc_par.maxi_Al,               $
+                  FIXPAR_AL:coinc_par.fixpar_Al,           $
+                  FIXBKG_AL:spectra.var.fixbkg_Al,         $
+                  FITMODEL_AL:spectra.var.fitmodel_Al,     $
+                  NOMPEAK_511:coinc_par.nompeak_511,       $
+                  NOMPEAK_1022:coinc_par.nompeak_1022,     $
+                  NOMSIGMA_511:coinc_par.nomsigma_511,     $
+                  NOMSIGMA_1022:coinc_par.nomsigma_1022,   $
+                  CH_AL26:coinc_par.ch_Al26,               $
+                  SIGMA_AL26:coinc_par.sigma_Al26,         $
+                  E_Al26_N:coinc_par.E_Al26_N,             $
+                  sE_Al26_N:coinc_par.sE_Al26_N            $
+                 }
+
+    var_Al_1022 = TITAN_COINCCYCLE_AL(spectra, coinc_info)
+
+    addinfo = STRTRIM(STRING(FIX(spectra.var.nsxsigma_1022)), 1)+ '-' + STRTRIM(STRING(FIX(spectra.var.ndxsigma_1022)), 1) + '-' + STRTRIM(STRING(FIX(spectra.var.steppeak_1022)), 1) + '-' +  $
+              STRTRIM(STRING(FIX(spectra.var.stepdelta_1022)), 1) + '_' + 'bkg' + STRTRIM(STRING(FIX(spectra.var.fixbkg_Al[0])), 1) + '-' + STRTRIM(STRING(FIX(spectra.var.fixbkg_Al[1])), 1)
+    nomefilesav = STRMID(spectra.var.nomefile, 0, STRLEN(spectra.var.nomefile)-4)+ '_Al_' + addinfo + '_' + coinc_info.coinc + '.sav'
+
+    SAVE, var_Al_1022, filename = nomefilesav
+    
+    RETURN, var_Al_1022
+  
+  ENDIF  
+  
+  ; the user selected the TOT Al coincidence
+  IF whichcoinc EQ 'TOT' THEN BEGIN
+    
+    dim_m_tot    = FLOOR((spectra.var.NaI_domain[1]-spectra.var.NaI_domain[0])/spectra.var.steppeak_tot)+1
+    dim_w_tot    = FLOOR((spectra.var.NaI_domain[1]-spectra.var.NaI_domain[0])/(2*spectra.var.stepdelta_tot))+2
+    mv_tot       = FINDGEN(dim_m_tot)
+    wv_tot       = FINDGEN(dim_w_tot)
+    mean_tot     = spectra.var.steppeak_tot*mv_tot
+    width_tot    = spectra.var.stepdelta_tot*wv_tot + spectra.var.stepdelta_tot
+    limitinf_tot = 0
+    limitsup_tot = spectra.var.NaI_domain[1]
+
+    coinc_info = {                                         $
+                  COINC:'tot',                             $
+                  DIM_M:dim_m_tot,                         $
+                  DIM_W:dim_w_tot,                         $
+                  MEAN:mean_tot,                           $
+                  WIDTH:width_tot,                         $
+                  LIMITINF:limitinf_tot,                   $
+                  LIMITSUP:limitsup_tot,                   $
+                  MINI_AL:coinc_par.mini_Al,               $
+                  MAXI_AL:coinc_par.maxi_Al,               $
+                  FIXPAR_AL:coinc_par.fixpar_Al,           $
+                  FIXBKG_AL:spectra.var.fixbkg_Al,         $
+                  FITMODEL_AL:spectra.var.fitmodel_Al,     $
+                  NOMPEAK_511:coinc_par.nompeak_511,       $
+                  NOMPEAK_1022:coinc_par.nompeak_1022,     $
+                  NOMSIGMA_511:coinc_par.nomsigma_511,     $
+                  NOMSIGMA_1022:coinc_par.nomsigma_1022,   $
+                  CH_AL26:coinc_par.ch_Al26,               $
+                  SIGMA_AL26:coinc_par.sigma_Al26,         $
+                  E_Al26_N:coinc_par.E_Al26_N,             $
+                  sE_Al26_N:coinc_par.sE_Al26_N            $
+                 }
+
+    var_Al_tot = TITAN_COINCCYCLE_AL(spectra, coinc_info)
+
+    addinfo = STRTRIM(STRING(FIX(spectra.var.steppeak_tot)), 1) + '-' + STRTRIM(STRING(FIX(spectra.var.stepdelta_tot)), 1) + '_' + 'bkg' + STRTRIM(STRING(FIX(spectra.var.fixbkg_TiBi[0])), 1) + '-' + $
+              STRTRIM(STRING(FIX(spectra.var.fixbkg_TiBi[1])), 1)
+    nomefilesav = STRMID(spectra.var.nomefile, 0, STRLEN(var.nomefile)-4)+ '_Al_' + addinfo + '_' + coinc_info.coinc + '.sav'
+
+    SAVE, var_Al_tot, filename = nomefilesav
+    
+    RETURN, var_Al_tot
+   
+  ENDIF
+
+ENDIF
+
+information = DIALOG_MESSAGE('Successful coincidence analysis. Remember to include the C/N normalization. Run WTITAN to display the results.', /information)
+
+END
\ No newline at end of file

TITAN/TITAN_COINCCYCLE_AL.pro

+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;; TITAN_COINCCYCLE_AL ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+
+; this function performs cycle over specified NaI channels 
+; and compute results over the Al Ge channel range
+
+FUNCTION TITAN_COINCCYCLE_AL, spectra,   $ ; spectra structure, passed from Titan main, that contains everything is needed   - INPUT
+                              coinc_info   ; information abaoutn coincidence - structure                                     - INPUT
+
+var_Al = TITAN_CREATE_VAR_AL(coinc_info)
+
+FOR m = 0, coinc_info.dim_m-1 DO BEGIN
+
+  FOR w = 0, coinc_info.dim_w-1 DO BEGIN
+
+    IF ((coinc_info.mean[m]-coinc_info.width[w]) GE coinc_info.limitinf) AND ((coinc_info.mean[m]+coinc_info.width[w]) LE coinc_info.limitsup) THEN BEGIN
+
+      down = coinc_info.mean[m] - coinc_info.width[w]
+      up   = coinc_info.mean[m] + coinc_info.width[w]
+      hge  = FLTARR(spectra.var.Ge_domain[1]-spectra.var.Ge_domain[0])
+      hge2 = TOTAL(spectra.h2d[coinc_info.mini_Al:coinc_info.maxi_Al,down+0:up+0], 2)
+      hge[coinc_info.mini_Al:coinc_info.maxi_Al] = hge2
+      ff   = hge(coinc_info.mini_Al:coinc_info.maxi_Al)
+      res  = TITAN_MPPEAKFIT(spectra.xGe[coinc_info.mini_Al:coinc_info.maxi_Al], ff, [coinc_info.ch_Al26], fixbkg=coinc_info.fixbkg_Al, fixpar=coinc_info.fixpar_Al, fitmodel=coinc_info.fitmodel_Al)
+
+      IF res.status EQ 0 THEN BEGIN
+        
+        var_Al.h_Al26[m,w]       = res.par[0]
+        var_Al.sh_Al26[m,w]      = res.spar[0]
+        var_Al.ch_Al26[m,w]      = res.par[1]
+        var_Al.sch_Al26[m,w]     = res.spar[1]
+        var_Al.sigma_Al26[m,w]   = res.par[2]
+        var_Al.ssigma_Al26[m,w]  = res.spar[2]
+        var_Al.tail_Al26[m,w]    = res.par[4]
+        var_Al.stail_Al26[m,w]   = res.spar[4]
+        var_Al.beta_Al26[m,w]    = res.par[5]
+        var_Al.sbeta_Al26[m,w]   = res.spar[5]
+        var_Al.bkg0_Al[m,w]      = res.par[6]
+        var_Al.sbkg0_Al[m,w]     = res.spar[6]
+        var_Al.bkg1_Al[m,w]      = res.par[7]
+        var_Al.sbkg1_Al[m,w]     = res.spar[7]
+        var_Al.bkg2_Al[m,w]      = res.par[8]
+        var_Al.sbkg2_Al[m,w]     = res.spar[8]
+        var_Al.step_Al[m,w]      = res.par[3]
+        var_Al.sstep_Al[m,w]     = res.spar[3]
+        var_Al.E_Al26[m,w]       = res.E[0]
+        var_Al.sE_Al26[m,w]      = res.sE[0]
+        var_Al.cn_Al26[m,w]      = res.E[0]/coinc_info.E_Al26_N
+        var_Al.scn_Al26[m,w]     = res.sE[0]/coinc_info.E_Al26_N
+        var_Al.status_Al[m,w]    = res.status
+        var_Al.chisq_Al[m,w]     = res.chisq
+        var_Al.spectra_Al[*,m,w] = res.hge
+        var_Al.yfit_Al[*,m,w]    = res.yfit
+
+      ENDIF ELSE BEGIN
+
+        var_Al.sE_Al26[m,w]      = 1
+        var_Al.status_Al[m,w]    = 1
+        var_Al.spectra_Al[*,m,w] = ff
+
+      ENDELSE
+
+    ENDIF
+
+  ENDFOR
+
+ENDFOR
+
+kk = WHERE(var_Al.sE_Al26 GT 0)
+var_Al.snr_Al26[kk] = var_Al.E_Al26/var_Al.sE_Al26[kk]
+
+RETURN, var_Al
+
+END
\ No newline at end of file

TITAN/TITAN_COINCCYCLE_TIBI.pro

+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;; TITAN_COINCCYCLE_TIBI ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+
+; this function performs cycle over specified NaI channels 
+; and compute results over the TiBi Ge channel range
+
+FUNCTION TITAN_COINCCYCLE_TIBI, spectra,   $ ; spectra structure, passed from Titan main, that contains everything is needed   - INPUT
+                                coinc_info   ; information abaoutn coincidence - structure                                     - INPUT
+ 
+var_tibi = TITAN_CREATE_VAR_TIBI(coinc_info)
+        
+FOR m = 0, coinc_info.dim_m-1 DO BEGIN
+
+  FOR w = 0, coinc_info.dim_w-1 DO BEGIN
+
+    IF ((coinc_info.mean[m]-coinc_info.width[w]) GE coinc_info.limitinf) AND ((coinc_info.mean[m]+coinc_info.width[w]) LE coinc_info.limitsup) THEN BEGIN
+
+      down = coinc_info.mean[m] - coinc_info.width[w]
+      up   = coinc_info.mean[m] + coinc_info.width[w]
+      hge  =FLTARR(spectra.var.Ge_domain[1]-spectra.var.Ge_domain[0])
+      hge2 = TOTAL(spectra.h2d(coinc_info.mini_TiBi:coinc_info.maxi_TiBi,down+0:up+0),2)
+      hge[coinc_info.mini_TiBi:coinc_info.maxi_TiBi] = hge2
+      ff   = hge[coinc_info.mini_TiBi:coinc_info.maxi_TiBi]
+      res  = TITAN_MPPEAKFIT(spectra.xGe[coinc_info.mini_TiBi:coinc_info.maxi_TiBi], ff, [coinc_info.ch_Bi214,coinc_info.ch_Ti44], fixbkg=coinc_info.fixbkg_TiBi, fixpar=coinc_info.fixpar_TiBi, fitmodel=coinc_info.fitmodel_TiBi)
+
+      IF res.status EQ 0 THEN BEGIN
+
+        var_tibi.h_Ti44[m,w] = res.par[6]
+        var_tibi.sh_Ti44[m,w] = res.spar[6]
+        var_tibi.h_Bi214[m,w] = res.par[0]
+        var_tibi.sh_Bi214[m,w] = res.spar[0]
+        var_tibi.ch_Ti44[m,w] = res.par[7]
+        var_tibi.sch_Ti44[m,w] = res.spar[7]
+        var_tibi.ch_Bi214[m,w] = res.par[1]
+        var_tibi.sch_Bi214[m,w] = res.spar[1]
+        var_tibi.sigma_Ti44[m,w] = res.par[8]
+        var_tibi.ssigma_Ti44[m,w] = res.spar[8]
+        var_tibi.sigma_Bi214[m,w] = res.par[2]
+        var_tibi.ssigma_Bi214[m,w] = res.spar[2]
+        var_tibi.tail_Ti44[m,w] = res.par[10]
+        var_tibi.stail_Ti44[m,w] = res.spar[10]
+        var_tibi.tail_Bi214[m,w] = res.par[4]
+        var_tibi.stail_Bi214[m,w] = res.spar[4]
+        var_tibi.beta_Ti44[m,w] = res.par[11]
+        var_tibi.sbeta_Ti44[m,w] = res.spar[11]
+        var_tibi.beta_Bi214[m,w] = res.par[5]
+        var_tibi.sbeta_Bi214[m,w] = res.spar[5]
+        var_tibi.bkg0_TiBi[m,w] = res.par[12]
+        var_tibi.sbkg0_TiBi[m,w] = res.spar[12]
+        var_tibi.bkg1_TiBi[m,w] = res.par[13]
+        var_tibi.sbkg1_TiBi[m,w] = res.spar[13]
+        var_tibi.bkg2_TiBi[m,w] = res.par[14]
+        var_tibi.sbkg2_TiBi[m,w] = res.spar[14]
+        var_tibi.step_TiBi[m,w] = res.par[15]
+        var_tibi.sstep_TiBi[m,w] = res.spar[15]
+        var_tibi.E_Ti44[m,w] = res.E[1]
+        var_tibi.sE_Ti44[m,w] = res.sE[1]
+        var_tibi.E_Bi214[m,w] = res.E[0]
+        var_tibi.sE_Bi214[m,w] = res.sE[0]
+        var_tibi.status_TiBi[m,w] = res.status
+        var_tibi.chisq_TiBi[m,w]  = res.chisq
+        var_tibi.spectra_TiBi[*,m,w] = res.hge
+        var_tibi.yfit_TiBi[*,m,w] = res.yfit
+
+
+      ENDIF ELSE BEGIN
+
+        var_tibi.sE_Bi214[m,w] = 1.
+        var_tibi.sE_Ti44[m,w] = 1.
+        var_tibi.status_TiBi[m,w] = 1.
+        var_tibi.spectra_TiBi[*,m,w] = ff
+
+      ENDELSE
+
+
+    ENDIF
+
+  ENDFOR
+
+ENDFOR
+
+kk = WHERE(var_tibi.sE_Ti44 GT 0)
+hh = WHERE(var_tibi.sE_Bi214 GT 0)
+var_tibi.snr_Ti44[kk]  = var_tibi.e_Ti44[kk]/var_tibi.sE_Ti44[kk]
+var_tibi.snr_Bi214[hh] = var_tibi.e_Bi214[hh]/var_tibi.sE_Bi214[hh]
+eTi    = var_tibi.E_Ti44/MAX(var_tibi.E_Ti44)
+seTi   = var_tibi.sE_Ti44/MAX(var_tibi.E_Ti44)
+eBi    = var_tibi.E_Bi214/MAX(var_tibi.E_Bi214)
+seBi   = var_tibi.sE_Bi214/MAX(var_tibi.E_Bi214)
+var_tibi.R_TiBi   = eTi/(eBi+1)
+var_tibi.sR_TiBi  = SQRT((1/(1+eBi)^2)*(seTi^2+(var_tibi.R_TiBi*seBi)^2))
+var_tibi.sR_TiBi  = SQRT((1/(1+eBi)^2)*(seTi^2+(var_tibi.R_TiBi*seBi)^2))
+var_tibi.Q_TiBi   = (eTi)*(1-eBi)
+var_tibi.sQ_TiBi  = SQRT(((1-eBi)*seTi)^2+(eTi*seBi)^2)
+var_tibi.R1_TiBi  = var_tibi.R_TiBi*var_tibi.snr_Ti44
+var_tibi.sR1_TiBi = var_tibi.sR_TiBi*var_tibi.snr_Ti44
+var_tibi.Q1_TiBi  = var_tibi.Q_TiBi*var_tibi.snr_Ti44
+var_tibi.sQ1_TiBi = var_tibi.sQ_TiBi*var_tibi.snr_Ti44
+ 
+RETURN, var_TiBi
+
+END
\ No newline at end of file

TITAN/TITAN_COINC_PAR.pro

+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;; TITAN_COINCPAR ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+
+; this function defines main parameters to be passed to the 
+; coincidence routines. Parameters are returned in a structure.
+
+FUNCTION TITAN_COINC_PAR, spectra ; spectra structure, passed from Titan main, that contains everything is needed
+
+isNaIps = STRLEN(FILE_SEARCH('*_NaI.ps'))
+if isNaIps gt 0 then isNaIps = 1
+spectratag  = tag_names(spectra)
+iispectra   = WHERE(spectratag EQ 'C511')
+if iispectra[0] EQ -1 then isNaI = 0 else isNaI = 1
+if isNaIps EQ 0 AND isNaI EQ 1 then isNaI = 0
+
+; call the Titan WindowsMet routine
+TITAN_WINDOWSMET, isNaI, isNaIps, spectra
+
+IF isNaIps NE 0 THEN BEGIN
+
+  coinc_par = {                   $
+              ISNAI:isNaI,        $
+              ISNAIPS:isNaIps,    $
+              NOMPEAK_511:0.,     $
+              NOMSIGMA_511:0.,    $
+              NOMPEAK_1022:0.,    $
+              NOMSIGMA_1022:0.,   $
+              E_AL26_N:0.,        $
+              SE_AL26_N:0.,       $
+              CH_BI214:0.,        $
+              CH_TI44:0.,         $
+              CH_AL26:0.,         $
+              SIGMA_BI214:0.,     $
+              SIGMA_TI44:0.,      $
+              SIGMA_AL26:0.,      $
+              MINI_TIBI:0.,       $
+              MAXI_TIBI:0.,       $
+              MINI_AL:0.,         $
+              MAXI_AL:0.          $
+              }
+  
+  coinc_par.nompeak_511   = spectra.c511.ch
+  coinc_par.nomsigma_511  = spectra.c511.sigma
+  coinc_par.nompeak_1022  = spectra.c1022.ch
+  coinc_par.nomsigma_1022 = spectra.c1022.sigma
+  
+  nomefiletxt = STRMID(spectra.var.nomefile, 0, STRLEN(spectra.var.nomefile)-4)+'_titan.txt'
+  READCOL,nomefiletxt, id, multi, status, t, peak_ch, A, sA, x0, sx0, sigma1, ssigma, step, sstep, tai, stai, bet, sbet, chisq, Energy, sigEnergy, CPM, sigCPM, Kev, sigKev, $
+    lineName, linekeV, format = '(I,I,I,I,I,F,F,F,F,F,F,F,F,F,F,F,F,F,F,F,F,F,F,F,A,A)', /silent
+  
+  iiAl = WHERE(lineName EQ 'Al26' AND linekeV EQ '(1808.65)')
+  coinc_par.E_AL26_N  = Energy[iiAl]
+  coinc_par.sE_Al26_N = sigEnergy[iiAl]
+  
+  FOR i = 0, spectra.kev2chcal.order DO BEGIN
+  
+    coinc_par.ch_Bi214 = coinc_par.ch_Bi214 + spectra.kev2chcal.calpar[i]*spectra.var.Bi214^i
+    coinc_par.ch_Ti44  = coinc_par.ch_Ti44 + spectra.kev2chcal.calpar[i]*spectra.var.Ti44^i
+    coinc_par.ch_Al26  = coinc_par.ch_Al26 + spectra.kev2chcal.calpar[i]*spectra.var.Al26^i
+  
+  ENDFOR
+  
+  FOR i = 0, spectra.ch2sigmacal.order DO BEGIN
+  
+    coinc_par.sigma_Bi214 = coinc_par.sigma_Bi214 + spectra.ch2sigmacal.calpar[i]*coinc_par.ch_Bi214^i
+    coinc_par.sigma_Ti44  = coinc_par.sigma_Ti44 + spectra.ch2sigmacal.calpar[i]*coinc_par.ch_Ti44^i
+    coinc_par.sigma_Al26  = coinc_par.sigma_Al26 + spectra.ch2sigmacal.calpar[i]*coinc_par.ch_Al26^i
+  
+  ENDFOR
+  
+  coinc_par.mini_TiBi = ROUND(coinc_par.ch_Bi214-spectra.var.fixbkg_TiBi[0])
+  coinc_par.maxi_TiBi = ROUND(coinc_par.ch_Ti44+spectra.var.fixbkg_TiBi[1])
+  coinc_par.mini_Al   = ROUND(coinc_par.ch_Al26-spectra.var.fixbkg_Al[0])
+  coinc_par.maxi_Al   = ROUND(coinc_par.ch_Al26+spectra.var.fixbkg_Al[1])
+  
+  fixpar_TiBi = [[1,coinc_par.ch_Bi214],[2,coinc_par.sigma_Bi214],[7,coinc_par.ch_Ti44],[8,coinc_par.sigma_Ti44],[13,0]]
+  fixpar_Al=[[1,coinc_par.ch_Al26],[2,coinc_par.sigma_Al26],[7,0.],[8,0.],[9,0.]]
+  
+  coinc_par = create_struct(coinc_par, 'FIXPAR_TIBI', fixpar_TiBi, 'FIXPAR_AL', fixpar_Al)
+  
+ENDIF ELSE coinc_par = {ISNAIPS:isNaIps, ISNAI:isNaI}
+ 
+RETURN, coinc_par
+
+END
\ No newline at end of file

TITAN/TITAN_CREATE_VAR_AL.pro

+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;; TITAN_CREATE_VAR_AL ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;
+; function that create variables structure for Al coincidence
+
+FUNCTION TITAN_CREATE_VAR_AL, coinc_info ; information about coincidence - structure     - INPUT
+
+
+a = FLTARR(coinc_info.dim_m, coinc_info.dim_w)
+b = FLTARR(coinc_info.maxi_Al-coinc_info.mini_Al+1, coinc_info.dim_m,coinc_info.dim_w)
+c = INTARR(coinc_info.maxi_Al-coinc_info.mini_Al+1, coinc_info.dim_m,coinc_info.dim_w)
+
+var_Al = {                      $
+          H_AL26:a,             $  
+          SH_AL26:a,            $
+          CH_AL26:a,            $
+          SCH_AL26:a,           $
+          SIGMA_AL26:a,         $
+          SSIGMA_AL26:a,        $
+          TAIL_AL26:a,          $
+          STAIL_AL26:a,         $
+          BETA_AL26:a,          $
+          SBETA_AL26:a,         $
+          BKG0_AL:a,            $
+          SBKG0_AL:a,           $
+          BKG1_AL:a,            $
+          SBKG1_AL:a,           $
+          BKG2_AL:a,            $
+          SBKG2_AL:a,           $
+          STEP_AL:a,            $
+          SSTEP_AL:a,           $
+          E_AL26:a,             $
+          SE_AL26:a,            $
+          CN_AL26:a,            $
+          SCN_AL26:a,           $
+          STATUS_AL:a,          $
+          CHISQ_AL :a,          $
+          SNR_AL26:a,           $
+          YFIT_AL:b,            $
+          SPECTRA_AL:c,         $
+          COINC_INFO:coinc_info $
+          }
+
+RETURN, var_al
+
+END
\ No newline at end of file

TITAN/TITAN_CREATE_VAR_TIBI.pro

+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;; TITAN_CREATE_VAR_TIBI ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+;
+; function that create variables structure for TiBi coincidence
+
+FUNCTION TITAN_CREATE_VAR_TIBI, coinc_info  ; information about coincidence - structure  - INPUT
+
+a = FLTARR(coinc_info.dim_m, coinc_info.dim_w)
+b = FLTARR(coinc_info.maxi_TiBi-coinc_info.mini_TiBi+1, coinc_info.dim_m,coinc_info.dim_w)
+c = INTARR(coinc_info.maxi_TiBi-coinc_info.mini_TiBi+1, coinc_info.dim_m,coinc_info.dim_w)
+
+var_tibi = {                      $
+            H_TI44:a,             $
+            SH_TI44:a,            $
+            H_BI214:a,            $
+            SH_BI214:a,           $
+            CH_TI44:a,            $
+            SCH_TI44:a,           $
+            CH_BI214:a,           $
+            SCH_BI214:a,          $
+            SIGMA_TI44:a,         $
+            SSIGMA_TI44:a,        $
+            SIGMA_BI214:a,        $
+            SSIGMA_BI214:a,       $
+            TAIL_TI44:a,          $
+            STAIL_TI44:a,         $
+            TAIL_BI214:a,         $
+            STAIL_BI214:a,        $
+            BETA_TI44:a,          $
+            SBETA_TI44:a,         $
+            BETA_BI214:a,         $
+            SBETA_BI214:a,        $
+            BKg0_TIBI:a,          $
+            SBKG0_TIBI:a,         $
+            BKG1_TIBI:a,          $
+            SBKG1_TIBI:a,         $
+            BKG2_TIBI:a,          $
+            SBKG2_TIBI:a,         $
+            STEP_TIBI:a,          $
+            SSTEP_TIBI:a,         $
+            E_TI44:a,             $
+            SE_TI44:a,            $
+            E_BI214:a,            $
+            SE_BI214:a,           $
+            STATUS_TIBI:a,        $
+            CHISQ_TIBI :a,        $
+            SNR_TI44:a,           $
+            SNR_BI214:a,          $
+            R_TIBI:a,             $
+            Q_TIBI:a,             $
+            R1_TIBI:a,            $
+            Q1_TIBI:a,            $
+            SR_TIBI:a,            $
+            SQ_TIBI:a,            $
+            SR1_TIBI:a,           $
+            SQ1_TIBI:a,           $
+            YFIT_TIBI:b,          $
+            SPECTRA_TIBI:c,       $
+            COINC_INFO:coinc_info $
+            }
+
+RETURN, var_tibi
+
+END
\ No newline at end of file