Add code and data dependencies

.gitmodules

@@ -0,0 +1,12 @@
+[submodule "src/structural"]
+	path = src/structural
+	url = https://github.com/psychoinformatics-de/studyforrest-data-structural.git
+[submodule "src/templatetransforms"]
+	path = src/templatetransforms
+	url = https://github.com/psychoinformatics-de/studyforrest-data-templatetransforms.git
+[submodule "src/aligned"]
+	path = src/aligned
+	url = https://github.com/psychoinformatics-de/studyforrest-data-aligned.git
+[submodule "src/freesurfer"]
+	path = src/freesurfer
+	url = https://github.com/psychoinformatics-de/studyforrest-data-freesurfer.git

README.rst

@@ -0,0 +1,69 @@
+studyforrest.org Dataset
+************************
+
+|license| |access| |doi|
+
+Some
+====
+
+
+For more information about the project visit: http://studyforrest.org
+
+
+How to obtain the data files
+----------------------------
+
+This repository contains metadata and information on the identity of all
+included files. However, the actual content of the (sometime large) data
+files is stored elsewhere. To obtain any dataset component, git-annex_ is
+required in addition to Git_.
+
+1. Clone this repository to the desired location.
+2. Enter the directory with the local clone and run::
+
+     git annex init
+
+   Older versions of git-annex may require you to run the following
+   command immediately afterwards::
+
+     git annex enableremote mddatasrc
+
+Now any desired dataset component can be obtained by using the ``git annex get``
+command. To obtain the entire dataset content run::
+
+     git annex get .
+
+Keep data up-to-date
+--------------------
+
+If updates to this dataset are made in the future, update any local clone by
+running::
+
+     git pull
+
+followed by::
+
+     git annex get .
+
+to fetch all new files.
+
+
+
+
+.. _Git: http://www.git-scm.com
+
+.. _git-annex: http://git-annex.branchable.com/
+
+.. |license|
+   image:: https://img.shields.io/badge/license-PDDL-blue.svg
+    :target: http://opendatacommons.org/licenses/pddl/summary
+    :alt: PDDL-licensed
+
+.. |access|
+   image:: https://img.shields.io/badge/data_access-unrestricted-green.svg
+    :alt: No registration or authentication required
+
+.. |doi|
+   image:: https://img.shields.io/badge/doi-missing-lightgrey.svg
+    :target: http://dx.doi.org/
+    :alt: DOI

code/RetMap_phaseshift

@@ -0,0 +1,71 @@
+#!/usr/bin/env python
+
+# Shift the phase of the input nifti.
+# It is useful to correct for different wedge or ring starting points,
+# ie. if the polar angle wedge doesn't start from the upper vertical
+# meridian and your further processing expects it to start at the upper
+# vertical meridian.
+
+def make_phaseshift(niftipath, shift):
+    import nibabel as nib
+    import numpy as np
+    import os
+
+
+    cond = os.path.basename(niftipath).split("_")[0]
+    print "condition: ", cond
+    
+    ## get condition
+    #dsback2niftisave = os.path.dirname(os.path.dirname(niftipath))+'/post_processing/'+cond+"_phShift.nii.gz"
+    dsback2niftisave = 'post_processing/'+cond+"_phShift.nii.gz"
+
+    ## set shift depending on condition
+    #~ if cond == 'ccw':
+        #~ shift = -281.25
+    #~ elif cond == 'clw':
+        #~ shift = 101.25
+    #~ elif cond == "con":
+        #~ shift = -337.5
+    #~ elif cond == "exp":
+        #~ shift = 0
+    shift=float(shift)
+    print "shifting by degree: ", shift
+    ## get nifti 
+    nifti= nib.load(niftipath)
+
+    data = nifti.get_data()
+
+    # do phase shift depending on condition
+    tmp0 = (data[:,:,:,0] + shift)%360
+    #tmp0[tmp0 >= 360] = tmp0[tmp0 >= 360] - 360
+    #tmp0[tmp0 < 0] = tmp0[tmp0 < 0] + 360
+    
+    #~ if cond == 'ccw' or cond == "con":
+        #~ tmp0 = 360-tmp0           
+    
+    data[:,:,:,0]=tmp0
+
+
+    dsback2nifti = nib.Nifti1Image(data, nifti.get_affine())
+    nib.save(dsback2nifti, dsback2niftisave)
+    
+
+if __name__ == "__main__":
+    import sys
+    niftipath = sys.argv[1]
+    shift = sys.argv[2]
+    print "processing file: \n", niftipath
+    make_phaseshift(niftipath, shift)
+
+
+
+
+
+
+
+
+
+
+
+
+

code/combine_volumes

@@ -0,0 +1,71 @@
+#!/usr/bin/env python
+
+# Combines the phase volumes given to generate teh pahse field maps.
+# The code is translated from AFNIs implementation:
+#
+#   for (i=0; i<nvox;++i) {
+#      dif[i] = (phi1[i]-phi2[i])/2.0;
+#      if (ABS(dif[i]) < 90/n || !rpud->fixsum) { /* should this be 90 / n ? */ 
+#         sum[i] = (phi1[i]+phi2[i])/2.0;
+#      } else { /* Too big a difference in phase, likely 
+#               summing about zero degrees where wrapping from
+#               360 occurs*/
+#         sum[i] = (phi1[i]+phi2[i])/2.0-180.0/n; 
+#         if (sum[i]<0) sum[i] = 360/n + sum[i];
+#      }
+#
+
+def combine_volumes(volume1_path, volume2_path, output_path):
+    import nibabel as nib
+    import numpy as np
+    import os
+    import scipy.stats as ss
+
+    volume1 = nib.load(volume1_path)
+    volume2 = nib.load(volume2_path)
+    real_volume1_data=volume1.get_data()[:,:,:,0]
+    real_volume2_data=volume2.get_data()[:,:,:,0]
+    combined_volume=ss.circmean([real_volume1_data, real_volume2_data],high=360,low=0,axis=0)
+
+#    combined_volume=np.zeros(real_volume1_data.shape)
+
+#    it = np.nditer(real_volume1_data, flags=['multi_index'])
+#    while not it.finished:
+#        diff=0.5*(real_volume1_data[it.multi_index]-real_volume2_data[it.multi_index])
+#        if abs(diff)< 90:
+#            combined_volume[it.multi_index]=0.5*(real_volume1_data[it.multi_index]+real_volume2_data[it.multi_index])
+#        else:
+#            sum_value=0.5*(real_volume1_data[it.multi_index]+real_volume2_data[it.multi_index])-180
+#            if sum_value<0:
+#                sum_value+=360
+#            combined_volume[it.multi_index]=sum_value
+#        it.iternext()
+     
+    volume_returned_data=(volume1.get_data()+volume2.get_data())*0.5
+    volume_returned_data[:,:,:,0]=combined_volume
+    nib.Nifti1Image(volume_returned_data, volume1.get_affine()).to_filename(output_path)
+    return volume_returned_data
+
+    
+
+if __name__ == "__main__":
+    import sys
+    volume1_path = sys.argv[1]
+    volume2_path = sys.argv[2]
+    output_path = sys.argv[3]
+    combined_volume=combine_volumes(volume1_path, volume2_path, output_path)
+    print ">> Combining volumes completed"
+
+
+
+
+
+
+
+
+
+
+
+
+
+

code/easyret_gui

@@ -0,0 +1,170 @@
+#!/usr/bin/env python
+
+import os
+import wx
+import wx.lib.filebrowsebutton as libtest
+import wx.lib.agw.floatspin as FS
+from os.path import join as opj
+
+MAIN_WINDOW_DEFAULT_SIZE = (650,700)
+
+class Frame(wx.Frame):
+    
+    def __init__(self, parent, id, title):
+        style=wx.DEFAULT_FRAME_STYLE ^ (wx.RESIZE_BORDER) # XOR to remove the resizeable border        
+        wx.Frame.__init__(self, parent, id, title=title, size=MAIN_WINDOW_DEFAULT_SIZE, style=style)
+        self.Center() # open in the centre of the screen
+        self.panel = wx.Panel(self)
+        self.subject_id_label = wx.StaticText(self.panel, pos=(50, 45), size=(140,-1), label="Enter Subject ID:")
+        self.subject_id = wx.TextCtrl(self.panel, pos=(200, 40), size=(60,-1))
+        self.contract = libtest.FileBrowseButton(self.panel, pos=(50, 80), initialValue='Enter path to Contracting volume', toolTip='Enter/Browse path to Contracting volume', fileMask='*.nii.gz', size=(500,30), labelWidth=100)
+        self.expand = libtest.FileBrowseButton(self.panel, pos=(50, 120), initialValue='Enter path to Expanding volume', toolTip='Enter/Browse path to Expanding volume', fileMask='*.nii.gz', size=(500,30), labelWidth=100)
+        self.clock = libtest.FileBrowseButton(self.panel, pos=(50, 160), initialValue='Enter path to Clockwise volume', toolTip='Enter/Browse path to Clockwise volume', fileMask='*.nii.gz', size=(500,30), labelWidth=100)
+        self.counter = libtest.FileBrowseButton(self.panel, labelWidth=10, pos=(50, 200),  initialValue='Enter path to Counter clockwise volume', toolTip='Enter/Browse path to Counter clockwise volume', fileMask='*.nii.gz', size=(500,30))
+        self.mask = libtest.FileBrowseButton(self.panel, labelWidth=10, pos=(50, 240),  initialValue='Enter path to mask', toolTip='Enter/Browse path to brain mask of the retmapping volumes', fileMask='*.nii.gz', size=(350,30))
+        self.FWHM_label = wx.StaticText(self.panel, pos=(50, 285), size=(250,-1), label="Spatial Smoothing FWHM (mm):")
+        self.FWHM = FS.FloatSpin(self.panel, value='4', pos=(300, 280), size=(45, -1), min_val=0, max_val=10, increment=1)
+        self.FWHM.SetDigits(0) 
+        self.output_folder_label = wx.StaticText(self.panel, pos=(50, 325), size=(250,-1), label="Select output folder:")
+        self.output = wx.DirPickerCtrl(self.panel, pos=(250, 320))
+        self.postprocessingshift = wx.CheckBox(self.panel, -1, 'Post Processing Shift', pos=(50, 385))
+        self.Bind(wx.EVT_CHECKBOX, self.OnPostProcessShift, self.postprocessingshift)
+        self.surfacemaps = wx.CheckBox(self.panel, -1, 'Surface Maps', pos=(50, 415))
+        self.Bind(wx.EVT_CHECKBOX, self.OnSurfaceMaps, self.surfacemaps)
+        self.debugmode = wx.CheckBox(self.panel, -1, 'Debug Mode', pos=(50, 580))      
+        self.process_button = wx.Button(self.panel, pos=(470, 580), size=(80, 25), label="EasyRet")
+        self.Bind(wx.EVT_BUTTON, self.OnProcess, self.process_button)
+        self.CreateMenuBar()
+        self.CreateStatusBar()
+    
+    
+    
+        
+    def CreateMenuBar(self):
+
+        menuBar = wx.MenuBar()
+        # Tell our Frame about this MenuBar
+        self.SetMenuBar(menuBar)
+        menuFile = wx.Menu()
+
+        # NOTE on wx ids - they're used everywhere, we don't care about them
+        # Used to handle events and other things
+        # An id can be -1 or wx.ID_ANY, wx.NewId(), your own id
+        # Get the id using object.GetId()
+        
+        # add a Help menu with an About item
+        menuHelp = wx.Menu()
+        menuBar.Append(menuHelp, '&Help')
+        helpMenuItem = menuHelp.Append(-1, '&About', 'About EasyRet')
+        self.Bind(wx.EVT_MENU, self.OnAbout, helpMenuItem)        
+        
+        # create a menu item for Exit and bind it to the OnExit function       
+        exitMenuItem = menuHelp.Append(-1, '&Exit', 'Exit EasyRet')        
+        self.Bind(wx.EVT_MENU, self.OnExit, exitMenuItem)
+        
+    def OnAbout(self,e):
+        # A message dialog box with an OK button. wx.OK is a standard ID in wxWidgets.
+        
+        help_text='I am working on it\n please bear with me'
+        dlg = wx.MessageDialog( self, help_text, "About EasyRet", wx.OK)
+        dlg.ShowModal() # Show it
+        dlg.Destroy() # finally destroy it when finished.
+        
+    def OnPostProcessShift(self, event):   
+        if self.postprocessingshift.GetValue(): 
+            self.con_shift_text = wx.StaticText(self.panel, pos=(250, 360), size=(200,-1), label="Con")        
+            self.con_shift = FS.FloatSpin(self.panel, value='22.5', pos=(250, 380), size=(70, -1), min_val=-360, max_val=360, increment=0.25)
+            self.con_shift.SetDigits(2)
+            self.exp_shift_text = wx.StaticText(self.panel, pos=(320, 360), size=(200,-1), label="Exp")        
+            self.exp_shift = FS.FloatSpin(self.panel, value='22.5', pos=(320, 380), size=(70, -1), min_val=-360, max_val=360, increment=0.25)
+            self.exp_shift.SetDigits(2)
+            self.clw_shift_text = wx.StaticText(self.panel, pos=(390, 360), size=(200,-1), label="Clw")        
+            self.clw_shift = FS.FloatSpin(self.panel, value='22.5', pos=(390, 380), size=(70, -1), min_val=-360, max_val=360, increment=0.25)
+            self.clw_shift.SetDigits(2)
+            self.ccw_shift_text = wx.StaticText(self.panel, pos=(460, 360), size=(200,-1), label="Ccw")        
+            self.ccw_shift = FS.FloatSpin(self.panel, value='22.5', pos=(460, 380), size=(70, -1), min_val=-360, max_val=360, increment=0.25)
+            self.ccw_shift.SetDigits(2)
+         
+        else:
+            self.con_shift.Destroy()
+            self.con_shift_text.Destroy()
+            self.exp_shift.Destroy()
+            self.exp_shift_text.Destroy()
+            self.clw_shift.Destroy()
+            self.clw_shift_text.Destroy()
+            self.ccw_shift.Destroy()
+            self.ccw_shift_text.Destroy()
+            
+    def OnSurfaceMaps(self, event):   
+        if self.surfacemaps.GetValue(): 
+            self.anatomy = libtest.FileBrowseButton(self.panel, labelWidth=10, pos=(150, 450),  initialValue='Enter path to anatomy', toolTip='Enter/Browse path to the anatomical volume that was used for Freesurfer recon-all analysis', fileMask='*.nii.gz', size=(400,30))
+            self.transformation = libtest.FileBrowseButton(self.panel, labelWidth=10, pos=(150, 490),  initialValue='Enter transformation matrix', toolTip='Enter/Browse path to the fsl type transformation file for co-registering BOLD retmapping volumes to anatomical volumes', fileMask='*.mat', size=(400,30))
+            self.ROI = libtest.FileBrowseButton(self.panel, labelWidth=10, pos=(150, 530),  initialValue='Enter path to ROI mask', fileMask='*.nii.gz',  toolTip='Enter/Browse path to the ROI mask if you want your surface maps to be restricted to Occipital lobe only', size=(400,30))
+
+        else:
+            self.anatomy.Destroy()
+            self.transformation.Destroy()
+            self.ROI.Destroy()
+            
+    def OnProcess(self, event):
+        
+        from subprocess import PIPE, CalledProcessError, check_call, Popen
+        bash_processing_script=opj(os.path.dirname(os.path.realpath(__file__)), 'process_retmap')
+        subject_id = str(self.subject_id.GetValue())
+        contract_volume = str(self.contract.GetValue())
+        expand_volume = str(self.expand.GetValue())
+        clock_volume = str(self.clock.GetValue())
+        counter_volume = str(self.counter.GetValue())
+        mask = str(self.mask.GetValue())
+        FWHM = str(self.FWHM.GetValue())
+        output_folder = str(self.output.GetPath())
+        contract_shift=expand_shift=clock_shift=counter_shift=anatomy_volume=transformation_matrix=ROI=''
+        post_processing_shift=str(self.postprocessingshift.GetValue())
+        if self.postprocessingshift.GetValue():
+            contract_shift=str(self.con_shift.GetValue())
+            expand_shift=str(self.exp_shift.GetValue())
+            clock_shift=str(self.clw_shift.GetValue())
+            counter_shift=str(self.ccw_shift.GetValue())
+        
+        surface_maps=str(self.surfacemaps.GetValue())  
+        if self.surfacemaps.GetValue():
+            anatomy_volume=str(self.anatomy.GetValue())
+            transformation_matrix=str(self.transformation.GetValue())
+            ROI=str(self.ROI.GetValue())
+        debug_mode=str(self.debugmode.GetValue())
+        error_log=opj(output_folder, subject_id+'_error.log')
+        #scripts/pyretmap/code/processing_pipeline/process_retmap sub-16 /home/data/psyinf/forrest_gump/anondata/sub016/BOLD/task005_run001/bold_bold3Tp2_to_subjbold3Tp2.nii.gz /home/data/psyinf/forrest_gump/anondata/sub016/BOLD/task005_run003/bold_bold3Tp2_to_subjbold3Tp2.nii.gz /home/data/psyinf/forrest_gump/anondata/sub016/BOLD/task005_run004/bold_bold3Tp2_to_subjbold3Tp2.nii.gz /home/data/psyinf/forrest_gump/anondata/sub016/BOLD/task005_run002/bold_bold3Tp2_to_subjbold3Tp2.nii.gz /home/data/psyinf/forrest_gump/anondata/sub016/templates/bold3Tp2/brain_mask.nii.gz /home/data/exppsy/spark/Study_Forrest checked -337.5 0.0 101.25 -281.25 checked /home/data/psyinf/forrest_gump/anondata/sub016/templates/t1w/head.nii.gz /home/data/psyinf/forrest_gump/anondata/sub016/templates/bold3Tp2/in_t1w/xfm_6dof.mat ''
+
+        arguments_list=[bash_processing_script, subject_id, contract_volume, expand_volume, clock_volume, counter_volume, mask, FWHM, output_folder, post_processing_shift, contract_shift, expand_shift, clock_shift, counter_shift, surface_maps, anatomy_volume, transformation_matrix, ROI, debug_mode]
+        print '\n\n'
+        print ' '.join(arguments_list)
+        print '\n\n'
+        print 'Processing Retmap'
+
+        with open(error_log, "w") as f:
+            try:
+                check_call(arguments_list, stderr=f)
+            except CalledProcessError as e:
+                print(e)
+                exit(1)
+        print 'Done'
+        
+    def OnExit(self, event):
+        "Close the application by Destroying the object"
+        # TRY add in your own print call here
+        self.Destroy() # SHOW HELP SIDEBAR
+        
+    
+class App(wx.App):
+    
+    def OnInit(self):
+        self.frame = Frame(parent=None, id=-1, title='EasyRet')
+        self.frame.Show()
+        self.SetTopWindow(self.frame)
+        return True
+    
+if __name__ == "__main__":       
+    # make an App object, set stdout to the console so we can see errors
+    app = App(redirect=False)
+        
+    app.MainLoop()

code/process_retmap

@@ -0,0 +1,312 @@
+#!/bin/bash
+set -e
+set -u
+
+# TODO: 
+# - add copyright stuff (Ayan Doesn't know how to do it, Michael will have address it) 
+# - make comments nice [http://mywiki.wooledge.org/BashGuide/CommandsAndArguments]
+# - check if rings start at zero or 360 deg -> adapt comments for phase shift 
+# "<post_processing> variable makes the decision of doing post processing or not.... it is required also by the GUI...cannot be ommitted"
+# - $subject: yes this is FS subject ID
+# - make it beautifullll -> What?
+# - adopt "Example" to generalize or remove it
+
+### Example Run & Usage
+# process_retmap - retinotopic mapping processing. 
+# Takes corrected niftis and generates the  phasemaps for polar angle and 
+# eccentricity, as nifti and freesurfer overlay for surface presentation.
+# 
+# Usage #
+# process_retmap <subject ID> <con> <exp> <clw> <ccw> <template_mask> <output_folder> <post_processing> <contract_shift> <expand_shift> <clock_shift> <counter_shift> <surface_maps> <anatomy_volume> <retmap2anat_mat> <ROI_mask>
+# 
+# Input definition #
+# subject        : Freesurfer subject ID
+# con            : nifti for your contracting ring
+# exp            : nifti for your expanding ring
+# clw            : nifti for your clockwise wedge
+# ccw            : nifti for your counter clockwise wedge 
+# template_mask  : mask to mask you phasemap, eg. brain mask, occ. lobe mask. 
+#                  Needs to be coregistered to your input tSeries!
+# output_folder  : where generate folders and save output
+#
+# post_processing : "True": the shift values if you didn't start as AFNI
+#                   intended and you want to shift them manually. 
+#                   Therefore give your values
+#
+# contract_shift  : range of phase shift values from 0-360
+# expand_shift    : range of phase shift values from 0-360
+# clock_shift     : range of phase shift values from 0-360
+# counter_shift   : range of phase shift values from 0-360
+#
+# surface_maps    : "True" if you want overlays.
+#
+# anatomy_volume  : your reference anatomy for the overlay. It should 
+#                   be your freesurfer anatomy, because your freesurfer 
+#                   surfaces are generated from it awhile recon-all
+# retmap2anat_mat : a '*.mat' file used for the overlays
+# ROI_mask        : if you want to get the overlay for your ROI, too
+# smooth          : for fslmaths: create a gauss kernel of sigma mm and perform mean filtering
+#
+# Example without post processing and phase shift#
+#scripts/pyretmap/code/processing_pipeline/process_retmap \
+#sub-16 \
+#/home/data/psyinf/forrest_gump/anondata/sub016/BOLD/task005_run001/bold_bold3Tp2_to_subjbold3Tp2.nii.gz \
+#/home/data/psyinf/forrest_gump/anondata/sub016/BOLD/task005_run003/bold_bold3Tp2_to_subjbold3Tp2.nii.gz \
+#/home/data/psyinf/forrest_gump/anondata/sub016/BOLD/task005_run004/bold_bold3Tp2_to_subjbold3Tp2.nii.gz \
+#/home/data/psyinf/forrest_gump/anondata/sub016/BOLD/task005_run002/bold_bold3Tp2_to_subjbold3Tp2.nii.gz \
+#/home/data/psyinf/forrest_gump/anondata/sub016/templates/bold3Tp2/brain_mask.nii.gz \
+#/home/data/exppsy/spark/Study_Forrest \
+#False \
+#'' \
+#'' \
+#'' \
+#'' \
+#False \
+#'' \
+#'' \
+#'' \
+#2
+ 
+#~/Documents/Auswertung/paper-forrest_phase2_data/code/retmap/processing_pipeline/process_retmap \
+#sub-16 \
+#/home/data/psyinf/forrest_gump/anondata/sub016/BOLD/task005_run001/bold_bold3Tp2_to_subjbold3Tp2.nii.gz \
+#/home/data/psyinf/forrest_gump/anondata/sub016/BOLD/task005_run003/bold_bold3Tp2_to_subjbold3Tp2.nii.gz \
+#/home/data/psyinf/forrest_gump/anondata/sub016/BOLD/task005_run004/bold_bold3Tp2_to_subjbold3Tp2.nii.gz \
+#/home/data/psyinf/forrest_gump/anondata/sub016/BOLD/task005_run002/bold_bold3Tp2_to_subjbold3Tp2.nii.gz \
+#/home/data/psyinf/forrest_gump/anondata/sub016/templates/bold3Tp2/brain_mask.nii.gz \
+#/home/fkaule/temp \
+#False \
+#'' \
+#'' \
+#'' \
+#'' \
+#False \
+#'' \
+#'' \
+#'' \
+#2
+
+
+### get necessary files and folders
+
+subject=${1}
+con=${2}
+exp=${3}
+clw=${4}
+ccw=${5}
+template_mask=${6}
+FWHM=${7}
+output_folder=${8}
+post_processing=${9}
+contract_shift=${10}
+expand_shift=${11}
+clock_shift=${12}
+counter_shift=${13}
+surface_maps=${14}  
+anatomy_volume=${15} 
+retmap2anat_mat=${16}  
+ROI_mask=${17} 
+debug_mode_flag=${18}
+
+
+stimulus_cycle_time=32 #Tstim
+highpass=`bc <<< "scale=2; 2/($stimulus_cycle_time*3)"` # should be around 1/Tstim * 2/3   
+lowpass=`bc <<< "scale=2; $highpass*3"`   # app. highpass*3
+req_TR=`fslinfo $con | grep pixdim4 | tr -s ' '| cut -d ' ' -f 2`
+code_path=$(dirname $(readlink -f "$0"))
+
+sigma="$(echo "scale=2; $FWHM / 2.35" | bc)"
+echo $sigma
+echo $code_path
+ 
+retmapping_folder=$output_folder'/'$subject
+echo 'create retmap folder'
+rm -rf $retmapping_folder
+mkdir -p $retmapping_folder
+cd $retmapping_folder
+
+echo 'create temporary processing folder'
+mkdir -p 'pre_processing'
+
+
+declare -a order_of_processing=('contract' 'expand' 'clock' 'counter');
+declare -A shift=( [con]=$contract_shift [exp]=$expand_shift [clw]=$clock_shift [ccw]=$counter_shift );
+
+retmap_volumes=`fslnvols $con`
+waver -TR $req_TR -GAM -numout ${retmap_volumes} -inline 2@0.0 1@1.0 15@0.0 1@1.0 15@0.0 1@1.0 15@0.0 1@1.0 15@0.0 1@1.0 25@0.0 > 'refts.1D'
+
+### processing
+
+index=0
+for f in $con $exp $clw $ccw; do
+
+    retmap_run_folder='pre_processing/'${order_of_processing[$index]}
+    mkdir -p $retmap_run_folder
+
+    if [ "$FWHM" -gt 0 ]; then
+        echo "Spatial Smoothing ->" 
+        fslmaths "$f" -s $sigma $retmap_run_folder'/bold_sm.nii.gz'
+    else
+        cp "$f" $retmap_run_folder'/bold_sm.nii.gz'
+    fi
+    
+    echo "Slice Time Correction ->" 
+    stcIn=$retmap_run_folder'/bold_sm.nii.gz'
+    stcOut=$retmap_run_folder'/bold_sm_STC.nii.gz'
+    3dTshift -TR $req_TR -prefix $stcOut $stcIn
+    
+    echo "De-obliquing ->" 
+    warpIn=$retmap_run_folder'/bold_sm_STC.nii.gz'
+    warpOut=$retmap_run_folder'/bold_sm_STC_deobl.nii.gz'
+    3dWarp -prefix $warpOut -deoblique $warpIn
+    
+    echo "Re-orienting to standard orientation after De-obliquing -> "
+    reorIn=$retmap_run_folder'/bold_sm_STC_deobl.nii.gz'
+    reorOut=$retmap_run_folder'/bold_sm_STC_deobl_reor.nii.gz'
+    fslreorient2std $reorIn $reorOut
+
+    echo "Band Pass Filtering -> " 
+    bpIn=$retmap_run_folder'/bold_sm_STC_deobl_reor.nii.gz'
+    bpOut=$retmap_run_folder'/bold_sm_STC_deobl_reor_filt.nii.gz'
+    3dBandpass -prefix $bpOut -norm $highpass $lowpass $bpIn
+    
+    echo "masking ->"
+    maskingIn=$retmap_run_folder'/bold_sm_STC_deobl_reor_filt.nii.gz'
+    maskingOut=$retmap_run_folder'/bold_sm_STC_deobl_reor_filt_masked.nii.gz'
+    if [ -f "$template_mask" ]; then
+		fslmaths $maskingIn -mas $template_mask $maskingOut
+    else
+		cp $maskingIn $maskingOut
+    fi
+    index="$((index + 1))"    
+done
+
+
+
+
+mkdir -p 'raw_outputs'
+# 3dRetinoPhase: 
+# the core processing. If you want to get the "rawest" version of AFNI retmap,
+# use this function and ignore the rest of this script.
+# Here you find help: [http://afni.nimh.nih.gov/pub/dist/doc/program_help/3dRetinoPhase.html]
+3dRetinoPhase -con 'pre_processing/'${order_of_processing[0]}'/bold_sm_STC_deobl_reor_filt_masked.nii.gz' \
+              -ccw 'pre_processing/'${order_of_processing[3]}'/bold_sm_STC_deobl_reor_filt_masked.nii.gz' \
+              -exp 'pre_processing/'${order_of_processing[1]}'/bold_sm_STC_deobl_reor_filt_masked.nii.gz' \
+              -clw 'pre_processing/'${order_of_processing[2]}'/bold_sm_STC_deobl_reor_filt_masked.nii.gz' \
+              -prefix 'BRIK_files/retmap' \
+              -Tstim 32 \
+              -nrings 1 \
+              -nwedges 1 \
+              -pre_stim 0 \
+              -ref_ts 'refts.1D' \
+              -phase_estimate DELAY
+
+# AFNI2nifti: 
+# since AFNI has its own file format we want to leave
+# the AFNI-universe here, transform to nifti
+mri_convert $(find 'BRIK_files' -name 'retmap.ecc.field*.BRIK' -type f) 'raw_outputs/combined_eccentricity_map_field.nii.gz' --in_type afni
+mri_convert $(find 'BRIK_files' -name 'retmap.pol.field*.BRIK' -type f) 'raw_outputs/combined_polar_map_field.nii.gz' --in_type afni
+mri_convert $(find 'BRIK_files' -name 'retmap.ecc-+*.BRIK' -type f) 'raw_outputs/con_eccentricity_map.nii.gz' --in_type afni
+mri_convert $(find 'BRIK_files' -name 'retmap.pol++*.BRIK' -type f) 'raw_outputs/clw_polar_map.nii.gz' --in_type afni
+mri_convert $(find 'BRIK_files' -name 'retmap.ecc++*.BRIK' -type f) 'raw_outputs/exp_eccentricity_map.nii.gz' --in_type afni
+mri_convert $(find 'BRIK_files' -name 'retmap.pol-+*.BRIK' -type f) 'raw_outputs/ccw_polar_map.nii.gz' --in_type afni
+
+# phaseshift and merge to efd and pfd: 
+if [ "$post_processing" == "True" ]; then
+    echo "Performing post processing"
+    mkdir -p 'post_processing'
+    
+    # phaseshift
+    # if there is a shift in stimulus onset, ie.: 
+    # - the polar wedge doesn't start at the upper vertical meridian. 
+    #   Both (ccw & clw) need to start at the upper vertical meridian
+    # - the exp & con rings don't start zero ..
+    # You can phase shift them to do so and shift the output phasemap 
+    # from 3dRetinoPhase to this position.
+    # Normally this shouldn't be necessary, unlike you do sth. special or
+    # mess with your stimulus.
+    for fileToShift in $(find 'raw_outputs' -name '*_map.nii.gz'); do
+        echo 'Post processing phase shift '$fileToShift   
+        IFS='_' read -a conds <<< $(basename $fileToShift)
+        echo ${shift[${conds[0]}]}
+        $code_path'/RetMap_phaseshift' $fileToShift ${shift[${conds[0]}]}
+    done
+
+
+
+    # merge opposite phasemaps into fieldmaps for ECC and POL:
+    # As we expect a phaseshift, we have to merge the opposite phasemaps 
+    # here. 
+    # The outcome is as the 'ecc.field*.BRIK' and 'pol.field*.BRIK' 
+    # output from 3dRetinoPhase. So if your stimulus is as AFNI intended,
+    # you can stick with the 3dRetinoPhase output.
+    $code_path'/combine_volumes' 'post_processing/con_phShift.nii.gz' 'post_processing/exp_phShift.nii.gz' 'post_processing/combined_pipeline_eccentricity.nii.gz'
+    $code_path'/combine_volumes' 'post_processing/clw_phShift.nii.gz' 'post_processing/ccw_phShift.nii.gz' 'post_processing/combined_pipeline_polar.nii.gz'
+
+    #~ #for fileToRT in $(find 'post_processing' -name '*pipeline*.nii.gz'); do
+    #~ #    echo 'Generating rtview'$fileToRT   
+    #~ #    $code_path'/deg2rtview' $fileToRT
+    #~ #done
+
+
+
+fi
+
+
+# get the phasemaps overlays for freesurfer surfaces
+if [ "$surface_maps" == "True" ]; then
+    mkdir -p 'alignments'       
+    mkdir -p 'surface_maps'
+
+    if [ "$post_processing" == "True" ]; then
+    
+        if [ "$debug_mode_flag" == "True" ]; then
+
+            files_to_process=`find 'post_processing' -name '*.nii.gz*' -type f`
+        else
+            files_to_process=`find 'post_processing' -name 'combined_*.nii.gz*' -type f`
+        fi
+    else
+        files_to_process=`find 'raw_outputs' -name 'combined_*.nii.gz*' -type f`
+    fi
+
+    for pm in $files_to_process; do
+            
+            pm_base_name=${pm##*/}
+            echo 'coregistering to anatomy -> '$pm_base_name
+            flirt \
+                -in $pm \
+                -ref "$anatomy_volume" \
+                -out 'alignments/'${pm_base_name%.*.*}'2FS_anat.nii.gz' \
+                -init $retmap2anat_mat \
+                -applyxfm
+            if [[ -f $ROI_mask ]]; then                
+                fslmaths 'alignments/'${pm_base_name%.*.*}'2FS_anat.nii.gz' -mas $ROI_mask 'alignments/'${pm_base_name%.*.*}'2FS_anat.nii.gz'
+            fi
+            ## left hemisphere
+            mri_vol2surf \
+                --hemi lh \
+                --mov 'alignments/'${pm_base_name%.*.*}'2FS_anat.nii.gz' \
+                --regheader $subject \
+                --projfrac-avg 0.0 0.75 0.001 \
+                --out 'surface_maps/'${pm_base_name%.*.*}'_lh.mgh' \
+                --surf white \
+                --out_type mgh
+
+            ## right hemisphere
+            mri_vol2surf \
+                --hemi rh \
+                --mov 'alignments/'${pm_base_name%.*.*}'2FS_anat.nii.gz' \
+                --regheader $subject \
+                --projfrac-avg 0.0 0.75 0.001 \
+                --out 'surface_maps/'${pm_base_name%.*.*}'_rh.mgh' \
+                --surf white \
+                --out_type mgh
+    done
+
+fi
+
+if [ "$debug_mode_flag" != "True" ]; then
+    echo 'delete intermediate steps'
+    rm -rf 'pre_processing'
+fi

src/aligned

@@ -0,0 +1 @@
+Subproject commit 72b577b301ad307f7f2d56e8218763b7a793597f

src/freesurfer

@@ -0,0 +1 @@
+Subproject commit 9bae995af237aed92437c4738c221176155e10df

src/structural

@@ -0,0 +1 @@
+Subproject commit 7f40a2fa64cadda2de1fe7d96b29d9f7d318e175

src/templatetransforms

@@ -0,0 +1 @@
+Subproject commit 98584d15796b62dce61ed88d06e72d76bff5511a