
proc CalculateLOCOFit { args } {
    
    global latticeIterFile latticeRMFile lteFile
    global dispersionFile variableFile
    global SVnumber badPointsLevel badPointsIterCounter bpmTiltIterCounter
    global rmDerivFile rmDerivInverseFile rmDifferenceFile twiFile
    global hrmMeasured vrmMeasured hrmIterFile vrmIterFile hrmDifference vrmDifference 
    global elementListVar elementListRM
    global hYes vYes nuxMeasured nuyMeasured corFraction solutionParamFileList
    global usePreviousSolution previousResultsFile resultsFile elementsUpdateFile variableArchiveFile
    global averageErrorFile allElementsFile specialElementsInputFile specialParamFile
    global svdProgramName useQsubLong recalcIterRMD recalcIterRMDErrorLevel
    global workDir tmpDir definitionFile verbose useQuadConstraints quadConstraintWeight adjustAverageGains
    global outputStatusDevice iterationsLogFile iterationsLogFileSDDS convergenceTolerance resErrorMax quadConstraintFile
    global convergenceTolerance2 iterations useParallelMultiplication matrixDivisions

    set time1 [exec date +%s]
    set requiredVarList [list latticeInputFile hrmMeasuredInput vrmMeasuredInput \
			     averLevel workDir beamlineName bigIter bigIterNumber \
			     elementListVar elementListRM optionsList nuxMeasured nuyMeasured useMeasuredDisp \
			     usePreviousSolution recalculateDerivative recalculateInverse \
			     continuePrevious nonlinLevel badPointsLevel useBpmWeight recombineDerivative \
			     SVnumber corFraction iterations resultsFile coupledMatrix filterDerivative \
			     definitionFile acceleratorCode calculateDerivOnly \
			     outputStatusDevice useKickFiles abortFile verbose]

    set resErrorList ""
    foreach var $requiredVarList {set $var ""}
    set args1 $args

    APSParseArguments $requiredVarList
    set missingArgs ""
    foreach var $requiredVarList {
	if ![string length [set $var]] {lappend missingArgs $var}
    }
    if [llength $missingArgs] {
	return -code error "There are missing arguments: $missingArgs"
    }
    if {![string length $optionsList]} {
        return -code error "No fit options given. \"options\" is empty."
    } else {
        global options
        array set options $optionsList
    }
    set rawLteFile $latticeInputFile

    #------ Not required variables:
    set nonrequiredVarList [list previousResultsFile kickFileX kickFileY \
				rawRMDeriv latticeParamFileList prepareMatricesOnly]
    foreach var $nonrequiredVarList {set $var ""}
    set prepareMatricesOnly 0
    APSParseArguments $nonrequiredVarList
    set lteFile $rawLteFile

    #------ definitionFile defines a lot of internal variables.
    if ![string length $definitionFile] {
	set definitionFile $workDir/initialDefinitions.tcl
	uplevel \#0 "set definitionFile $workDir/initialDefinitions.tcl"
    }
    if [catch {Fit_ReadDefinitionFile -definitionFile $definitionFile -workDir $workDir} result] {
	return -code error "Fit_ReadDefinitionFile: $result"
    }

    #------ check for possible conflicts in variables:
    global closedOrbit directRMCalc directDispersionCalc keepTunesOnTarget
    if {!$closedOrbit && $options(useEnergy)} {
	OutputStatusMessage "Warning: Energy fitting must be disabled when using non-closed trajectories. Disabling..."
	set options(useEnergy) 0
    }
    if !$options(useQuads) {
	if {$useQuadConstraints != 0} {
	    OutputStatusMessage "Warning: useQuadConstraints defined in $definitionFile is forced to 0 because there are no quads to fit."
	    set useQuadConstraints 0
	}
    }
    if !$options(fitDispersion) {
	if {[string compare $adjustAverageGains dispersion] == 0} {
	    OutputStatusMessage "Warning: adjustAverageGains is set to \"dispersion\" while no dispersion fitting is chosen. Will force to \"none\"."
	    set adjustAverageGains none
	}
    }
    if {$directRMCalc != 1} {
	if $directDispersionCalc {
	    OutputStatusMessage "Warning: directDispersionCalc should be disabled if directCalc!=1. Disabling..."
	    set directDispersionCalc 0
	}
    }
    if $keepTunesOnTarget {
	if !$options(useQuads) {
	    OutputStatusMessage "Warning: keepTunesOnTarget=1 while useQuads=0. Will force keepTunesOnTarget to 0."
	    set keepTunesOnTarget 0
	}
    }

    #------ Read the number of iterations in the previous runs or delete the log file if it's the first bigIter:
    if {$bigIter == 0} {
	catch {file delete $iterationsLogFile $iterationsLogFileSDDS $variableArchiveFile}
    } else {
	if [catch {exec sddscombine $iterationsLogFileSDDS -pipe=out -merge \
		       | sddsprocess -pipe -clip=0,1,inv \
		       | sdds2stream -pipe=in -col=Iteration} result] {
	    return -code error "Error reading $iterationsLogFileSDDS: $result"
	}
	set iterations0 [format %1.0f $result]
    }

    OutputStatusMessage "Calculations started..."

    if [catch {InitialPreparations -lteFile $lteFile \
		   -latticeParamFileList $latticeParamFileList \
		   -hrmMeasuredInput $hrmMeasuredInput -vrmMeasuredInput $vrmMeasuredInput \
		   -usePreviousSolution $usePreviousSolution -previousResultsFile $previousResultsFile \
		   -useMeasuredDisp $useMeasuredDisp -coupledMatrix $coupledMatrix \
		   -useBpmWeight $useBpmWeight -beamlineName $beamlineName \
		   -nonlinLevel $nonlinLevel -acceleratorCode $acceleratorCode \
		   -useKickFiles $useKickFiles -kickFileX $kickFileX -kickFileY $kickFileY \
		   -specialElementsInputFile $specialElementsInputFile -specialParamFile $specialParamFile \
		   -prepareMatricesOnly $prepareMatricesOnly -bigIter $bigIter} result] {
	return -code error "InitialPreparations: $result"
    }

### I think I do not need to ensure that recalculateDerivative recombineDerivative recalculateInverse are all to set to 1 when recalcIterRMD=1
###    if {$recalcIterRMD && (!$options(useQuads) && !$options(useSkew))} {
###	OutputStatusMessage "Warning: recalcIterRMD = 1 while no quad nor skew fit is chosen. Will reset recalcIterRMD to 0."
###	set recalcIterRMD
###    }
###    if $recalcIterRMD {foreach {recalculateDerivative recombineDerivative recalculateInverse} {1 1 1} {}}

    if !$prepareMatricesOnly  {
	if {[string length $nuxMeasured] && [string length $nuyMeasured]} {
	    OutputStatusMessage "[format "%s %10.4f %10.4f" "Measured Tunes are   " $nuxMeasured $nuyMeasured]"
	} else {
	    OutputStatusMessage "Warining: Measured betatron tunes are missing."
	}

	set badPointsIterCounter 0
	set bpmTiltIterCounter 0
	set residualError 10.0
	if $usePreviousSolution {
	    if ![catch {exec sdds2stream $previousResultsFile -para=ResidualError} result] {
		#------ Need lindex here because the file is multipage.
		set residualError [lindex $result 0]
	    }
	    if ![catch {exec sdds2stream $previousResultsFile -para=BadPointsIterCounter} result] {
		set badPointsIterCounter [expr int([lindex $result 0])]
	    }
	}
	if [catch {CalculateRMDifference -latticeFile $latticeIterFile -acceleratorCode $acceleratorCode \
		       -xyColumnFile $rmDifferenceFile -dispersionFile $dispersionFile -allElementsFile $allElementsFile \
		       -hrmMeasured $hrmMeasured -hrmIterFile $hrmIterFile -hrmDifference $hrmDifference \
		       -vrmMeasured $vrmMeasured -vrmIterFile $vrmIterFile -vrmDifference $vrmDifference \
		       -residualError $residualError -badPointsLevel $badPointsLevel -nonlinLevel $nonlinLevel \
		       -currentIteration -1 \
		       -averLevel $averLevel -averageErrorFile $averageErrorFile -coupledMatrix $coupledMatrix \
		       -useQuadConstraints $useQuadConstraints -quadConstraintWeight $quadConstraintWeight \
		       -useKickFiles $useKickFiles -kickFileX $kickFileX -kickFileY $kickFileY} rmsList] {
	    return -code error "CalculateRMDifference: $rmsList"
	}
	array set rmsArray $rmsList
	set residualError $rmsArray(rmsTotalWQ)
	set residualErrorPrev $residualError
	lappend resErrorList $residualError
	if {$bigIter == 0} {
	    #------ Copy only if it is the first bigIter
	    #------ 2plane file has the same BPMs in X and Y planes. It will be used in BuildResponseMatrixDeriv
	    foreach filename [list $hrmIterFile $vrmIterFile] {
		file copy -force $filename $filename.before
		if {$options(useBPMsTilt) || $options(useBPMsCoef)} {
		    file copy -force [file root $filename].2plane[file ext $filename] [file root $filename].2plane[file ext $filename].before
		}
	    }
	}
	if [catch {IterationResultsOutput -rmsList $rmsList -coupledMatrix $coupledMatrix -bigIter $bigIter \
		       -fitDispersion $options(fitDispersion) -iteration 0 -sddsFile $iterationsLogFileSDDS} result] {
	    return -code error "IterationResultsOutput: $result"
	}
	
	if [file exists $abortFile] {
	    return -code error "Iterrupted by user."
	}

	#------ Build the Response Matrix Derivative
	set rmdChanged 0
	if {$recalculateDerivative && $iterations} {
	    if [catch {BuildResponseMatrixDeriv -continuePrevious $continuePrevious \
			   -rmDerivFile $rmDerivFile.original -useQsubLong $useQsubLong -acceleratorCode $acceleratorCode \
			   -hrmFile $hrmIterFile.before -vrmFile $vrmIterFile.before \
			   -hrmFile2plane [file root $hrmIterFile].2plane[file ext $hrmIterFile].before \
			   -vrmFile2plane [file root $vrmIterFile].2plane[file ext $vrmIterFile].before \
			   -coupledMatrix $coupledMatrix \
			   -lteFile $lteFile -beamlineName $beamlineName -solutionParamFileList $solutionParamFileList \
			   -useKickFiles $useKickFiles -kickFileX $kickFileX -kickFileY $kickFileY \
			   -abortFile $abortFile} result] {
		return -code error "BuildResponseMatrixDeriv: $result"
	    }
	    file delete -force $rmDerivFile
	    exec ln -s $rmDerivFile.original $rmDerivFile
	    set rmdChanged 1
	}
	if $calculateDerivOnly {
	    return -code ok
	}
    }
    set currentRMDFile $rmDerivFile.original

    #------ Filter out the Response Matrix Derivative
    #------ Filtering is always done on .original or rawRMDeriv
    if {$filterDerivative && $iterations} {
	if {[string length $rawRMDeriv] && [file exists $rawRMDeriv]} {
	    set currentRMDFile $rawRMDeriv
	}
        if [catch {FilterResponseMatrixDeriv -inputFile $currentRMDFile -outputFile $rmDerivFile.filtered \
                       -fitDispersion $options(fitDispersion) \
                       -elementListVar $elementListVar -elementListRM $elementListRM -optionsList $optionsList \
                       -coupledMatrix $coupledMatrix -useTune $options(useTune) \
		       -useKickFiles $useKickFiles -kickFileX $kickFileX -kickFileY $kickFileY} result] {
            return -code error "FilterResponseMatrixDeriv: $result"
        }
	file delete -force $rmDerivFile
	exec ln -s $rmDerivFile.filtered $rmDerivFile
	set currentRMDFile $rmDerivFile.filtered
	set rmdChanged 1
    }
    
    #------ Recombine the Response Matrix Derivative with quad constraints
    #------ If rmdChanged, used currentRMDFile. If not, use source file from previous calculations.
    if {$recombineDerivative && $iterations} {
	if $useQuadConstraints {
	    if !$rmdChanged {set currentRMDFile [exec sdds2stream $rmDerivFile.QC -para=SourceFile]}
	    if [catch {RecombineResponseMatrixDeriv -input $currentRMDFile -output $rmDerivFile.QC \
			   -constraintWeight $quadConstraintWeight -constraintFile $quadConstraintFile \
			   -variableFile $variableFile -coupledMatrix $coupledMatrix \
			   -useQuadConstraints $useQuadConstraints} result] {
		return -code error "RecombineResponseMatrixDeriv: $result"
	    }
	    file delete -force $rmDerivFile
	    exec ln -s $rmDerivFile.QC $rmDerivFile
	} else {
	    #------ Recombining with useQuadConstraints=0 - removing quad constraints:
	    #------ Done only when the matrix was not recalculated:
	    if !$rmdChanged {
		set currentRMDFile [exec sdds2stream $rmDerivFile.QC -para=SourceFile]
		file delete -force $rmDerivFile
		exec ln -s $currentRMDFile $rmDerivFile
	    }
	}
    }
    
    if [file exists $abortFile] {return -code error "Iterrupted by user."}

    if {$recalculateInverse && $iterations} {
	if [catch {CalculateRMDerivInverse -programName $svdProgramName -SVnumber $SVnumber \
		       -rmDerivFile $rmDerivFile -rmDerivInverseFile $rmDerivInverseFile \
		       -rmDifferenceFile $rmDifferenceFile} result] {
	    return -code error "CalculateRMDerivInverse: $result"
	}
    }
    
    #------ Split Inverse matrix for parallel multiplication
    if $prepareMatricesOnly {
	if [catch {Fit_ReadDefinitionFile -definitionFile $definitionFile -workDir $workDir \
		       -variableList "useParallelMultiplication matrixDivisions" -abortOnError 0} result] {
	    return -code error "Fit_ReadDefinitionFile: $result"
	}
    }
    if [info exists useParallelMultiplication] {
	if $useParallelMultiplication {if ![info exists matrixDivisions] {set matrixDivisions 4}}
    } else {
	set useParallelMultiplication 0
    }
    if {$useParallelMultiplication && $bigIter == 0 && $iterations != 0} {
	set doSplit 0
	if $recalculateInverse {
	    set doSplit 1
	} else {
	    for {set i 1} {$i <= $matrixDivisions} {incr i} {if ![file exists ${rmDerivInverseFile}$i] {set doSplit 1}}
	}
	if $doSplit {
	    OutputStatusMessage "[exec date +%D--%T]: Splitting derivative matrix..."
	    set rmdRows [llength [join $elementListVar ]]
	    if $options(fitDispersion) {set rmdRows [expr $rmdRows + 1]}
	    set rmdRowLimit [expr round($rmdRows / $matrixDivisions) + 1]
	    exec sddsbreak $rmDerivInverseFile -pipe=out -rowlimit=$rmdRowLimit \
		| sddssplit -pipe=in -rootname=$rmDerivInverseFile -digits=1 -ext=
	}
    }
    set time2 [exec date +%s]
    OutputStatusMessage "Preparation took [expr $time2 - $time1] sec." -logFileOnly 1
    if $prepareMatricesOnly {
	OutputStatusMessage "Done: Matrices prepared for fit."
	return
    }
#
# ------- Starting iterations: ------------------------------------------------------------------------------------------
#

    if {$recalcIterRMD && (!$options(useQuads) && !$options(useSkew))} {
	OutputStatusMessage "Warning: recalcIterRMD = 1 while no quad nor skew fit is chosen. Will reset recalcIterRMD to 0."
	set recalcIterRMD 0
    }
    set pwd [pwd]
    file copy -force $hrmMeasured $hrmMeasured.original
    file copy -force $vrmMeasured $vrmMeasured.original
    if [catch {GetTolerance -bigIter $bigIter} convergenceTol] {
	return -code error "GetTolerance: $convergenceTol"
    }
    
    set convergenceCounter 0
    set divergenceCounter 0
    catch {file delete $workDir/vector_out}
    for {set I 0} {$I < $iterations} {incr I 1} {
        
	set time1 [exec date +%s]
	if {$recalcIterRMD && $I != 0} {
	    #------ We look at both residualError and vectorStdDev because when bpm noise is large, residualError may never get small, while
	    #------ vectorStdDev will get smaller
	    if {$residualError > $recalcIterRMDErrorLevel && $vectorStdDev > $recalcIterRMDErrorLevel} {
		if [catch {CalculateAndInvertRMD -lteFile $lteFile -beamlineName $beamlineName -solutionParamFileList $solutionParamFileList \
			       -coupledMatrix $coupledMatrix -continuePrevious $continuePrevious -acceleratorCode $acceleratorCode \
			       -useKickFiles $useKickFiles -kickFileX $kickFileX -kickFileY $kickFileY -abortFile $abortFile} result] {
		    return -code error "CalculateAndInvertRMD: $result"
		}
	    }
	}

	if {$bigIter != 0} {
	    OutputStatusMessage "-------> Loop $bigIter: Iteration number $I ($iterations0)"
	} else {
	    OutputStatusMessage "-------> Loop $bigIter: Iteration number $I "
	}

        #------ Matrix multiplication ---
	set time1 [exec date +%s]
	if $useParallelMultiplication {
	    if [catch {ParallelMultiplication -matrix1 $rmDerivInverseFile -matrix2 $rmDifferenceFile \
			   -matrixDivisions $matrixDivisions -output $workDir/vector_out} result] {
		return -code error "ParallelMultiplication: $result"
	    }
	} else {
	    if [catch {exec sddsmatrixmult $rmDerivInverseFile $rmDifferenceFile $workDir/vector_out} result] {
		return -code error "Error running sddsmatrixmult ($rmDerivInverseFile $rmDifferenceFile): $result"
	    }
	}
	set time2 [exec date +%s]
	OutputStatusMessage "Main Loop:: sddsmatrixmult took [expr $time2 - $time1] sec." -logFileOnly 1

        set vectorStdDev [exec sddsprocess $workDir/vector_out -pipe=out -process=MatrixDifference,standardDeviation,Std \
			      | sdds2stream -pipe=in -para=Std]
        
        #------ Updating the model after the iteration ---
        if [catch {UpdateModel -elementsUpdateFile $elementsUpdateFile -vectorOutFile $workDir/vector_out \
                       -corFraction $corFraction -variableFile $variableFile -allElementsFile $allElementsFile} result] {
            return -code error "UpdateModel: $result"
        }
        
        #------ Calculate response matrix difference ---
	set time3 [exec date +%s]
        if [catch {CalculateRMDifference -acceleratorCode $acceleratorCode \
                       -dispersionFile $dispersionFile -averageErrorFile $averageErrorFile \
                       -xyColumnFile $rmDifferenceFile -allElementsFile $allElementsFile \
                       -hrmMeasured $hrmMeasured -vrmMeasured $vrmMeasured \
                       -hrmDifference $hrmDifference -vrmDifference $vrmDifference \
                       -hrmIterFile $hrmIterFile -vrmIterFile $vrmIterFile \
                       -latticeFile $latticeIterFile -coupledMatrix $coupledMatrix \
                       -residualError $residualError \
                       -badPointsLevel $badPointsLevel -nonlinLevel $nonlinLevel \
                       -averLevel $averLevel -currentIteration $I \
		       -useQuadConstraints $useQuadConstraints -quadConstraintWeight $quadConstraintWeight \
		       -useKickFiles $useKickFiles -kickFileX $kickFileX -kickFileY $kickFileY} rmsList] {
	    OutputStatusMessage "Last variable vector std: $vectorStdDev."	    
	    OutputStatusMessage "CalculateRMDifference: $rmsList"
 	    if [catch {StopExecution -resultsFile $resultsFile -allElementsFile $allElementsFile \
                           -hrmIterFile $hrmIterFile -vrmIterFile $vrmIterFile \
			   -rmsList $rmsList -badPointsIterCounter $badPointsIterCounter} result] {
                return -code error "StopExecution(1): $result"
            }
	    OutputStatusMessage "Last variable vector std: $vectorStdDev."	    
            if [catch {StopExecution -resultsFile $resultsFile -allElementsFile $allElementsFile \
                           -hrmIterFile $hrmIterFile -vrmIterFile $vrmIterFile \
			   -rmsList $rmsList -badPointsIterCounter $badPointsIterCounter} result] {
                return -code error "StopExecution(2): $result"
            }
            return -code error "CalculateRMDifference: $rmsList"
        }
	set time4 [exec date +%s]
	OutputStatusMessage "Main Loop:: Matrix difference took [expr $time4 - $time3] sec." -logFileOnly 1
	array set rmsArray $rmsList
	set residualError $rmsArray(rmsTotalWQ)
	lappend resErrorList $residualError
        
        if [catch {IterationResultsOutput -rmsList $rmsList -vectorStdDev $vectorStdDev -iteration [expr $I + 1] \
                       -fitDispersion $options(fitDispersion) -coupledMatrix $coupledMatrix \
		       -sddsFile $iterationsLogFileSDDS -bigIter $bigIter} result] {
            return -code error "IterationResultsOutput: $result"
        }
        
	if [catch {ArchiveVariableResults -currentFile $allElementsFile.sdds -archiveFile $variableArchiveFile -iteration [expr $I + 1] \
		       -bigIter $bigIter -residualError $rmsArray(rmsTotalWQ)} result] {
	    return -code error "ArchiveVariableResults: $result"
	}

        if [file exists $abortFile] {
            if [catch {StopExecution -resultsFile $resultsFile -allElementsFile $allElementsFile \
                           -hrmIterFile $hrmIterFile -vrmIterFile $vrmIterFile \
			   -rmsList $rmsList -badPointsIterCounter $badPointsIterCounter} result] {
                return -code error "StopExecution(2): $result"
            }
            return -code error "Interrupted by user."
        }
	set time2 [exec date +%s]
	OutputStatusMessage "Main Loop:: Iteration $I took [expr $time2 - $time1] sec." -logFileOnly 1

	#------ Checking for convergence and for number diverging iterations:
	if [expr abs($residualError - $residualErrorPrev) < $convergenceTol] {incr convergenceCounter}
	if {$residualError > $residualErrorPrev} {incr divergenceCounter}
	if {$convergenceCounter > 1} {
	    OutputStatusMessage "Iterations converged."
	    break
	}
	if {$divergenceCounter > 5} {
	    OutputStatusMessage "Iterations are diverging - iterations are interrupted."
	    break
	}

	set residualErrorPrev $residualError
    }
    
    if [catch {StopExecution -resultsFile $resultsFile -allElementsFile $allElementsFile \
                   -hrmIterFile $hrmIterFile -vrmIterFile $vrmIterFile \
		   -rmsList $rmsList -badPointsIterCounter $badPointsIterCounter} result] {
        return -code error "StopExecution(3): $result"
    }
    OutputStatusMessage "Iterations are done."
    return $resErrorList
}

#------------------------------------------------------------------------------------------------------------------------
proc ArchiveVariableResults {args} {
    global tmpDir hrmIterFile vrmIterFile twiFile iterationsLogFileSDDS
    APSParseArguments {currentFile archiveFile iteration bigIter residualError}
    set tmpRoot $tmpDir/[APSTmpString]-archive
    set currentFileList [list $currentFile $hrmIterFile $vrmIterFile $twiFile]
    set archiveFileList [list $archiveFile $hrmIterFile.iterations $vrmIterFile.iterations $twiFile.iterations]

    #--- Find total iteration:
    if {$bigIter == 0} {
	set outputIter $iteration
    } else {
	if [catch {exec sddscombine $iterationsLogFileSDDS -pipe=out -merge \
		       | sddsprocess -pipe -filter=col,BigIteration,0,[expr $bigIter - 1] -clip=0,1,inv \
		       | sdds2stream -pipe=out -col=Iteration} iterations0] {
	    return -code error "Error reading iterations0: $iterations0"
	}
	set outputIter [expr $iteration + $iterations0]
    }

    if ![file exists $archiveFile] {
	foreach filename1 $currentFileList filename2 $archiveFileList {
	    exec sddsprocess $filename1 $filename2 -redef=para,BigIteration,$bigIter,type=long -redef=para,Iteration,$outputIter,type=long \
		-redef=para,ResidualError,$residualError -nowarning
	}
    } else {
	foreach filename1 $currentFileList filename2 $archiveFileList {
	    exec sddsprocess $filename1 $tmpRoot.1 -redef=para,BigIteration,$bigIter,type=long -redef=para,Iteration,$outputIter,type=long \
		-redef=para,ResidualError,$residualError -nowarning
	    exec sddscombine $filename2 $tmpRoot.1 $tmpRoot.2 -overWrite
	    file copy -force $tmpRoot.2 $filename2
	    file delete $tmpRoot.1 $tmpRoot.2
	}
    }
}
#------------------------------------------------------------------------------------------------------------------------
proc ParallelMultiplication {args} {
    global tmpDir qsubCommandShort qsubRespProcCommandShort queueSystemNameShort submissionPause
    global qsubCommand qsubRespProcCommand queueSystemName workDir
    set qsubCommand $qsubCommandShort
    set qsubRespProcCommand $qsubRespProcCommandShort
    set queueSystemName $queueSystemNameShort
    set submissionPause 0
    APSParseArguments {matrix1 matrix2 matrixDivisions output}
    set tmpRoot $workDir/[APSTmpString]-matrMult
    set deleteFiles ""
    set commandList ""
    set jobList ""
    set workDir [file dirname $output]
    for {set i 1} {$i <= $matrixDivisions} {incr i} {
	set jobName matrixMult$i
	set doneFile $workDir/matrixMult$i.done
	catch {file delete $doneFile}
	set command "mkdir -p $tmpDir; sddsmatrixmult $matrix1$i $matrix2 $tmpRoot.$i; touch $doneFile"
	lappend commandList $command
	lappend jobNameList $jobName
	lappend doneFileList $doneFile
	lappend outputFiles $tmpRoot.$i
	lappend deleteFiles $doneFile $tmpRoot.$i
    }
    if [catch {Fit_SubmitJobs -commandList $commandList -jobNameList $jobNameList \
		   -doneFileList $doneFileList -verbose 0} returnList] {
	return -code error "Fit_SubmitJobs: $returnList"
    }
    set pidList [lindex $returnList 0]
    set logFileList [lindex $returnList 1]
    set commandFileList [lindex $returnList 2]
    set dotOFileList [lindex $returnList 3]
    set waitTime 25
    set waitInterval 1.0
    if [catch {Fit_WaitForTasks -pidList $pidList -jobNameList $jobNameList -commandList $commandList \
		   -logFileList $logFileList -waitTime $waitTime \
		   -doneFileList $doneFileList -waitInterval $waitInterval -usePopupWindow 0 -verbose 0 \
		   -useQsub 1 -abortFile dummyMult} result] {
	return -code error "Fit_WaitForTasks: $result"
    }
    #------ Sometimes files don't get written fast enough, need to check the size:
    set keepWaiting 1
    set counter 0
    while {$keepWaiting == 1} {
	set keepWaiting 0
	foreach filename $outputFiles {if {[file size $filename] == 0} {set keepWaiting 1}}
	if $keepWaiting {after 1000; puts stdout "----------->>> Extra waiting for matrix multiplication files to appear..."}
	incr counter
	if {$counter > 20} {return -code error "Error: Parallel matrix multiplication failed to create files."}
    }
    if [catch {eval exec sddscombine $outputFiles $output -merge -overWrite} result] {
	return -code error "Error combining matrix parts ($outputFiles): $result"
    }
    eval file delete [concat $deleteFiles $commandFileList $dotOFileList]
}
#------------------------------------------------------------------------------------------------------------------------
#------------------------------------------------------------------------------------------------------------------------
#------ Procedures related to initial preparation of the problem (mainly creation of files) --------------------------------
#------------------------------------------------------------------------------------------------------------------------
#------------------------------------------------------------------------------------------------------------------------

proc CalculateAndInvertRMD {args} {
    global rmDerivFile useQsubLong latticeIterFile hrmIterFile vrmIterFile useQuadConstraints quadConstraintWeight
    global svdProgramName workDir SVnumber rmDerivInverseFile rmDifferenceFile variableFile svdProgramName SVnumber 
    global elementListVar useParallelMultiplication matrixDivisions options allElementsFile

    APSParseArguments {lteFile beamlineName solutionParamFileList coupledMatrix continuePrevious acceleratorCode abortFile useKickFiles kickFileX kickFileY}

    set localParamFileList $solutionParamFileList
    if $options(useQuads) {lappend $localParamFileList $allElementsFile.Quad}
    if $options(useSkew)  {lappend $localParamFileList $allElementsFile.Skew}
    if [catch {BuildResponseMatrixDeriv -continuePrevious $continuePrevious \
		   -rmDerivFile $rmDerivFile.current -useQsubLong $useQsubLong -acceleratorCode $acceleratorCode \
		   -hrmFile $hrmIterFile.before -vrmFile $vrmIterFile.before \
		   -hrmFile2plane [file root $hrmIterFile].2plane[file ext $hrmIterFile].before \
		   -vrmFile2plane [file root $vrmIterFile].2plane[file ext $vrmIterFile].before \
		   -coupledMatrix $coupledMatrix \
		   -lteFile $lteFile -beamlineName $beamlineName -solutionParamFileList $solutionParamFileList \
		   -useKickFiles $useKickFiles -kickFileX $kickFileX -kickFileY $kickFileY \
		   -abortFile dummyRMD -verboseLocal 0} result] {
	return -code error "BuildResponseMatrixDeriv: $result"
    }
    file delete -force $rmDerivFile
    exec ln -s $rmDerivFile.current $rmDerivFile

    if $useQuadConstraints {
	if [catch {RecombineResponseMatrixDeriv -input $rmDerivFile.current -output $rmDerivFile.QC \
		       -constraintWeight $quadConstraintWeight -constraintFile $quadConstraintFile \
		       -variableFile $variableFile -coupledMatrix $coupledMatrix \
		       -useQuadConstraints $useQuadConstraints} result] {
	    return -code error "RecombineResponseMatrixDeriv: $result"
	}
	file delete -force $rmDerivFile
	exec ln -s $rmDerivFile.QC $rmDerivFile
    }

    if [catch {CalculateRMDerivInverse -programName $svdProgramName -SVnumber $SVnumber \
		   -rmDerivFile $rmDerivFile -rmDerivInverseFile $rmDerivInverseFile \
		   -rmDifferenceFile $rmDifferenceFile -verbose 0} result] {
	return -code error "CalculateRMDerivInverse: $result"
    }

    if $useParallelMultiplication {
	set rmdRows [llength [join $elementListVar ]]
	if $options(fitDispersion) {set rmdRows [expr $rmdRows + 1]}
	set rmdRowLimit [expr round($rmdRows / $matrixDivisions) + 1]
	exec sddsbreak $rmDerivInverseFile -pipe=out -rowlimit=$rmdRowLimit \
	    | sddssplit -pipe=in -rootname=$rmDerivInverseFile -digits=1 -ext=
    }

}

#------------------------------------------------------------------------------------------------------------------------

proc InitialPreparations { args } {

    APSParseArguments {lteFile previousResultsFile usePreviousSolution coupledMatrix \
			   useBpmWeight useMeasuredDisp hrmMeasuredInput vrmMeasuredInput \
			   nonlinLevel latticeParamFileList beamlineName acceleratorCode \
			   useKickFiles kickFileX kickFileY \
			   specialElementsInputFile specialParamFile prepareMatricesOnly bigIter}

    global tmpDir useSigmaWeightFile
    global twiFile dispersionFile elementsUpdateFile variableFile elementOrderFile allElementsFile
    global hrmMeasured vrmMeasured workDir weightFile manualWeightFile weightCorrectorsFile averageErrorFile 
    global latticeIterFile latticeRMFile dispColumn dispMeasuredFile nonlinFile sigmaWeightFile
    global hYes vYes elementListRM options indexList variableTypeList 
    global zeroBpmTol solutionParamFileList bRho

    set hYes [llength [lindex $elementListRM 0]]
    set vYes [llength [lindex $elementListRM 1]]

    global hrmIterFile vrmIterFile
    if {$bigIter == 0} {
	catch {file delete $hrmIterFile.iterations $vrmIterFile.iterations $twiFile.iterations \
		   $hrmIterFile.before $hrmIterFile.after $vrmIterFile.before $vrmIterFile.after $weightFile}
    }
    catch {eval file delete $hrmIterFile.flhrm $hrmIterFile.flhrm.XRM $hrmIterFile.flvrm $hrmIterFile.flvrm.XRM}

    if [catch {CreateElementOrderFile -elementOrderFile $elementOrderFile -acceleratorCode $acceleratorCode \
                   -twiFile $twiFile -lteFile $lteFile -latticeParamFileList $latticeParamFileList \
                   -beamlineName $beamlineName -dispColumn $dispColumn -useKickFiles $useKickFiles \
		   -kickFileX $kickFileX -kickFileY $kickFileY} result] {
        return -code error "CreateElementOrderFile: $result"
    }
    
    #------ Create parameter files: 
    if [catch {CreateParameterFiles -twiFile $twiFile -dispersionFile $dispersionFile \
                   -elementsUpdateFile $elementsUpdateFile -allElementsFile $allElementsFile \
                   -usePreviousSolution $usePreviousSolution -previousResultsFile $previousResultsFile \
                   -useMeasuredDisp $useMeasuredDisp -variableFile $variableFile \
                   -averageErrorFile $averageErrorFile -elementOrderFile $elementOrderFile \
                   -lteFile $lteFile -latticeParamFileList $latticeParamFileList \
		   -beamlineName $beamlineName -coupledMatrix $coupledMatrix -dispColumn $dispColumn \
                   -nonlinFile $nonlinFile -nonlinLevel $nonlinLevel \
		   -acceleratorCode $acceleratorCode \
		   -useKickFiles $useKickFiles -kickFileX $kickFileX -kickFileY $kickFileY \
		   -specialElementsInputFile $specialElementsInputFile -specialParamFile $specialParamFile} result] {
        return -code error "CreateParameterFiles: $result"
    }
    if $prepareMatricesOnly {return}

    #------ Prepare measured response files:
    if [catch {PrepareMeasuredResponseFiles -hrmMeasuredInput $hrmMeasuredInput -vrmMeasuredInput $vrmMeasuredInput \
		   -hrmMeasured $hrmMeasured -vrmMeasured $vrmMeasured -elementOrderFile $elementOrderFile \
		   -coupledMatrix $coupledMatrix -useKickFiles $useKickFiles -kickFileX $kickFileX -kickFileY $kickFileY \
		   -dispMeasuredFile $dispMeasuredFile -dispColumn $dispColumn} result] {
	return -code error "PrepareMeasuredResponseFiles: $result"
    }

    set latticeIterFile $workDir/iterations.ele
    #------ Iteration ele file --- Include all files in case the variables exist in the prevSolFile
    set solutionParamFileList $latticeParamFileList
    lappend solutionParamFileList $allElementsFile.Quad
    lappend solutionParamFileList $allElementsFile.Skew
    lappend solutionParamFileList $allElementsFile.hCorr $allElementsFile.vCorr
    lappend solutionParamFileList $allElementsFile.hCorrTilt $allElementsFile.vCorrTilt
    if $options(useSpecElem) {lappend solutionParamFileList $specialParamFile}
    #------ Making iterations file...
    set energy [format %10.4e [expr $bRho * 299792458e-6]]
    if [catch {Fit_GenerateElegantFileFromLTE1 -eleOutputFile $latticeIterFile \
		   -run_setup [list "lattice $lteFile use_beamline $beamlineName default_order 2 p_central_mev $energy"] \
		   -load_parameters [list "filename_list \"$solutionParamFileList\""] \
		   -twiss_output [list "filename $twiFile"] \
		   -run_control [list "n_steps 1"] \
		   -bunched_beam [list "n_particles_per_bunch 1"] \
	       } result] {
	return -code error "Fit_GenerateElegantFileFromLTE1: $result"
    }
    
    #------ Takes manually created file with user-entered bpm wights and creates a file for the program's use.
    if [catch {MakeWeightFile -useBpmWeight $useBpmWeight -manualWeightFile $manualWeightFile -weightFile $weightFile \
                   -hrmFile $hrmMeasured -vrmFile $vrmMeasured -weightCorrectorsFile $weightCorrectorsFile \
                   -coupledMatrix $coupledMatrix} result] {
        return -code error "makeWeightFiles: $result"
    }
}


#-------------------------------------------------------------------------------------------------------------------------
proc PrepareMeasuredResponseFiles {args} {
    global zeroBpmTol options
    APSParseArguments {hrmMeasuredInput vrmMeasuredInput hrmMeasured vrmMeasured elementOrderFile coupledMatrix \
			   useKickFiles kickFileX kickFileY dispColumn}

    #------ Check for dispColumn presence:
    if $options(fitDispersion) {
	if {[lsearch [exec sddsquery -col $hrmMeasuredInput] $dispColumn] == -1} {
	    return -code error "Error: file $hrmMeasuredInput does not have $dispColumn column while fitDispersion is requested."
	}
    }

    file copy -force $hrmMeasuredInput $hrmMeasured
    file copy -force $vrmMeasuredInput $vrmMeasured
    #------ Filter according to element lists:
    if [catch {FilterMeasuredResponseFiles -hrmMeasured $hrmMeasured -vrmMeasured $vrmMeasured -elementOrderFile $elementOrderFile \
		   -coupledMatrix $coupledMatrix -useKickFiles $useKickFiles -kickFileX $kickFileX -kickFileY $kickFileY \
		   -dispColumn $dispColumn} result] {
        return -code error "FilterMeasuredResponseFiles: $result"
    }

    #------ Sort response files according to twiFile order and look for zero reading BPMs:
    if [catch {SortMeasuredResponseFiles -hrmMeasured $hrmMeasured -vrmMeasured $vrmMeasured -elementOrderFile $elementOrderFile \
		   -coupledMatrix $coupledMatrix} result] {
        return -code error "SortMeasuredResponseFiles ($hrmMeasured -- $vrmMeasured): $result"
    }
   
    #------ Detect zero reading BPMs to issue warning only:
    if [catch {DetectZeroBPMs -hrmMeasured $hrmMeasured -vrmMeasured $vrmMeasured -zeroBpmTol $zeroBpmTol -dispColumn $dispColumn} result] {
	return -code error "DetectZeroBPMs: $result"
    }

}
#-------------------------------------------------------------------------------------------------------------------------

proc CreateElementOrderFile {args} {

    global tmpDir options elementListRM
    APSParseArguments {elementOrderFile twiFile lteFile latticeParamFileList beamlineName dispColumn \
			    acceleratorCode useKickFiles kickFileX kickFileY}
    
    set tmpRoot $tmpDir/[APSTmpString]-elementOrder
    set sddsLatticeFile $tmpRoot.lattice.sdds
    set deleteFiles ""
    set latticeOptions "-lteFile $lteFile -beamlineName $beamlineName -latticeParamFileList \"$latticeParamFileList\""
    if [catch {eval CalculateResponseMatrix -getSDDSLattice 1 -workDir $tmpDir \
                   -sddsLatticeFile $sddsLatticeFile -acceleratorCode $acceleratorCode \
		   $latticeOptions} result] {
        return -code error "CalculateResponseMatrix: $result"
    }
    exec sddsconvert $sddsLatticeFile -pipe=out -retain=col,ElementName,ElementType,ElementOccurence,TagName \
	| sddsprocess -pipe=in $elementOrderFile "-match=col,ElementType=*KICK*,ElementType=*MON*,|,ElementType=*QUAD,|,ElementType=*BEND,|,ElementType=MULT,|" 
    lappend deleteFiles $sddsLatticeFile
    set combineFiles $elementOrderFile

    if $useKickFiles {
	#--- Remove names from kick files that are already in the order file:
	exec sddscollapse $kickFileX -pipe=out \
	    | sddsconvert -pipe -retain=col,ColumnName -rename=col,ColumnName=ElementName \
	    | sddsprocess -pipe -print=col,ElementType,MEASUREMENTX -edit=col,TagName,ElementName, -redef=col,ElementOccurence,1,type=long \
	    | sddsselect -pipe=in $elementOrderFile $tmpRoot.x -match=TagName -inv
	exec sddscollapse $kickFileY -pipe=out \
	    | sddsconvert -pipe -retain=col,ColumnName -rename=col,ColumnName=ElementName \
	    | sddsprocess -pipe -print=col,ElementType,MEASUREMENTY -edit=col,TagName,ElementName, -redef=col,ElementOccurence,1,type=long \
	    | sddsselect -pipe=in $elementOrderFile $tmpRoot.y -match=TagName -inv
	set combineFiles "$tmpRoot.x $tmpRoot.y $elementOrderFile"
	lappend deleteFiles $tmpRoot.x $tmpRoot.y $tmpRoot.order
    }
    #--- Dispersion needs to go after kickFiles so that dispersion column was after trajectory columns:
    if $options(fitDispersion) {
	exec sddsmakedataset $tmpRoot.D -col=ElementName,type=string -data=$dispColumn \
	    -col=ElementType,type=string -data=DISPERSION \
	    -col=TagName,type=string -data=$dispColumn -col=ElementOccurence,type=long -data=1
	lappend combineFiles $tmpRoot.D
	lappend deleteFiles $tmpRoot.D
    }
    if [llength $combineFiles] {
	eval exec sddscombine $combineFiles $tmpRoot -merge
	file copy -force $tmpRoot $elementOrderFile
	lappend deleteFiles $tmpRoot
    }
    eval file delete $deleteFiles
}

#-------------------------------------------------------------------------------------------------------------------------
proc CompareLists {args} {
    APSParseArguments {list0 subList}
    set missingElements ""
    foreach element $subList {if {[lsearch -exact $list0 $element] == -1} {lappend missingElements $element}}
    return $missingElements
}
#-------------------------------------------------------------------------------------------------------------------------
proc CheckElementListsSanity {args} {
    global options ignoreBpmListDiffs
    APSParseArguments {listRM listVar}

    #------ Make sure there are no correctors or bpms that are in variables list but not in RM list:
    set hcorrRMIndex 0
    set vcorrRMIndex 1
    set hcorrVarIndex 1
    set vcorrVarIndex 2
    set hcorrTiltRMIndex 0
    set vcorrTiltRMIndex 1
    set hcorrTiltVarIndex 7
    set vcorrTiltVarIndex 8
    set xBpmRMIndex 2
    set yBpmRMIndex 3
    set xBpmVarIndex 3
    set yBpmVarIndex 4
    set bpmTiltVarIndex 9
    set bpmCoefVarIndex 10
    set energyVarIndex 5
    set rootList ""
    set list0List ""
    set subListList ""
    if $options(useCorrs) {
	lappend rootList hCorr vCorr
	lappend list0List [lindex $listRM $hcorrRMIndex] [lindex $listRM $vcorrRMIndex]
	lappend subListList [lindex $listVar $hcorrVarIndex] [lindex $listVar $vcorrVarIndex]
    }
    if $options(useCorrsTilt) {
	lappend rootList hCorrTilt vCorrTilt
	lappend list0List [lindex $listRM $hcorrTiltRMIndex] [lindex $listRM $vcorrTiltRMIndex]
	lappend subListList [lindex $listVar $hcorrTiltVarIndex] [lindex $listVar $vcorrTiltVarIndex]
    }
    if $options(useBPMs) {
	lappend rootList xBpm yBpm
	lappend list0List [lindex $listRM $xBpmRMIndex] [lindex $listRM $yBpmRMIndex]
	lappend subListList [lindex $listVar $xBpmVarIndex] [lindex $listVar $yBpmVarIndex]
    }
    if $options(useEnergy) {
	lappend rootList Energy
	lappend list0List [lindex $listRM $hcorrRMIndex]
	lappend subListList [lindex $listVar $energyVarIndex]
    }
    set error 0
    foreach list0 $list0List subList $subListList root $rootList {
	set missingElements [CompareLists -list0 $list0 -subList $subList]
	if [string length $missingElements] {
	    OutputStatusMessage "Error: There are elements in Variables list that are not in RM list ($root): $missingElements"
	    set error 1
	}
    }

    #--- Check if xBpm and yBpm lists are the same:
    foreach index1 [list $xBpmRMIndex $yBpmRMIndex] index2 [list $yBpmRMIndex $xBpmRMIndex] {
	set list0 [lindex $listRM $index1]
	set list1 [lindex $listRM $index2]
	set list0_NoPlane ""
	set list1_NoPlane ""
	foreach bpm $list0 {regsub #\[XY\] $bpm "" bpm_NoPlane; lappend list0_NoPlane $bpm_NoPlane}
	foreach bpm $list1 {regsub #\[XY\] $bpm "" bpm_NoPlane; lappend list1_NoPlane $bpm_NoPlane}
	set missingElements [CompareLists -list0 $list0_NoPlane -subList $list1_NoPlane]
	if [string length $missingElements] {
	    if $ignoreBpmListDiffs {
		OutputStatusMessage "Warning: X and Y BPM lists are not equal: $missingElements"
	    } else {
		OutputStatusMessage "Error: X and Y BPM lists are not equal: $missingElements"
		set error 1
	    }
	}
    }

    #--- Check if bpmTilt/bpmCoef are equal to xBpm
    set rootList ""
    set list0List ""
    set subListList ""
    set bpmCompList [concat [lindex $listVar $xBpmVarIndex] [lindex $listVar $yBpmVarIndex]]
    if $options(useBPMsTilt) {
	lappend rootList bpmTilt
	lappend list0List $bpmCompList
	lappend subListList [lindex $listVar $bpmTiltVarIndex]
    }
    if $options(useBPMsCoef) {
	lappend rootList bpmCoef
	lappend list0List $bpmCompList
	lappend subListList [lindex $listVar $bpmCoefVarIndex]
    }
    foreach list0 $list0List list1 $subListList root $rootList {
	set list0_NoPlane ""
	set list1_NoPlane ""
	foreach bpm $list0 {regsub #\[XY\] $bpm "" bpm_NoPlane; lappend list0_NoPlane $bpm_NoPlane}
	foreach bpm $list1 {regsub #\[XY\] $bpm "" bpm_NoPlane; lappend list1_NoPlane $bpm_NoPlane}
	set missingElements [CompareLists -list0 $list0_NoPlane -subList $list1_NoPlane]
	if [string length $missingElements] {
	    if $ignoreBpmListDiffs {
		OutputStatusMessage "Warning: There are elements in Variables list that are not in RM list ($root): $missingElements"
	    } else {
		OutputStatusMessage "Error: There are elements in Variables list that are not in RM list ($root): $missingElements"
		set error 1
	    }
	}
    }

    if $error {return -code error "There are missing elements"}
}
#-------------------------------------------------------------------------------------------------------------------------
#------ Creates elements files and energy.sdds and noElements.sdds and dispersion.sdds files

proc CreateParameterFiles { args } {

    APSParseArguments { twiFile dispersionFile usePreviousSolution previousResultsFile elementsUpdateFile \
                            useMeasuredDisp variableFile averageErrorFile elementOrderFile \
                            allElementsFile coupledMatrix dispColumn nonlinLevel nonlinFile \
                            lteFile latticeParamFileList beamlineName acceleratorCode \
			    useKickFiles kickFileX kickFileY specialElementsInputFile specialParamFile}
    
    global options
    global workDir tmpDir
    global elementListRM elementListVar indexList variableTypeList
    set tmpRoot $tmpDir/[APSTmpString]-createParamFiles
    set deleteFiles ""

    catch {eval file delete [glob -nocomplain -- $allElementsFile.*] [glob -nocomplain -- $variableFile.*]}

    if [catch {CheckElementListsSanity -listRM $elementListRM -listVar $elementListVar} result] {
	return -code error "CheckElementListsSanity: $result"
    }

    set NVars 11
    set indexList [list useQuads useCorrs useCorrs useBPMs useBPMs useEnergy useSkew useCorrsTilt useCorrsTilt \
                       useBPMsTilt useBPMsCoef]
    set variableTypeList [list Quad hCorr vCorr xBpm yBpm Energy Skew hCorrTilt vCorrTilt bpmTilt bpmCoef]
    set parameterList [list K1 CALIBRATION CALIBRATION XCALIBRATION YCALIBRATION ENERGY TILT TILT TILT TILT COEF]
    set parameterModeList [list differential fractional fractional fractional fractional fractional differential \
                               differential differential differential differential]
    
    if $options(fitDispersion) {
        if [llength [lindex $elementListVar 1]] {
            set elementListVar [lreplace $elementListVar 1 1 [concat [lindex $elementListVar 1] $dispColumn]]
        }
    }

    set allElementsList $elementListVar
    if $usePreviousSolution {
        #------ Make allElementsList which contains all present variables plus any variables from $previousResultsFile
	#------ that are not used presently
	for {set i 0} {$i < $NVars} {incr i} {
	    set pageID [lindex $variableTypeList $i]
	    if {[exec sddsprocess $previousResultsFile -pipe=out -nowarning -match=para,PageID=$pageID \
		     | sdds2stream -pipe=in -rows=bare] != 0} {
		set tmpfile $tmpRoot.$pageID.1
		set variableList [lindex $allElementsList $i]
		exec sddsmakedataset $tmpfile -col=TagName,type=string -data=[join $variableList ,]
		lappend deleteFiles $tmpfile
		set additionalElements [join [exec sddsprocess $previousResultsFile -pipe=out -nowarning -match=para,PageID=$pageID \
						  | sddsselect -pipe $tmpfile -match=TagName -inv \
						  | sdds2stream -pipe=in -col=TagName] ]
		if [llength $additionalElements] {
		    set elementList [concat $variableList $additionalElements]
		    set allElementsList [lreplace $allElementsList $i $i $elementList]
		}
	    }
	}
    }

    #------ Creating zero allElements and variables files for everything regardless of options
    foreach filename [list $variableFile $allElementsFile] bigElementList [list $elementListVar $allElementsList] \
        fileList [list varFileList elemFileList] {
            set $fileList ""
            for {set i 0} {$i < $NVars} {incr i} {
                set index         [lindex $indexList $i]
                set parameter     [lindex $parameterList $i]
                set parameterMode [lindex $parameterModeList $i]
                set pageID        [lindex $variableTypeList $i]
                set elementList   [lindex $bigElementList $i]
		#------ "sddsmakedataset -data=" creates 1-row file while "sddsmakedataset -data" creates 0-row file
		if [llength $elementList] {
		    exec sddsmakedataset -pipe=out -col=TagName,type=string -data=[join $elementList ,] \
			| sddsprocess -pipe=in $filename.$pageID -print=para,PageID,$pageID \
			-def=col,ParameterValue,0.0 \
			-print=col,ElementParameter,$parameter \
			-print=col,ParameterMode,$parameterMode \
			-print=col,VariableType,$pageID
		    if [catch {Fit_MakeElementNameFromTagName -input $filename.$pageID} result] {
			return -code error "Fit_MakeElementNameFromTagName: $result"
		    }
		} else {
		    exec sddsmakedataset $filename.$pageID -para=PageID,type=string -data=$pageID \
			-col=TagName,type=string -data \
			-col=ParameterValue,type=double -data \
			-col=ElementParameter,type=string -data \
			-col=ParameterMode,type=string -data \
			-col=VariableType,type=string -data \
			-col=ElementName,type=string -data
		}
                #------ don't sort quads and skews because the issues with the OPTIM(quad name and quad variable are different!)
                if {!(![string compare $index useQuads] || ![string compare $index useSkew])} {
                    if [exec sdds2stream $filename.$pageID -rows=bare] {
                        if [catch {ResortSingleFile -orderFile $elementOrderFile -orderColumn TagName \
                                       -file2sort $filename.$pageID -column2sort TagName} result] {
                            return -code error "ResortSingleFile($pageID): $result"
                        }
                    }
                }
                lappend $fileList $filename.$pageID
            }
        }

    #------ Different types of elements have different parameters for tilt errors (*BEND -- ETILT; *QUAD -- TILT)
    #------ Correct elements.Skew and variables.Skew files:
    exec sddsmakedataset $tmpRoot.k1Types -col=ElementType,type=string -data=CCBEND,CSBEND,RBEN,SBEN,KQUAD,KQUSE,QUAD,MULT \
	-col=ElementParameter,type=string -data=K1,K1,K1,K1,K1,K1,K1,KNL
    exec sddsmakedataset $tmpRoot.tiltTypes -col=ElementType,type=string -data=CCBEND,CSBEND,RBEN,SBEN,KQUAD,KQUSE,QUAD \
	-col=ElementParameter,type=string -data=ETILT,ETILT,ETILT,ETILT,TILT,TILT,TILT
    foreach filename [list $variableFile $allElementsFile] {
	exec sddsxref $filename.Quad $elementOrderFile -pipe=out -take=ElementType -match=TagName -nowarning \
	    | sddsxref -pipe $tmpRoot.k1Types -replace=col,ElementParameter -match=ElementType -reuse -leave=* \
	    | sddsconvert -pipe=in $tmpRoot.quadNew -del=col,ElementType
	file copy -force $tmpRoot.quadNew $filename.Quad
	exec sddsxref $filename.Skew $elementOrderFile -pipe=out -take=ElementType -match=TagName -nowarning \
	    | sddsxref -pipe $tmpRoot.tiltTypes -replace=col,ElementParameter -match=ElementType -reuse -leave=* \
	    | sddsconvert -pipe=in $tmpRoot.skewNew -del=col,ElementType
	file copy -force $tmpRoot.skewNew $filename.Skew
    }
    lappend deleteFiles $tmpRoot.k1Types $tmpRoot.quadNew $tmpRoot.tiltTypes $tmpRoot.skewNew 

    #------ Creating file for special variables (for example, average gains)...
    set elementNameColumn ""
    set parameterValueColumn ""
    set elementParameterColumn ""
    set parameterModeColumn ""
    set variableTypeColumn ""
    if $options(useAverGains) {
	lappend elementNameColumn AverGainX AverGainY
	lappend parameterValueColumn 0 0
	lappend elementParameterColumn AVERGAIN AVERGAIN
	lappend parameterModeColumn differential differential
	lappend variableTypeColumn special special
	lappend elementOccurenceColumn 1 1
    }
    if $options(fitPathLength) {
	lappend elementNameColumn MALIN
	lappend tagNameColumn MALIN#1
	lappend parameterValueColumn 0
	lappend elementParameterColumn DZ
	lappend parameterModeColumn differential
	lappend variableTypeColumn special
	lappend elementOccurenceColumn 1
    }
    if $options(useSpecElem) {
	set specialNameList [join [exec sdds2stream $specialElementsInputFile -col=TagName] ]
	set specialParameterList [join [exec sdds2stream $specialElementsInputFile -col=ElementParameter] ]
	foreach specialName $specialNameList specialParameter $specialParameterList {
	    foreach {elementName elementOccurence} [Fit_MakeElementNameFromTagNameSingle -tagName $specialName] {}
	    lappend elementNameColumn $elementName
	    lappend tagNameColumn $specialName
	    lappend parameterValueColumn 0
	    lappend elementParameterColumn $specialParameter
	    lappend parameterModeColumn differential
	    lappend variableTypeColumn special
	    lappend elementOccurenceColumn $elementOccurence
	}
    }
    if [llength $elementNameColumn] {
	exec sddsmakedataset $variableFile.special \
	    -col=ElementName,type=string -data=[join $elementNameColumn ,] \
	    -col=TagName,type=string -data=[join $tagNameColumn ,] \
	    -col=ParameterValue,type=double -data=[join $parameterValueColumn ,] \
	    -col=ElementParameter,type=string -data=[join $elementParameterColumn ,] \
	    -col=ParameterMode,type=string -data=[join $parameterModeColumn ,] \
	    -col=VariableType,type=string -data=[join $variableTypeColumn ,] \
	    -col=ElementOccurence,type=long -data=[join $elementOccurenceColumn ,] \
	    -para=PageID,type=string -data=special
    } else {
	exec sddsmakedataset $variableFile.special \
	    -col=ElementName,type=string -data \
	    -col=ElementParameter,type=string -data \
	    -col=ParameterMode,type=string -data \
	    -col=ParameterValue,type=double -data \
	    -col=VariableType,type=string -data \
	    -col=ElementOccurence,type=long -data \
	    -para=PageID,type=string -data=special
    }

    file copy $variableFile.special $allElementsFile.special
    lappend varFileList $variableFile.special
    lappend elemFileList $allElementsFile.special

    #------ Including previous solution values into initial files...
    if $usePreviousSolution {
        foreach pageID [linsert $variableTypeList 100 special] {
            set tmpFile ${tmpRoot}-prevSol.$pageID
            exec sddsprocess $previousResultsFile $tmpFile -match=par,PageID=$pageID -nowarning 
	    lappend deleteFiles $tmpFile
            if [exec sdds2stream -rows=bare $tmpFile] {
                exec sddsxref $allElementsFile.$pageID $tmpFile -pipe=out -nowarning -take=ParameterValue \
                    -rename=col,ParameterValue=ParameterValueOld -match=TagName -fillIn \
                    | sddsprocess -pipe "-redef=col,ParameterValue,ParameterValue ParameterValueOld +" -nowarning \
                    | sddsconvert -pipe=in $tmpFile.1 -del=col,ParameterValueOld
                file copy -force $tmpFile.1 $allElementsFile.$pageID
                lappend deleteFiles $tmpFile.1
            }
        }
    }
    #------ Combining final files...
    eval exec sddscombine $varFileList $variableFile.sdds -overWrite
    eval exec sddscombine $varFileList $elementsUpdateFile -merge -overWrite
    eval exec sddscombine $elemFileList -pipe=out | sddsconvert -pipe=in $allElementsFile.sdds -del=col,ElementOccurence

    #--- Make real parameter file from elements.special file
    if $options(useSpecElem) {
	if [catch {MakeRealParamFileFromSpecialVariables -specialElementsInputFile $specialElementsInputFile -specialElementsFile $allElementsFile.special \
		       -specialParamFile $specialParamFile} result] {
	    return -code error "MakeRealParamFileFromSpecialVariables: $result"
	}
    }

    #------ Creating nonlinFiles from allElements bpm files ----
    file delete $nonlinFile.X $nonlinFile.Y
    if {$options(useBPMs) && $nonlinLevel > 0} {
        #------ This will give error if some bpms are not in fitting but in the response matrix.
        exec sddsconvert $allElementsFile.xBpm -pipe=out -retain=col,ElementName,TagName -rename=col,ElementName=BPMName \
            | sddsprocess -pipe=in $tmpRoot.nonlinX -def=col,pfitIntercept,0 -def=col,pfitSlope,0 -def=col,pfitCubic,0 \
	    -reprint=para,PageID,xBpmNonlin
        exec sddsconvert $allElementsFile.yBpm -pipe=out -retain=col,ElementName,TagName -rename=col,ElementName=BPMName \
            | sddsprocess -pipe=in $tmpRoot.nonlinY -def=col,pfitIntercept,0 -def=col,pfitSlope,0 -def=col,pfitCubic,0 \
	    -reprint=para,PageID,yBpmNonlin
        switch -exact -- $coupledMatrix {
            0 {
                file copy -force $tmpRoot.nonlinX $nonlinFile.X
                file copy -force $tmpRoot.nonlinY $nonlinFile.Y
            }
            1 {
                return -code error "Nonlin correction for coupledMatrix=1 is not done yet."
            }
            2 {
                exec sddscombine $tmpRoot.nonlinX $tmpRoot.nonlinY $nonlinFile.X -overWrite
                file copy -force $nonlinFile.X $nonlinFile.Y
            }
        }
        lappend deleteFiles $tmpRoot.nonlinX $tmpRoot.nonlinY
    }
    
    #------ Creating averageErrorFile --------------------------
    switch -exact -- $coupledMatrix {
        0 {
            set hResponseBPMs [lindex $elementListRM 2]
            set vResponseBPMs [lindex $elementListRM 3]
        }
        1 {
            set hResponseBPMs [lindex $elementListRM 3]
            set vResponseBPMs [lindex $elementListRM 2]
        }
        2 {
            set hResponseBPMs [concat [lindex $elementListRM 2] [lindex $elementListRM 3]]
            set vResponseBPMs [concat [lindex $elementListRM 2] [lindex $elementListRM 3]]
        }
    }
    foreach plane [list x y] bpmList [list $hResponseBPMs $vResponseBPMs] {
        if [llength $bpmList] {
            exec sddsmakedataset $averageErrorFile.$plane -col=TagName,type=string -data=[join $bpmList ,]
	    if [catch {Fit_MakeElementNameFromTagName -input $averageErrorFile.$plane} result] {
		return -code error "Fit_MakeElementNameFromTagName: $result"
	    }
        } else {
            exec sddsmakedataset $averageErrorFile.$plane -col=BPMName,type=string -data \
		-col=TagName,type=string -data -col=ElementOccurence,type=long -data
        }
        if [catch {ResortSingleFile -orderFile $elementOrderFile -orderColumn TagName \
                       -file2sort $averageErrorFile.$plane -column2sort TagName} result] {
            return -code error "ResortSingleFile (aver BPMs): $result"
        }
    }
    exec sddscombine $averageErrorFile.x $averageErrorFile.y -pipe=out \
        | sddsprocess -pipe=in $averageErrorFile -def=col,AverageError,0 -nowarning
    lappend deleteFiles $averageErrorFile.x $averageErrorFile.y
    
    #------ Creating dispersionFile ----------------------------------------------------
    if $options(useEnergy) {
	set latticeOptions "-lteFile $lteFile -beamlineName $beamlineName -latticeParamFileList \"$latticeParamFileList\""
        if [catch {eval CalculateResponseMatrix -getBetaFunctions 1 -workDir $workDir \
                       $latticeOptions -twiFile $twiFile \
		       -acceleratorCode $acceleratorCode } result] {
            return -code error "CalculateResponseMatrix: $result"
        }
        #------ Making files with bpms used in fitting
        set fileList ""
        foreach etaColumn [list etax etay] plane [list X Y] bpmFile [list $allElementsFile.xBpm $allElementsFile.yBpm] {
            exec sddsconvert $twiFile -pipe=out -retain=col,ElementName,ElementOccurence,$etaColumn -rename=col,$etaColumn=$dispColumn \
                | sddsprocess -pipe -print=para,Plane,$plane "-redef=col,$dispColumn,$dispColumn 1000 *,units=mm" \
		-print=col,TagName,%s#%1d#$plane,ElementName,ElementOccurence \
		| sddsselect -pipe=in $bpmFile $tmpRoot-energy.$plane -match=ElementName
            lappend fileList $tmpRoot-energy.$plane
        }
        eval exec sddscombine $fileList $dispersionFile -overWrite
        set deleteFiles [concat $deleteFiles $fileList]
    }
    eval file delete $deleteFiles
}

#-------------------------------------------------------------------------------------------------------------------------
# This procedures resort variableFile according to orderFile.

proc ResortSingleFile {args} {
    global tmpDir
    APSParseArguments {orderFile orderColumn file2sort column2sort}
    set tmpRoot $tmpDir/[APSTmpString]-resort
    set rowsBefore [exec sdds2stream $file2sort -rows=bare]
    if {[lsearch -exact [exec sddsquery -col $orderFile] Index] == -1 } {
        exec sddsprocess $orderFile $tmpRoot.index -def=col,Index,i_row
    } else {
        file copy -force $orderFile $tmpRoot.index
    }
    exec sddsxref $file2sort $tmpRoot.index -pipe=out -take=Index -match=$column2sort=$orderColumn -nowarning \
        | sddssort -pipe -col=Index \
        | sddsconvert -pipe=in $tmpRoot.sorted -del=col,Index
    set rowsAfter [exec sdds2stream $tmpRoot.sorted -rows=bare]
    if {$rowsBefore != $rowsAfter} {
        file copy $file2sort $tmpRoot.beforeSorting
	if [expr $rowsBefore-$rowsAfter < 10] {
	    set extraElements [join [exec sddsselect $tmpRoot.beforeSorting $tmpRoot.sorted -pipe=out -match=TagName -inv \
					 | sdds2stream -pipe=in -col=TagName] ]
	    OutputStatusMessage "Warning: ResortSingleFile: The following elements were removed from [file tail $file2sort] during sorting: $extraElements"
	} else {
	    OutputStatusMessage "Warning: ResortSingleFile: Many elements were removed from [file tail $file2sort] during sorting (see loco.log)."
	}
        OutputStatusMessage "Warning: ResortSingleFile: Number of rows in file $file2sort before the sorting \
                                       ($rowsBefore) is different from after the sorting ($rowsAfter). The original file \
                                       $file2sort was copied to $tmpRoot.beforeSorting." -logFileOnly 1
    }
    file copy -force $tmpRoot.sorted $file2sort
    file delete $tmpRoot.index $tmpRoot.sorted
}

#-------------------------------------------------------------------------------------------------------------------------
proc FilterMeasuredResponseFiles {args} {
    global tmpDir elementListRM options
    APSParseArguments {hrmMeasured vrmMeasured elementOrderFile useKickFiles kickFileX kickFileY coupledMatrix dispColumn}
    set tmpRoot $tmpDir/[APSTmpString]-filterRM
    set deleteFiles ""

    if $useKickFiles {
	set corrXList [join [exec sddscollapse $kickFileX -pipe=out | sdds2stream -pipe=in -col=ColumnName] ]
	set corrYList [join [exec sddscollapse $kickFileY -pipe=out | sdds2stream -pipe=in -col=ColumnName] ]
    } else {
	set corrXList [lindex $elementListRM 0]
	set corrYList [lindex $elementListRM 1]
    }
    if $options(fitDispersion) {set corrXList [concat $corrXList $dispColumn]}
    set bpmXList  [lindex $elementListRM 2]
    set bpmYList  [lindex $elementListRM 3]
    
    exec sddsmakedataset $tmpRoot.bpmsx -col=TagName,type=string -data=[join $bpmXList ,]
    exec sddsmakedataset $tmpRoot.bpmsy -col=TagName,type=string -data=[join $bpmYList ,]
    exec sddscombine $tmpRoot.bpmsx $tmpRoot.bpmsy $tmpRoot.bpms -overWrite
    lappend deleteFiles $tmpRoot.bpmsx $tmpRoot.bpmsy $tmpRoot.bpms
    exec sddsconvert $hrmMeasured -pipe=out -retain=column,TagName,BPMName,[join $corrXList ,] \
	| sddsselect -pipe=in $tmpRoot.bpms $hrmMeasured.filtered -match=TagName
    exec sddsconvert $vrmMeasured -pipe=out -retain=column,TagName,BPMName,[join $corrYList ,] \
	| sddsselect -pipe=in $tmpRoot.bpms $vrmMeasured.filtered -match=TagName
    lappend deleteFiles $hrmMeasured.filtered $vrmMeasured.filtered
    #------ Keeping required pages in calculated files:
    switch -regexp -- $coupledMatrix {
	0 {
	    exec sddsprocess $hrmMeasured.filtered $hrmMeasured -match=para,Plane=X
	    exec sddsprocess $vrmMeasured.filtered $vrmMeasured -match=para,Plane=Y
	}
	2 {
	    file copy -force $hrmMeasured.filtered $hrmMeasured
	    file copy -force $vrmMeasured.filtered $vrmMeasured
	}
    }
    eval file delete $deleteFiles
}
#-------------------------------------------------------------------------------------------------------------------------
proc AddDispersionColumn {args} {
    global tmpDir dispFitWeight
    APSParseArguments {hrmMeasured dispMeasuredFile dispColumn}
    set tmpRoot $tmpDir/[APSTmpString]-addD
    #------ See if any BPMs are missing from dispersion file:
    if [catch {exec sddsselect $hrmMeasured $dispMeasuredFile -pipe=out -match=TagName -inv \
		   | sddscombine -pipe -merge \
		   | sdds2stream -pipe=in -col=TagName} misingBpms] {
	return -code error "Error reading missing BPMs in $dispMeasuredFile file $misingBpms"
    }
    if [llength $misingBpms] {
	puts stdout "Error: the following BPMs are missing from dispersion file dispMeasuredFile: $misingBpms"
	exit
    }
    exec sddsxref $hrmMeasured $dispMeasuredFile -pipe=out -take=$dispColumn -match=TagName -nowarning \
	| sddsprocess -pipe=in $tmpRoot "-redef=col,$dispColumn,$dispColumn $dispFitWeight *,symbol=,units="
    file copy -force $tmpRoot $hrmMeasured
    file delete $tmpRoot
}
#-------------------------------------------------------------------------------------------------------------------------
proc SortMeasuredResponseFiles {args} {
    global tmpDir
    APSParseArguments {hrmMeasured vrmMeasured elementOrderFile coupledMatrix}
    set tmpRoot $tmpDir/[APSTmpString]-reorder
    set deleteFiles ""
    foreach rmFile [list $hrmMeasured $vrmMeasured] {
	if {$coupledMatrix == 2} {
	    exec sddsprocess $rmFile $tmpRoot.pagex -match=para,Plane=X
	    exec sddsprocess $rmFile $tmpRoot.pagey -match=para,Plane=Y
	    set singlePageFileList [list $tmpRoot.pagex $tmpRoot.pagey]
	    lappend deleteFiles $tmpRoot.pagex $tmpRoot.pagey
	} else {
	    set singlePageFileList $rmFile
	}
	foreach singlePageFile $singlePageFileList {
	    if [catch {ResortSingleFile -orderFile $elementOrderFile -orderColumn TagName \
			   -file2sort $singlePageFile -column2sort TagName} result] {
		return -code error "ResortSingleFile (file2sort: $file2sort): $result"
	    }
	    exec sddstranspose $singlePageFile $singlePageFile.cols -newcolumn=TagName -oldcolumn=Corrector
	    if [catch {ResortSingleFile -orderFile $elementOrderFile -orderColumn TagName \
			   -file2sort $singlePageFile.cols -column2sort Corrector} result] {
		return -code error "ResortSingleFile (file2sort: $singlePageFile.cols): $result"
	    }
	    exec sddstranspose $singlePageFile.cols $singlePageFile.rows -newcolumn=Corrector -oldcolumn=TagName
	    exec sddsconvert $singlePageFile -pipe=out -retain=col,*Name \
		| sddsxref -pipe=in $singlePageFile.rows $singlePageFile.result -take=* -leave=*Name
	    file copy -force $singlePageFile.result $singlePageFile
	    lappend deleteFiles $singlePageFile.cols $singlePageFile.rows $singlePageFile.result
	}
	if {$coupledMatrix == 2} {
	    eval exec sddscombine $singlePageFileList $rmFile -overWrite
	}
    }
    eval file delete $deleteFiles
}
#-------------------------------------------------------------------------------------------------------------------------
proc DetectZeroBPMs {args} {
    global tmpDir
    APSParseArguments {hrmMeasured vrmMeasured zeroBpmTol dispColumn}
    set tmpRoot $tmpDir/[APSTmpString]-detectZeroBpms
    foreach measuredFile [list $hrmMeasured $vrmMeasured] plane [list H V] {
        if {[string compare "H" $plane] == 0} {set page 1} else {set page 2}
        exec sddsconvert $measuredFile $tmpRoot -keepPages=$page
	if [catch {Fit_GetZeroBpms -filename $tmpRoot -threshold $zeroBpmTol -column TagName -dispColumn $dispColumn} zeroBpms] {
	    return -code error "Fit_GetZeroBpms (measuredFile: $measuredFile): $zeroBpms"
	}
	if [llength $zeroBpms] {
	    if {[llength $zeroBpms] < 10} {
		OutputStatusMessage "DetectZeroBPMs: Warning: the following bpms in $plane \
                have zero readings only: [join $zeroBpms ,]"
	    } else {
		OutputStatusMessage "DetectZeroBPMs: Warning: the following bpms in $plane \
                have zero readings only: [join $zeroBpms ,]" -logFileOnly 1
		OutputStatusMessage "DetectZeroBPMs: Warning: some bpms in $plane have zero readings only (see log file)"
	    }
	}
    }
    file delete $tmpRoot
}
#-------------------------------------------------------------------------------------------------------------------------

proc MakeWeightFile { args } {
    global tmpDir useBpmWeight
    APSParseArguments { useBpmWeight weightFile manualWeightFile hrmFile vrmFile weightCorrectorsFile coupledMatrix }
    if $useBpmWeight {
        if ![file exists $manualWeightFile] {
            return -code error "useBpmWeight is on, but the weight file $manualWeightFile is missing."
        }
        set tmpRoot $tmpDir/[APSTmpString]-weight
        switch -exact -- $coupledMatrix {
            0 {
                set rmFileList [list $hrmFile $vrmFile]
                set planeList [list X Y]
            }
            1 {
                set rmFileList [list $vrmFile $hrmFile]
                set planeList [list Y X]
            }
            2 {
                exec sddsconvert $hrmFile $tmpRoot.page1 -keepPage=1
                exec sddsconvert $hrmFile $tmpRoot.page2 -keepPage=2
                set rmFileList [list $tmpRoot.page1 $tmpRoot.page2]
                set planeList [list X Y]
            }
        }
        foreach rmFile $rmFileList plane $planeList {
            exec sddsprocess $manualWeightFile $tmpRoot.manual -nowarning -match=col,Plane=$plane
            if [exec sdds2stream $tmpRoot.manual -rows=bare] {
                exec sddsconvert $rmFile -pipe=out -retain=col,TagName \
                    | sddsprocess -pipe=in $tmpRoot.$plane -redef=col,Index,i_row -del=para,*
                set missingBpms ""
                if [catch {exec sddsselect $tmpRoot.manual $tmpRoot.$plane -pipe=out -match=TagName -invert \
                               | sdds2stream -pipe=in -col=TagName} missingBpms] {
                    return -code error "$missingBpms"
                }
                if [llength $missingBpms] {
                    OutputStatusMessage "makeWeightFile: Warning processing $plane: The following bpms in \
                         the weight file: [join $missingBpms ] are missing from the $rmFile file." -logFileOnly 1
                }
                exec sddsselect $tmpRoot.$plane $tmpRoot.manual -pipe=out -match=TagName \
                    | sddsxref -pipe=in $tmpRoot.manual $tmpRoot.$plane.weight -take=Weight -match=TagName
                exec sddsselect $tmpRoot.$plane $tmpRoot.manual -pipe=out -match=TagName -invert \
                    | sddsprocess -pipe=in $tmpRoot.$plane.noweight -def=col,Weight,1.0
                exec sddscombine $tmpRoot.$plane.noweight $tmpRoot.$plane.weight -pipe=out -merge \
                    | sddssort -pipe -col=Index,increasing \
                    | sddsprocess -pipe -reprint=para,Plane,$plane \
                    | sddsconvert -pipe=in $tmpRoot.$plane -del=col,Index
                file delete $tmpRoot.$plane.weight $tmpRoot.$plane.noweight
            } else {
                exec sddsconvert $rmFile -pipe=out -retain=col,TagName \
                    | sddsprocess -pipe=in $tmpRoot.$plane -def=col,Weight,1.0 -reprint=para,Plane,$plane
            }
            file delete $tmpRoot.manual
        }
        exec sddscombine -overWrite $tmpRoot.X $tmpRoot.Y $weightFile
        file delete $tmpRoot.X $tmpRoot.Y $tmpRoot.page1 $tmpRoot.page2
    }
    if [file exists $weightCorrectorsFile] {
        OutputStatusMessage "makeWeightFile: Warning, corrector weight are not working now!!!"
    } 
}

#-------------------------------------------------------------------------------------------------------------------------

proc RemoveAllSpaces {line} {
    set strLength [string length $line]
    set line1 ""
    for {set I 0} {$I < $strLength} {incr I} {
	if {[string compare " " [string range $line $I $I]] != 0} {
	    append line1 [string range $line $I $I]
	}
    }
    return $line1
}

#-------------------------------------------------------------------------------------------------------------------------

proc MakeSigmaWeightFiles {args} {
    global tmpDir dispColumn options sigmaWeightFileX sigmaWeightFileY
    APSParseArguments {dispMeasuredFile hrmMeasured vrmMeasured}
    set tmpRoot $tmpDir/[APSTmpString]-sigmaWeight
    set deleteFiles ""
    #------ This part makes matrix that has 1.0 everywhere except in the dispColumn, where it has weights
    if $options(fitDispersion) {
	#------ dispMeasuredFile has DispersionSigma column.
	if [catch {exec sddsoutlier $dispMeasuredFile -pipe=out -col=DispersionSigma -stdev=3 \
		       | sddsoutlier -pipe -col=DispersionSigma -stdev=3 \
		       | sddsprocess -pipe=in $tmpRoot.1 -process=DispersionSigma,average,SigmaAver} result] {
	    return -code error "Error getting average sigma: $result"
	}
	set sigmaWeightCutoff 5
	exec sddsxref $dispMeasuredFile $tmpRoot.1 -pipe=out -leave=* -transfer=para,SigmaAver \
	    | sddsprocess -pipe=in $tmpRoot.2 "-redef=col,SigmaWeight,DispersionSigma SigmaAver /" \
	    "-redef=col,SigmaWeight,SigmaWeight $sigmaWeightCutoff < ? 1.0 : SigmaWeight \$"
	lappend deleteFiles $tmpRoot.1 $tmpRoot.2
	exec sddsconvert $hrmMeasured -pipe=out -del=col,$dispColumn \
	    | sddsprocess -pipe -redef=col,%s,1.0,select=*,exclude=*Name \
	    | sddsxref -pipe=in $tmpRoot.2 $sigmaWeightFileX -take=SigmaWeight -match=TagName \
	    -rename=col,SigmaWeight=$dispColumn -fillin -nowarning
	exec sddsprocess $vrmMeasured $sigmaWeightFileY -redef=col,%s,1.0,select=*,exclude=*Name
    } else {
	exec sddsprocess $hrmMeasured $sigmaWeightFileX -redef=col,%s,1.0,select=*,exclude=*Name 
	exec sddsprocess $vrmMeasured $sigmaWeightFileY -redef=col,%s,1.0,select=*,exclude=*Name
    }
    #------ Here do processing of RM sigmas if they ever become available...
    catch {eval file delete $deleteFiles}
}

#-------------------------------------------------------------------------------------------------------------------------

proc GetTolerance {args} {
    global bigIterNumber
    APSParseArguments {bigIter}
    set varList [list convergenceTolerance convergenceTolerance0 convergenceTolerance1 convergenceTolerance2]
    set tolList ""
    foreach var $varList {
	global $var
	if [info exists $var] {lappend tolList [set $var]}
    }
    set max 0
    foreach tol $tolList {if {$tol > $max} {set max $tol}}
    set min 1e10
    foreach tol $tolList {if {$tol < $min} {set min $tol}}
    #--- Set tol to max only if it is the first big iteration and bigIterNumber is larger than 1
    if {$bigIterNumber > 1 && $bigIter == 0} {set tolerance $max} else { set tolerance $min}
    return $tolerance
}

#-------------------------------------------------------------------------------------------------------------------------
#-------------------------------------------------------------------------------------------------------------------------
#------ Procedures related to calculation of response matrices during iterations --------------------------------------------
#-------------------------------------------------------------------------------------------------------------------------
#-------------------------------------------------------------------------------------------------------------------------

#-------------------------------------------------------------------------------------------------------------------------
# Corrects measured matrix by removing bad points and average non-zero error (proc CorrectMeasuredResponse),
# then calculates response matrices (proc CalculateMatrices). 
# Then calculates difference with the measured ones, then applies BPM weight correction, then makes one column.

proc CalculateRMDifference { args } {
    
    set residualError 10.0
    set badPointsLevel 0
    set nonlinLevel 0
    set averLevel 0
    set averageErrorFile ""
    set coupledMatrix 0
    APSParseArguments { xyColumnFile hrmMeasured hrmIterFile hrmDifference acceleratorCode \
                                  vrmMeasured vrmDifference vrmIterFile latticeFile \
                                  dispersionFile allElementsFile averageErrorFile coupledMatrix \
                                  residualError badPointsLevel nonlinLevel averLevel currentIteration \
				  useQuadConstraints quadConstraintWeight useKickFiles kickFileX kickFileY}
    
    global elementListRM tmpDir options keepTunesOnTarget
    
    set hYes [llength [lindex $elementListRM 0]]
    set vYes [llength [lindex $elementListRM 1]]
    
    #------ Calculate updated response matrix and make nonlinear correction.
    set time1 [exec date +%s]
    if [catch {CalculateMatrices -residualError $residualError \
                   -hrmFile $hrmIterFile -vrmFile $vrmIterFile -nonlinLevel $nonlinLevel \
                   -hrmMeasured $hrmMeasured -vrmMeasured $vrmMeasured \
                   -dispersionFile $dispersionFile -allElementsFile $allElementsFile \
                   -coupledMatrix $coupledMatrix -acceleratorCode $acceleratorCode -currentIteration $currentIteration \
		   -useKickFiles $useKickFiles -kickFileX $kickFileX -kickFileY $kickFileY} result] {
        return -code error "CalculateMatrices: $result"
    }
    set time2 [exec date +%s]
    OutputStatusMessage "CalculateMatrices took [expr $time2 - $time1] sec." -logFileOnly 1
    
    #------ Correct measured response matrix: Removing bad points and subtracting average.
    set time1 [exec date +%s]
    if [catch {CorrectMeasuredResponse -hYes $hYes -vYes $vYes -hrmIterFile $hrmIterFile -vrmIterFile $vrmIterFile \
                   -hrmMeasured $hrmMeasured -vrmMeasured $vrmMeasured -badPointsLevel $badPointsLevel \
                   -averLevel $averLevel -averageErrorFile $averageErrorFile \
                   -residualError $residualError} result] {
        return -code error "CorrectMeasuredResponse: $result"
    }
    set time2 [exec date +%s]
    OutputStatusMessage "CorrectMeasuredResponse took [expr $time2 - $time1] sec." -logFileOnly 1
    
    #------ Calculate difference between measured and calculated matrices and combine into a single column.
    set time1 [exec date +%s]
    if [catch {SubtractMatricesAndMakeOneColumn -hYes $hYes -vYes $vYes -coupledMatrix $coupledMatrix \
                   -hrmMeasured $hrmMeasured -vrmMeasured $vrmMeasured \
                   -hrmIterFile $hrmIterFile -vrmIterFile $vrmIterFile -hrmDifference $hrmDifference \
                   -vrmDifference $vrmDifference -xyColumnFile $xyColumnFile \
		   -useQuadConstraints $useQuadConstraints -quadConstraintWeight $quadConstraintWeight} rmsList] {
        return -code error "SubtractMatricesAndMakeOneColumn: $rmsList"
    }
    set time2 [exec date +%s]
    OutputStatusMessage "SubtractMatricesAndMakeOneColumn took [expr $time2 - $time1] sec." -logFileOnly 1
    exec sddsconvert $xyColumnFile -nowarning -retain=para,ConstraintLines
    file delete ${xyColumnFile}~

    #------ Calculate in-iteration tune adjustment to keep them on target:
    if $keepTunesOnTarget {
	if [catch {KeepTunesOnTarget} result] {
	    return -code error "KeepTunesOnTarget: $result"
	}
    }
    return $rmsList
}

#-------------------------------------------------------------------------------------------------------------------------
proc KeepTunesOnTarget {args} {
    global tmpDir nuxMeasured nuyMeasured twiFile nuxUpFile nuyUpFile tuneCorrectionFile
    set corrFraction 0.75
    set tuneTol 0.01
    set tmpRoot $tmpDir/[APSTmpString]-KeepTunesOnTarget
    set calculatedTunes [exec sdds2stream -para=nux,nuy $twiFile]
    if {[lsearch -glob $calculatedTunes *NaN*] != -1 && [lsearch -glob $calculatedTunes *nan*] != -1} {
        return -code error "Error: Tunes are unstable!"
    }
    set nuxDelta0 [exec sdds2stream $nuxUpFile -para=DeltaNu]
    set nuyDelta0 [exec sdds2stream $nuyUpFile -para=DeltaNu]
    set nuxDelta [expr [lindex $calculatedTunes 0] - $nuxMeasured]
    set nuyDelta [expr [lindex $calculatedTunes 1] - $nuyMeasured]
    if [expr abs($nuxDelta) > $tuneTol || abs($nuyDelta) > $tuneTol] {set applyCorr 1} else {set applyCorr 0}
    exec sddsprocess $nuxUpFile $tmpRoot.nux "-redef=col,ParameterValue,ParameterValue $nuxDelta0 / $nuxDelta * $corrFraction * -1 *"
    exec sddsprocess $nuyUpFile $tmpRoot.nuy "-redef=col,ParameterValue,ParameterValue $nuyDelta0 / $nuyDelta * $corrFraction * -1 *"
    exec sddscombine $tmpRoot.nux $tmpRoot.nuy -pipe=out -merge \
	| sddssort -pipe -col=TagName \
	| sddsbreak -pipe -change=TagName \
	| sddsprocess -pipe -process=ParameterValue,sum,Sum -clip=1,0,inv -redef=col,ParameterValue,Sum \
	| sddscombine -pipe -merge \
	| sddsxref -pipe $twiFile -take=s -match=TagName -nowarning \
	| sddssort -pipe -col=s \
	| sddsprocess -pipe -redef=para,dQx,$nuxDelta -redef=para,dQy,$nuyDelta -redef=para,ApplyCorrection,$applyCorr,type=long \
	| sddsconvert -pipe=in $tuneCorrectionFile -del=col,s -del=para,Sum,DeltaNu
    file delete $tmpRoot.nux $tmpRoot.nuy
}
#-------------------------------------------------------------------------------------------------------------------------
# This procedure removes bad points and subtracts averages.

proc CorrectMeasuredResponse {args} {
    
    APSParseArguments { hYes vYes hrmIterFile vrmIterFile hrmMeasured vrmMeasured badPointsLevel \
                            averLevel averageErrorFile residualError } 
    
    global badPointsIterCounter
    if ![info exists badPointsIterCounter] {set badPointsIterCounter 0}
    
    #------ Removing bad points ------------------------------
    if {$residualError < $badPointsLevel} {
        incr badPointsIterCounter
        switch -exact -- $badPointsIterCounter {
            1       {set NSigma 5}
            2       {set NSigma 4}
            default {set NSigma 3}
        }
        
        set xPoints ""
        set yPoints ""
        if $hYes {
            if [catch {RemoveBadPoints -computedResponseFile $hrmIterFile \
                           -measuredResponseFile $hrmMeasured -NSigma $NSigma \
                           -outputFile $hrmMeasured} xPoints] {
                return -code error "RemoveBadPoints: $xPoints"
            }
        }
        if $vYes {
            if [catch {RemoveBadPoints -computedResponseFile $vrmIterFile \
                           -measuredResponseFile $vrmMeasured -NSigma $NSigma \
                           -outputFile $vrmMeasured} yPoints] {
                return -code error "RemoveBadPoints: $yPoints"
            }
        }
	set pointsRemoved [expr [lindex $xPoints 0] + [lindex $yPoints 0]]
	set pointsRemovedTotal [expr [lindex $xPoints 1] + [lindex $yPoints 1]]
        OutputStatusMessage \
	    "RemoveBadPoints ($NSigma sigmas): [format %.0f $pointsRemoved] bad points removed (total: [format %.0f $pointsRemovedTotal])."
    }
    
    #------ Subtracting average -------------------------------
    if {$residualError < $averLevel} {
        if $hYes {
            if [catch {SubtractAverageError -computedResponseFile $hrmIterFile \
                           -measuredResponseFile $hrmMeasured -averageErrorFile $averageErrorFile \
                           -plane X} result] {
                return -code error "SubtractAverageError: $result"
            }
        }
        if $vYes {
            if [catch {SubtractAverageError -computedResponseFile $vrmIterFile \
                           -measuredResponseFile $vrmMeasured -averageErrorFile $averageErrorFile \
                           -plane Y} result] {
                return -code error "SubtractAverageError: $result"
            }
        }
    }
}

#-------------------------------------------------------------------------------------------------------------------------
#------ This procedure changes hrmMeasured and vrmMeasured files

proc RemoveBadPoints { args } {

    global tmpDir workDir dispColumn
    
    APSParseArguments { measuredResponseFile computedResponseFile outputFile NSigma }
    
    set tmpRoot $tmpDir/[APSTmpString]-removeBadPoints
    set deleteFiles ""
    #------ $dispColumn is removed if exists and merge if two pages...
    exec sddsconvert $measuredResponseFile $tmpRoot.measRM -del=col,$dispColumn
    exec sddsconvert $computedResponseFile $tmpRoot.compRM -del=col,$dispColumn
    lappend deleteFiles $tmpRoot.measRM $tmpRoot.compRM
    set twoPages 0
    if {[exec sdds2stream $computedResponseFile -npages=bare] == 2} {
        set twoPages 1
        exec sddscombine $tmpRoot.measRM $tmpRoot.measRM.tmp -merge
        exec sddscombine $tmpRoot.compRM $tmpRoot.compRM.tmp -merge
        file copy -force $tmpRoot.measRM.tmp $tmpRoot.measRM
        file copy -force $tmpRoot.compRM.tmp $tmpRoot.compRM
        file delete $tmpRoot.measRM.tmp $tmpRoot.compRM.tmp
        exec sddsprocess $computedResponseFile -pipe=out "-def=col,PageNumber,i_page" \
            | sddscombine -pipe -merge \
            | sddsconvert -pipe=in $tmpRoot.pageNumber -retain=col,PageNumber
    }
    
    #------ Create a mask file...
    set maskFile $tmpRoot.mask
    exec sddschanges $tmpRoot.measRM -pipe=out -base=$tmpRoot.compRM -copy=TagName -change=exclude=*Name,* \
        | sddsconvert -pipe "-edit=col,ChangeIn*,a 8d" \
        | sddstranspose -pipe -newColumn=TagName -oldColumn=Corrector \
        | tee $maskFile.transp \
        | sddsconvert -pipe -del=col,Corrector \
        | sddsprocess -pipe "-proc=*,rms,%sRms" \
        | sddsprocess -pipe "-redef=col,%s,%s abs %sRms $NSigma * > ? 1 : 0 \$,select=*" \
        | sddsxref -pipe $maskFile.transp -take=Corrector \
        | sddstranspose -pipe=in $maskFile -newColumn=Corrector -oldColumn=TagName
    
    #------ Get the number of bad points and the number of previously removed points:
    set nPoints [exec sddsrowstats $maskFile -pipe=out -sum=ColSum,* \
                     | sddsprocess -pipe -proc=ColSum,sum,Sum \
                     | sdds2stream -pipe=in -para=Sum] 
    if [catch {exec sddscombine $measuredResponseFile -pipe=out -merge \
		   | sdds2stream -pipe=in -para=BadPointsRemoved} nPointsTotal] {
	set nPointsTotal 0
    }
    set nPointsTotal [expr $nPointsTotal + $nPoints]

    #------ Substitute computed elements using the mask and redefining BadPointsRemoved parameter:
    file copy -force $tmpRoot.measRM $tmpRoot.measRM.copy
    eval exec sddsmatrixop $tmpRoot.measRM.copy -pipe=out -push=$tmpRoot.measRM -push=$tmpRoot.compRM \
        -subtract -push=$maskFile -multiply=hadamard -subtract -columnNames=file=$maskFile.transp,column=Corrector \
        | sddsxref -pipe $tmpRoot.measRM.copy -take=TagName \
	| sddsprocess -pipe=in $tmpRoot.measRM -redef=para,BadPointsRemoved,$nPointsTotal
    
    #------ Split back in 2 pages...
    if $twoPages {
        exec sddsxref $tmpRoot.measRM $tmpRoot.pageNumber -pipe=out -take=PageNumber \
            | sddsbreak -pipe -change=PageNumber \
	    | sddsxref -pipe $measuredResponseFile -leave=* -transfer=para,Plane \
            | sddsconvert -pipe=in $tmpRoot.measRM.tmp2 -del=col,PageNumber
        file copy -force $tmpRoot.measRM.tmp2 $tmpRoot.measRM
        lappend deleteFiles $tmpRoot.measRM.tmp2 $tmpRoot.pageNumber
    }
    
    #------ Reorder columns to put string columns first (for aesthetic purposes):
    exec sddsconvert $tmpRoot.measRM -pipe=out -retain=col,*Name \
	| sddsxref -pipe=in $tmpRoot.measRM $tmpRoot.measRM.reorder -take=* -leave=*Name
    lappend deleteFiles $tmpRoot.measRM.reorder

    #------ Return dispColumn back:
    if {[lsearch -exact [exec sddsquery -col $measuredResponseFile] $dispColumn] != -1 } {
        exec sddsxref $tmpRoot.measRM.reorder $measuredResponseFile $tmpRoot.measRM.1 -take=$dispColumn
        file copy -force $tmpRoot.measRM.1 $outputFile
        lappend deleteFiles $tmpRoot.measRM.1
    } else {
        file copy -force $tmpRoot.measRM.reorder $outputFile
    }
    eval file delete $deleteFiles $maskFile $tmpRoot.measRM.copy $maskFile.transp

    return [list $nPoints $nPointsTotal]
}

#-------------------------------------------------------------------------------------------------------------------------
# Subtracts nonzero average from each bpm. This is supposed to take care of dispersion.

proc SubtractAverageError {args} {
    
    global tmpDir dispColumn
    
    set tmpFile $tmpDir/[APSTmpString].avError
    APSParseArguments {computedResponseFile measuredResponseFile averageErrorFile plane}
    if [catch {exec sddschanges $computedResponseFile -pipe=out -base=$measuredResponseFile \
                   -change=exclude=*Name,* -copy=*Name -parallelPages \
                   | sddsconvert -pipe -del=col,ChangeIn$dispColumn \
                   | sddsrowstats -pipe -mean=AverageError,ChangeIn* \
                   | sddsconvert -pipe -retain=col,*Name -retain=col,AverageError \
                   | tee $tmpFile \
                   | sddscombine -pipe=in $tmpFile.merged -merge -overWrite} result] {
        return -code error "Error calculating average errors: $result"
    }
    if [catch {exec sddsxref $measuredResponseFile $tmpFile -pipe=out -take=AverageError \
                   | sddsconvert -pipe -del=col,$dispColumn \
                   | sddsprocess -pipe "-redef=col,%s,%s AverageError +,select=*,exclude=*Name" \
                   | sddsconvert -pipe=in  $tmpFile.1 -del=col,AverageError} result] {
        return -code error "Error subtracting average: $result"
    }
    if {[lsearch -exact [exec sddsquery -col $measuredResponseFile] $dispColumn] != -1 } {
        exec sddsxref $tmpFile.1 $measuredResponseFile $tmpFile.2 -take=$dispColumn
        file copy -force $tmpFile.2 $measuredResponseFile
        file delete $tmpFile.2
    } else {
        file copy -force $tmpFile.1 $measuredResponseFile
    }
    
    exec sddsconvert $averageErrorFile $tmpFile.x -keepPages=1
    exec sddsconvert $averageErrorFile $tmpFile.y -keepPages=2
    switch -exact -- $plane {
        X {set inputFile $tmpFile.x}
        Y {set inputFile $tmpFile.y}
    }
    if [catch {exec sddsxref $inputFile $tmpFile.merged -pipe=out -take=AverageError -match=TagName \
                   -rename=col,AverageError=AverageError1 \
                   | sddsprocess -pipe "-redef=col,AverageError,AverageError AverageError1 +" \
                   | sddsconvert -pipe=in  $tmpFile.2 -del=col,AverageError1} result] {
        return -code error "Error updating average error file: $result"
    }
    switch -exact -- $plane {
        X {exec sddscombine $tmpFile.2 $tmpFile.y $averageErrorFile -overWrite}
        Y {exec sddscombine $tmpFile.x $tmpFile.2 $averageErrorFile -overWrite}
    }
    file delete $tmpFile $tmpFile.1 $tmpFile.2 $tmpFile.x $tmpFile.y $tmpFile.merged
}

#-------------------------------------------------------------------------------------------------------------------------
# This procedure calculates response matrices:
#    1. Runs acceleratorCode to calculate Response Matrix
#    2. Removes unused BPMs from the RM
#    3. Applies calibration correction to correctors (corrector tilts are in calculated matrix)
#    4. Applies gain and tilt correction to BPMs
#    5. Adds energy change to horizontal RM
#    6. Applies nonlinearity correction for BPMs

proc CalculateMatrices { args } {
    
    APSParseArguments { hrmFile vrmFile residualError nonlinLevel acceleratorCode \
                                  dispersionFile hrmMeasured vrmMeasured allElementsFile \
                                  coupledMatrix useKickFiles kickFileX kickFileY currentIteration}
    
    global workDir tmpDir options twiFile dispColumn nonlinFile directRMCalc tiltCorFraction
    global solutionParamFileList verbose

    if [catch {Fit_DeleteElementsFromList -elementList [join [exec sddsquery -col $hrmMeasured] ] \
                   -deleteElements "BPMName TagName $dispColumn"} hCorrList] {
        return -code error "Fit_DeleteElementsFromList: $hCorrList"
    }
    if [catch {Fit_DeleteElementsFromList -elementList [join [exec sddsquery -col $vrmMeasured] ] \
                   -deleteElements "BPMName TagName $dispColumn"} vCorrList] {
        return -code error "Fit_DeleteElementsFromList: $vCorrList"
    }
    set time1 [exec date +%s]
    if [catch {CalculateResponseMatrix \
		   -acceleratorCode $acceleratorCode -directCalc $directRMCalc \
		   -latticeParamFileList $solutionParamFileList \
                   -hrmFile $hrmFile -vrmFile $vrmFile -hCorrList $hCorrList -vCorrList $vCorrList \
                   -workDir $workDir -twiFile $twiFile -coupledMatrix $coupledMatrix \
                   -fitDispersion $options(fitDispersion) -useKickFiles $useKickFiles \
		   -kickFileX $kickFileX -kickFileY $kickFileY} result] {
        return -code error "CalculateResponseMatrix: $result"
    }
    set time2 [exec date +%s]
    OutputStatusMessage "calculateResponseMatrix took [expr $time2 - $time1] sec." -logFileOnly 1
    set hYes [llength $hCorrList]
    set vYes [llength $vCorrList]
    set tmpRoot $tmpDir/[APSTmpString]-calcMatr
    set planeFlagList [list $hYes $vYes]
    set computedFileList [list $hrmFile $vrmFile]
    set measuredFileList [list $hrmMeasured $vrmMeasured]

    if [catch {Filter2planeRM -hrmFile $hrmFile -vrmFile $vrmFile -hrmMeasured $hrmMeasured -vrmMeasured $vrmMeasured} result] {
	return -code error "Filter2planeRM: $result"
    }

    #------ Filtering out not used bpms by comparing with measuredFile...
    #------ This is where we lose double-plane readings of single-plane BPMs
    foreach planeFlag $planeFlagList computedFile $computedFileList measuredFile $measuredFileList {
        if $planeFlag {
	    exec sddsconvert $computedFile -pipe=out -delete=column,s \
		| sddsselect -pipe=in $measuredFile $tmpRoot -match=TagName
	    if $verbose {
		if {[exec sddscombine $tmpRoot -pipe=out -merge | sdds2stream -pipe=in -rows=bare] == 0} {
		    OutputStatusMessage "Error: zero rows after filtering (computedFile: $computedFile; measuredFile: $measuredFile)" \
			-logFileOnly 1
		    return -code error "Error: zero rows after filtering."
		}
	    }
	    file copy -force $tmpRoot $computedFile
	}
    }
    file delete $tmpRoot

    set time1 [exec date +%s]
    #------ Corrector calibration is now in parameter files and is included automatically in the RM calculation
    #------ Here we do dispersion only:
    if $options(useCorrs) {
	if [catch {ApplyCorrectorCalibrationCorrection -hrmFile $hrmFile -vrmFile $vrmFile \
		       -hYes $hYes -vYes $vYes -allElementsFile $allElementsFile -dispersionOnly 1} result] {
	    return -code error "ApplyCorrectorCalibrationCorrection: $result"
	}
    }
    set time2 [exec date +%s]
    OutputStatusMessage "Corrector calibration took [expr $time2 - $time1] sec." -logFileOnly 1

    if $options(useEnergy) {
        if [catch {AddEnergy2ResponseMatrix -responseFile $hrmFile -allElementsFile $allElementsFile \
                       -dispersionFile $dispersionFile} result] {
            return -code error "AddEnergy2ResponseMatrix: $result"
        }
    }
    #------ Nonlin correction is not applied to cross-response
    if {$residualError < $nonlinLevel} {
	set time1 [exec date +%s]
        set nonlinFileList [list $nonlinFile.X $nonlinFile.Y]
	set planeFlagList [list $hYes $vYes]
	set computedFileList [list $hrmFile $vrmFile]
	set measuredFileList [list $hrmMeasured $vrmMeasured]
	set pageList [list 1 2]
        foreach planeFlag $planeFlagList computedFile $computedFileList measuredFile $measuredFileList \
            nonlFile $nonlinFileList usePage $pageList {
		if [file exists $computedFile] { 
		    if [catch {ApplyNonlinCorrection -computedFile $computedFile -measuredFile $measuredFile \
				   -nonlinFile $nonlFile -coupledMatrix $coupledMatrix -usePage $usePage} result] {
			return -code error "ApplyNonlinCorrection ($nonlFile): $result"
		    }
		}
            }
	set time2 [exec date +%s]
	OutputStatusMessage "Nonlin correction took [expr $time2 - $time1] sec." -logFileOnly 1
    }
    set nuList [exec sdds2stream $twiFile -para=nux,nuy]
#    exec sddsprocess $hrmFile -nowarning -def=para,nux,[lindex $nuList 0] -def=para,nuy,[lindex $nuList 1]
#    exec sddsprocess $vrmFile -nowarning -def=para,nux,[lindex $nuList 0] -def=para,nuy,[lindex $nuList 1]
    file delete ${hrmFile}~ ${vrmFile}~
}

#-------------------------------------------------------------------------------------------------------------------------
# 2-plane RM is needed for BPM tilt/coef calculations. Regardless of possible HMON/VMON, it has all BPMs in both planes.
# Page 1 has #X suffix only and page 2 has #Y suffix only
proc Filter2planeRM {args} {
    global options allElementsFile tmpDir
    APSParseArguments {hrmFile vrmFile hrmMeasured vrmMeasured}
    set tmpRoot $tmpDir/[APSTmpString]-2plane
    if {$options(useBPMsTilt) || $options(useBPMsCoef)} {
	#--- Get list of all BPMs (mentioned in X or Y planes):
	exec sddsprocess $hrmMeasured -pipe=out "-edit=col,TagNameNoPlane,TagName,%@#X@@ %@#Y@@" \
	    | sddscombine -pipe -merge \
	    | sddssort -pipe -col=TagNameNoPlane \
	    | sddsbreak -pipe -change=TagNameNoPlane \
	    | sddsprocess -pipe -clip=1,0,inv \
	    | sddscombine -pipe=in $tmpRoot.bmps -merge -overwrite
	foreach rmFile [list $hrmFile $vrmFile] {
	    exec sddsconvert $rmFile -pipe=out -keep=1 \
		| sddsprocess -pipe -match=col,TagName=*#X* "-edit=col,TagNameNoPlane,TagName,%@#X@@ %@#Y@@" \
		| sddsselect -pipe=in $tmpRoot.bmps $tmpRoot.x -reuse=page -match=TagNameNoPlane
	    exec sddsconvert $rmFile -pipe=out -keep=2 \
		| sddsprocess -pipe -match=col,TagName=*#Y* "-edit=col,TagNameNoPlane,TagName,%@#X@@ %@#Y@@" \
		| sddsselect -pipe=in $tmpRoot.bmps $tmpRoot.y -reuse=page -match=TagNameNoPlane
	    exec sddscombine $tmpRoot.x $tmpRoot.y -pipe=out \
		| sddsconvert -pipe=in [file root $rmFile].2plane[file ext $rmFile] -del=col,TagNameNoPlane
	}
	file delete $tmpRoot.bmps $tmpRoot.x $tmpRoot.y
    }
}
#-------------------------------------------------------------------------------------------------------------------------

proc CalculateResponseMatrix { args } {    
    global locoBinDir verbose bRho elemByElemSR allElementsFile

    set getSDDSLattice 0
    set getBetaFunctions 0
    set lteFile ""
    set beamlineName ""
    set latticeParamFileList ""
    set sddsLatticeFile ""
    set twiFile ""
    set workDir ""
    set hrmFile ""
    set vrmFile ""
    set hCorrList ""
    set vCorrList ""
    set useQsub ""
    set dispFitWeight ""
    set directCalc ""
    set kickX ""
    set kickY ""
    set fixedOrbitLength ""
    set fitDispersion ""
    set dispColumn ""
    set coupledMatrix ""
    set closedOrbit ""
    set useDoubleOrbits ""
    set acceleratorCode ""
    set numberOfQsubTasks ""
    set waitTime ""
    set waitInterval ""
    set qsubCommand ""
    set queueSystemName ""
    set qsubRespProcCommand ""
    set submissionPause ""
    set twissRefFile ""
    set twissRefElement ""
    set malignElement ""
    set useKickFiles 0
    set kickFileX ""
    set kickFileY ""

    APSParseArguments { getSDDSLattice getBetaFunctions lteFile beamlineName latticeParamFileList \
			    sddsLatticeFile twiFile workDir \
                            hrmFile vrmFile hCorrList vCorrList acceleratorCode closedOrbit useDoubleOrbits \
                            useQsub bRho dispFitWeight directCalc kickX kickY numberOfQsubTasks \
                            fixedOrbitLength fitDispersion dispColumn coupledMatrix \
			    waitTime waitInterval qsubCommand queueSystemName qsubRespProcCommand submissionPause \
			    twissRefFile twissRefElement malignElement useKickFiles kickFileX kickFileY}

    if $getSDDSLattice {
        if [catch {exec $locoBinDir/locoCalculateResponseMatrix -runMode getSDDSLattice -acceleratorCode $acceleratorCode \
                       -workDir $workDir -lteFile $lteFile -beamlineName $beamlineName -bRho $bRho \
		       -latticeParamFileList $latticeParamFileList \
		       -sddsLatticeFile $sddsLatticeFile >@ stdout} result] {
            return -code error "locoCalculateResponseMatrix (mode getSDDSLattice): $result"
        }
        return 0
    }
    if $getBetaFunctions {
	if ![string length $closedOrbit] {set closedOrbit 1}
        if [catch {exec $locoBinDir/locoCalculateResponseMatrix -runMode getBetaFunctions -acceleratorCode $acceleratorCode \
                       -workDir $workDir -lteFile $lteFile -beamlineName $beamlineName -bRho $bRho \
		       -latticeParamFileList $latticeParamFileList \
		       -twiFile $twiFile \
		       -closedOrbit $closedOrbit -twissRefFile $twissRefFile \
		       -twissRefElement $twissRefElement -malignElement $malignElement >@ stdout} result] {
            return -code error "locoCalculateResponseMatrix (mode getBetaFunctions): $result"
        }
        return 0
    }
    
    #------ Main call.
    #------ List of variables which has the same name as arguments of the locoCalculateResponseMatrix:
    set globalVarList [list bRho dispFitWeight fixedOrbitLength dispColumn closedOrbit \
			   coupledMatrix \
			   lteFile beamlineName latticeParamFileList useDoubleOrbits \
			   twissRefFile twissRefElement malignElement]
    foreach var $globalVarList {
        if ![string length [set $var]] {
            unset $var
            global $var
        }
    }
    set possibleEmptyVars [list twissRefFile twissRefElement malignElement kickFileX kickFileY]
    foreach var $possibleEmptyVars {
	if ![string length [set $var]] {set $var dummy}
    }
    #------ List of the variables, which are different from the arguments of the locoCalculateResponseMatrix:
    set optionsVariablesList [list {useQsubShort useQsub} {qsubCommandShort qsubCommand} {waitTimeShort waitTime} \
				  {numberOfQsubTasksShort numberOfQsubTasks} {waitIntervalShort waitInterval} \
				  {directRMCalc directCalc} {thetaKickX kickX} {thetaKickY kickY} \
				  {submissionPauseShort submissionPause} {queueSystemNameShort queueSystemName} \
				  {qsubRespProcCommandShort qsubRespProcCommand}]
    foreach optionList $optionsVariablesList {
	set globalVar [lindex $optionList 0]
	set optionName [lindex $optionList 1]
        if ![string length [set $optionName]] {
            global $globalVar
            set $optionName [set $globalVar]
        }
    }
    if {$useQsub && ![string length $numberOfQsubTasks]} {
	return -code error "Variable numberOfQsubTasks in not defined."
    }
    if ![string length $fitDispersion] {
        global options
        set fitDispersion $options(fitDispersion)
    }
    set directDispCalcString ""
    if $fitDispersion {
	global directDispersionCalc dispDirectParamList
	if $directDispersionCalc {set directDispCalcString " -directDispersionCalc $directDispersionCalc -dispDirectParamList \"$dispDirectParamList\""}
    }
    #------ calculateResponseMatrix runs appropriate accelerator code to calculate response matrix
    set runCommand "exec $locoBinDir/locoCalculateResponseMatrix \
                   -runMode getResponseMatrix \
		   -acceleratorCode $acceleratorCode \
                   -lteFile $lteFile \
		   -beamlineName $beamlineName \
		   -latticeParamFileList \"$latticeParamFileList\" \
                   -hrmFile $hrmFile \
                   -vrmFile $vrmFile \
                   -hCorrList \"$hCorrList\" \
                   -vCorrList \"$vCorrList\" \
                   -workDir $workDir \
                   -useQsub $useQsub \
                   -bRho $bRho \
                   -twiFile $twiFile \
                   -dispFitWeight $dispFitWeight \
                   -directCalc $directCalc \
                   -kickX $kickX \
                   -kickY $kickY \
	           -useKickFiles $useKickFiles \
	           -kickFileX $kickFileX \
	           -kickFileY $kickFileY \
                   -fixedOrbitLength $fixedOrbitLength \
                   -fitDispersion $fitDispersion \
                   -dispColumn $dispColumn \
                   -coupledMatrix $coupledMatrix \
                   -elemByElemSR $elemByElemSR \
		   -closedOrbit $closedOrbit \
		   -twissRefFile $twissRefFile \
		   -twissRefElement $twissRefElement \
		   -malignElement $malignElement \
                   -useDoubleOrbits $useDoubleOrbits \
		   -numberOfQsubTasks $numberOfQsubTasks \
		   -waitTime $waitTime \
		   -waitInterval $waitInterval \
		   -qsubCommand \"[ProtectDollarSigns $qsubCommand]\" \
		   -queueSystemName $queueSystemName \
		   -submissionPause $submissionPause \
		   -qsubRespProcCommand \"[ProtectDollarSigns $qsubRespProcCommand]\"\
                   -bpmGainXFile $allElementsFile.xBpm \
                   -bpmGainYFile $allElementsFile.yBpm \
                   -bpmTiltFile $allElementsFile.bpmTilt \
                   -bpmCoefFile $allElementsFile.bpmCoef \
                   $directDispCalcString \
                   -verbose $verbose"
    OutputStatusMessage $runCommand -logFileOnly 1
    if [catch {eval $runCommand  >@ stdout} result] {
        return -code error "locoCalculateResponseMatrix (mode getResponseMatrix): $result"
    }
}

#-------------------------------------------------------------------------------------------------------------------------

proc ApplyCorrectorCalibrationCorrection {args } {
    global tmpDir dispColumn options
    APSParseArguments {hrmFile vrmFile allElementsFile hYes vYes dispersionOnly}

    set tmpRoot $tmpDir/[APSTmpString]-calibCorr
    if $dispersionOnly {
	if $options(fitDispersion) {
	    set dispGain [exec sddsprocess $allElementsFile.hCorr -pipe=out -match=col,TagName=${dispColumn}* \
			      | sdds2stream -pipe=in -col=ParameterValue]
	    exec sddsprocess $hrmFile $tmpRoot "-redef=col,$dispColumn,$dispColumn 1 $dispGain + *"
	    file copy -force $tmpRoot $hrmFile
	    file delete $tmpRoot
	}
    } else {
	set hList [list $hrmFile]
	set vList [list $vrmFile]
	foreach planeFlag [list $hYes $vYes] calcFile [list $hrmFile $vrmFile] corrFile [list $allElementsFile.hCorr $allElementsFile.vCorr] {
	    if {$planeFlag && [exec sdds2stream $corrFile -rows=bare]} {
		if [file exists $calcFile] {
		    set rows [lindex [exec sdds2stream $calcFile -rows=bare] 0]
		    if [catch {exec sddscombine $calcFile -pipe=out -merge \
				   | sddstranspose -pipe -oldColumnNames=CorrectorNameValue -newColumnNames=TagName \
				   | sddsxref -pipe $corrFile -fillIn -take=ParameterValue -match=CorrectorNameValue=ElementName \
				   | sddsprocess -pipe "-redef=col,%s,%s 1 ParameterValue + *,select=*,exclude=*Value" \
				   | sddsconvert -pipe -del=col,ParameterValue -del=para,* \
				   | sddstranspose -pipe -oldColumnNames=TagName -newColumnNames=CorrectorNameValue \
				   | sddsbreak -pipe -rowlimit=$rows \
				   | sddsxref -pipe=in $calcFile $tmpRoot.out -leave=* -transfer=para,* } result] {
			OutputStatusMessage \
			    "Error applying calibration correction (calcFile: $calcFile; corrFile: $corrFile): $result" \
			    -logFileOnly 1
			return -code error "Error applying calibration correction: $result"
		    }
		    file copy -force $tmpRoot.out $calcFile
		}
	    }
	}
	file delete $tmpRoot.out
    }
}

#-------------------------------------------------------------------------------------------------------------------------

proc AddEnergy2ResponseMatrix { args } {
    
    global tmpDir options dispColumn
    APSParseArguments { responseFile dispersionFile allElementsFile }
    
    set tmpRoot $tmpDir/[APSTmpString]-addEnergy
    set energyFile $tmpRoot.energy
    lappend deleteFiles $energyFile
    exec sddsconvert $allElementsFile.Energy $energyFile -retain=col,TagName,ParameterValue
    set twoPageFile [exec sdds2stream $responseFile -npages=bare]
    if {$twoPageFile == 2} {
        exec sddsprocess $responseFile -pipe=out -define=col,PageNumber,i_page \
            | sddsconvert -pipe -retain=col,PageNumber \
            | sddscombine -pipe=in $tmpRoot.pageNumber -merge -overWrite
        exec sddscombine $dispersionFile $tmpRoot.disp -merge
        exec sddscombine $responseFile $tmpRoot.response -merge
	lappend deleteFiles $tmpRoot.pageNumber $tmpRoot.disp $tmpRoot.response
    } else {
        exec sddsprocess $dispersionFile $tmpRoot.disp -match=para,Plane=X
        file copy -force $responseFile $tmpRoot.response
	lappend deleteFiles $tmpRoot.disp $tmpRoot.response
    }
    #------ Add zero row at the end of energy file to match the dimensions of the response file with dispersion:
    if {[lsearch -exact [join [exec sddsquery $responseFile -col] ] $dispColumn] != -1} {
	exec sddsmakedataset $tmpRoot.lastrow -col=TagName,type=string -data=$dispColumn \
	    -col=ParameterValue,type=double -data=0
	exec sddscombine $energyFile $tmpRoot.lastrow $tmpRoot.energy1 -merge
	set energyFile $tmpRoot.energy1
	lappend deleteFiles $tmpRoot.lastrow $tmpRoot.energy1
    }
    if [catch {exec sddsmatrixop $tmpRoot.disp -pipe=out -push=$energyFile -transpose -multiply -push=$tmpRoot.response -add \
		   -columnNames=file=$allElementsFile.hCorr,columnName=TagName \
		   | sddsxref -pipe=in $tmpRoot.response $tmpRoot -take=TagName -transfer=para,nux,nuy,Kick} result] {
	return -code error "Error adding matrices ($responseFile, $tmpRoot.disp, $energyFile, and $tmpRoot.response): $result"
    }
    lappend deleteFiles $tmpRoot
    if {$twoPageFile == 2} {
        exec sddsxref $tmpRoot $tmpRoot.pageNumber -pipe=out -take=PageNumber \
            | sddsbreak -pipe -changeOf=PageNumber \
            | sddsconvert -pipe=in $responseFile -del=col,PageNumber
    } else {
        file copy -force $tmpRoot $responseFile
    }
    eval file delete $deleteFiles
}

#-------------------------------------------------------------------------------------------------------------------------

proc ApplyNonlinCorrection { args } {

    APSParseArguments { computedFile measuredFile nonlinFile coupledMatrix usePage }

    global tmpDir dispColumn
    set tmpRoot $tmpDir/[APSTmpString]-nonlin

    #------ Previously useSlope was 1, but that moves part of BPM gain to nonlinearity file.
    global useIntercept useSlope useCubic
    if ![info exists useIntercept] {set useIntercept 0}
    if ![info exists useSlope] {set useSlope 0}
    if ![info exists useCubic] {set useCubic 1}
    OutputStatusMessage "NonlinCorrection: $useIntercept -- $useSlope -- $useCubic" -logFileOnly 1

    set deleteFiles ""

    #------ We don't store previous nonlinear terms - calculate full on each iteration.

    #------ Don't understand cubic term in X-Y and Y-X responses, therefore don't use nonlinear correction for
    #------ corresponsing pages.

    #------ Calculate BPM nonlinearities:
    if {$coupledMatrix == 2} {
	set page $usePage
	if {$page == 1} {set otherPage 2} else {set otherPage 1}
	exec sddsconvert $computedFile $tmpRoot.comp -keepPage=$page
	exec sddsconvert $measuredFile $tmpRoot.meas -keepPage=$page
	exec sddsconvert $nonlinFile $tmpRoot.nonlin.$otherPage -keepPage=$otherPage
	if [catch {GetBPMNonlinearity -computedResponseFile $tmpRoot.comp -measuredResponseFile $tmpRoot.meas \
		       -outputFile $tmpRoot.nonlin.$page} result] {
	    return -code error "GetBPMNonlinearity: $result"
	}
	lappend deleteFiles $tmpRoot.comp $tmpRoot.meas $tmpRoot.nonlin.$page $tmpRoot.nonlin.$otherPage
        exec sddscombine $tmpRoot.nonlin.1 $tmpRoot.nonlin.2 $tmpRoot.pfit -overWrite
    } else {
        if [catch {GetBPMNonlinearity -computedResponseFile $computedFile -measuredResponseFile $measuredFile \
                       -outputFile $tmpRoot.pfit} result] {
            return -code error "GetBPMNonlinearity: $result"
        }
    }
    lappend deleteFiles $tmpRoot.pfit
    exec sddsprocess $tmpRoot.pfit $nonlinFile \
	"-redef=col,pfitIntercept,pfitIntercept $useIntercept *" \
	"-redef=col,pfitSlope,pfitSlope $useSlope *" \
	"-redef=col,pfitCubic,pfitCubic $useCubic *"

    #------ Subtract nonlinear terms from the matrix:

    if {[lsearch -exact [exec sddsquery -col $computedFile] $dispColumn] != -1 } {
        exec sddsconvert $computedFile $tmpRoot.dispVector -retain=col,TagName,$dispColumn
        exec sddsconvert $computedFile $tmpRoot.noDisp -del=col,$dispColumn
        set addDisp 1
	lappend deleteFiles $tmpRoot.dispVector $tmpRoot.noDisp
    } else {
        file copy -force $computedFile $tmpRoot.noDisp
        set addDisp 0
	lappend deleteFiles $tmpRoot.noDisp
    }
    if [catch {exec sddsxref $tmpRoot.noDisp $nonlinFile -pipe=out -take=pfitIntercept,pfitSlope,pfitCubic -match=TagName \
		   | sddsconvert -pipe -del=col,TagName,BPMName \
		   | sddsprocess -pipe \
		   "-redef=col,%s,%s pfitIntercept $useIntercept * + %s pfitSlope $useSlope * * + %s 3 pow pfitCubic $useCubic * * +,select=*,exclude=pfit*" \
		   | sddsxref -pipe $computedFile -take=TagName \
		   | sddsconvert -pipe=in $tmpRoot -del=col,pfit*} result] {
	return -code error "Error completing nonlin calculations ($tmpRoot.noDisp; $nonlinFile; $computedFile): $result"
    }
    lappend deleteFiles $tmpRoot
    if $addDisp {
        exec sddsxref $tmpRoot $tmpRoot.dispVector $computedFile -take=$dispColumn
    } else {
        file copy -force $tmpRoot $computedFile
    }
    eval file delete $deleteFiles
}

#-------------------------------------------------------------------------------------------------------------------------

proc GetBPMNonlinearity { args } {

    global tmpDir dispColumn
    APSParseArguments {measuredResponseFile computedResponseFile outputFile}

    set tmpRoot $tmpDir/[APSTmpString]-nonlin
    set measuredFile $tmpRoot.measured
    set computedFile $tmpRoot.computed
    if {[lsearch -exact [exec sddsquery -col $measuredResponseFile] $dispColumn] != -1 } {
        exec sddsconvert $measuredResponseFile $measuredFile -del=col,$dispColumn
        exec sddsconvert $computedResponseFile $computedFile -del=col,$dispColumn
    } else {
        file copy -force $measuredResponseFile $measuredFile
        file copy -force $computedResponseFile $computedFile
    }
    exec sddscollect $computedFile $tmpRoot.collect -nowarning -collect=match=*,exclude=*Name,column=Orbit,edit=
    exec sddsmakedataset $tmpRoot.col -col=CoefficientName,type=string -data=[join [list pfitIntercept pfitSlope pfitCubic] ,]
    if [catch {exec sddschanges $measuredFile -pipe=out -base=$computedFile -copy=TagName -changesIn=exclude=*Name,* \
		   | sddscollect -pipe -nowarning -collect=prefix=ChangeIn \
		   | sddsxref -pipe $tmpRoot.collect -take=Orbit -match=Rootname \
		   | sddsconvert -pipe -rename=col,ChangeIn=Error \
		   | tee $tmpRoot.errors \
		   | sddspfit -pipe -col=Orbit,Error -orders=0,1,3 -generate \
		   | sddsxref -pipe $tmpRoot.errors -transfer=para,TagName -rename=para,TagName=TagNameParam \
		   | sddsconvert -pipe -del=col,* \
		   | sddsarray2column -pipe -convert=Coefficient \
		   | sddsxref -pipe $tmpRoot.col -take=CoefficientName -reuse=page \
		   | sddstranspose -pipe -newColumnName=CoefficientName \
		   | sddsprocess -pipe -print=col,TagName,%s,TagNameParam \
		   | sddsconvert -pipe -del=col,OldColumnNames -del=para,* \
		   | sddscombine -pipe=in $outputFile -merge} result] {
	#------ Try again:
	if [catch {exec sddschanges $measuredFile -pipe=out -base=$computedFile -copy=TagName -changesIn=exclude=*Name,* \
		       | sddscollect -pipe -nowarning -collect=prefix=ChangeIn \
		       | sddsxref -pipe $tmpRoot.collect -take=Orbit -match=Rootname \
		       | sddsconvert -pipe -rename=col,ChangeIn=Error \
		       | tee $tmpRoot.errors \
		       | sddspfit -pipe -col=Orbit,Error -orders=0,1,3 -generate \
		       | sddsxref -pipe $tmpRoot.errors -transfer=para,TagName -rename=para,TagName=TagNameParam \
		       | sddsconvert -pipe -del=col,* \
		       | sddsarray2column -pipe -convert=Coefficient \
		       | sddsxref -pipe $tmpRoot.col -take=CoefficientName -reuse=page \
		       | sddstranspose -pipe -newColumnName=CoefficientName \
		       | sddsprocess -pipe -print=col,TagName,%s,TagNameParam \
		       | sddsconvert -pipe -del=col,OldColumnNames -del=para,* \
		       | sddscombine -pipe=in $outputFile -merge} result] {
	    return -code error "Error performing nonlinear fitting - some BPMs might have zero readings only \
                            (measuredFile: $measuredFile; computedFile: $computedFile; tmpRoot: $tmpRoot): $result"
	}
    }
    file delete $tmpRoot.collect $tmpRoot.col $tmpRoot.errors $measuredFile $computedFile
}

#-------------------------------------------------------------------------------------------------------------------------

proc SubtractMatricesAndMakeOneColumn { args } {
    
    APSParseArguments {hYes vYes hrmMeasured vrmMeasured hrmIterFile vrmIterFile hrmDifference vrmDifference \
                           xyColumnFile useQuadConstraints quadConstraintWeight coupledMatrix}
    
    global tmpDir dispColumn weightFile offDiagWeight options
    
    set tmpRoot $tmpDir/[APSTmpString]-make1col
    set planeFlagList [list $hYes $vYes]
    set measuredFileList [list $hrmMeasured $vrmMeasured]
    set iterFileList [list $hrmIterFile $vrmIterFile]
    set planeList [list X Y]
    set differenceFileList [list $hrmDifference $vrmDifference]
    set deleteFiles ""
    foreach planeFlag $planeFlagList measuredFile $measuredFileList iterFile $iterFileList plane $planeList \
        differenceFile $differenceFileList {
            if $planeFlag {
                if [catch {exec sddschanges $measuredFile -pipe=out -base=$iterFile -copy=TagName -changes=exclude=*Name,* \
			       -parallelPages \
			       | sddsconvert -pipe=in $differenceFile "-editNames=col,Change*,a 8d"} result] {
		    return -code error "Error calculating difference (measuredFile - $measuredFile; iterFile - $iterFile): $result"
		}
                if [catch {ApplyWeightCorrection -differenceFile $differenceFile -weightFile $weightFile \
                               -plane $plane} result] {
                    return -code error "ApplyWeightCorrection: $result"
                }
                if {[exec sdds2stream $differenceFile -npages=bare] == 2} {
                    #------ Applying offDiagWeight...
                    if {$plane == "X"} { set pageNum 2 } else { set pageNum 1 }
                    exec sddsprocess $differenceFile $tmpRoot.dif \
                        "-redef=col,%s,i_page $pageNum == pop pop ? %s $offDiagWeight * : %s \$,select=*,exclude=*Name"
                } else {
                    file copy -force $differenceFile $tmpRoot.dif
                }
		lappend deleteFiles $tmpRoot.dif
                if [catch {Fit_MakeColumnFromMatrix -matrixFile $tmpRoot.dif -columnFile $tmpRoot.$plane \
                               -columnName MatrixDifference -rowNameColumn TagName -deleteColumns BPMName} result] {
                    return -code error "Fit_MakeColumnFromMatrix: $result"
                }
            } else {
                exec sddsmakedataset $tmpRoot.$plane -col=Rootname,type=string -data \
                    -col=MatrixDifference,type=double -data
            }
        }
    set combineFiles "$tmpRoot.X $tmpRoot.Y"

    #------ Vector of tune differences:
    if $options(useTune) {
	if [catch {MakeTunesColumn -filename $tmpRoot.tunes} result] {
	    return -code error "MakeTunesColumn: $result"
	}
	lappend combineFiles $tmpRoot.tunes
    }

    #------ Calculating vector of quad constraint additions...
    if $useQuadConstraints {
	global variableFile allElementsFile dontUsePrevQuadConstr
	set varFile $tmpRoot.QuadSkew
	lappend deleteFiles $tmpRoot.QuadSkew
	#------ allElementsFile contains data from previous run and might contain more elements
	if $dontUsePrevQuadConstr {set quadValueFile $variableFile} else {set quadValueFile $allElementsFile}
	if ![exec sdds2stream $quadValueFile.Quad -rows=bare] {
	    return -code error "Quad constraints should not be used if there are no quads to fit."
	}
	if {$coupledMatrix == 2 && $options(useSkew) && $useQuadConstraints == 2} {
	    exec sddsselect $quadValueFile.Quad $variableFile.Quad $tmpRoot.Quad -match=TagName
	    exec sddsselect $quadValueFile.Skew $variableFile.Skew $tmpRoot.Skew -match=TagName
	    exec sddscombine $tmpRoot.Quad $tmpRoot.Skew $varFile -merge -overWrite
	    lappend deleteFiles $tmpRoot.Quad $tmpRoot.Skew
	} else {
	    exec sddsselect $quadValueFile.Quad $variableFile.Quad $varFile -match=TagName
	}
	if [catch {exec sddsconvert $varFile -pipe=out \
		       -retain=col,TagName,ParameterValue -rename=col,ElementName=Rootname,ParameterValue=MatrixDifference \
		       | sddsprocess -pipe=in $tmpRoot.quads \
		       "-redef=col,MatrixDifference,MatrixDifference $quadConstraintWeight * -1 *"} $result] {
	    return -code error "Error making constraint vector (varFile - $varFile): $result"
	}
	lappend combineFiles $tmpRoot.quads
	set constraintLines [exec sdds2stream $tmpRoot.quads -rows=bare]
    } else {
	set constraintLines 0
    }
    eval exec sddscombine $combineFiles -pipe=out -merge \
	| sddsprocess -pipe=in $xyColumnFile -redef=para,ConstraintLines,$constraintLines,type=long

    if [catch {CalculateResidualRms -vectorX $tmpRoot.X -vectorY $tmpRoot.Y \
                   -measuredFileX $hrmMeasured -measuredFileY $vrmMeasured -coupledMatrix $coupledMatrix \
		   -vectorQuads $tmpRoot.quads -quadConstraintWeight $quadConstraintWeight} rmsList] {
        return -code error "CalculateResidualRms: $rmsList"
    }
    eval file delete $combineFiles $deleteFiles
    return $rmsList
}

#-------------------------------------------------------------------------------------------------------------------------

proc ApplyWeightCorrection {args} {

    APSParseArguments {differenceFile weightFile plane}

    global tmpDir useBpmWeight weightCorrectorsFile sigmaWeightFileX sigmaWeightFileY useSigmaWeightFile
    
    set deleteFiles ""
    if $useBpmWeight {
        set tmpFile $tmpDir/[APSTmpString].weight
        if {[exec sdds2stream $differenceFile -npages=bare] != 2} {
            exec sddsprocess $weightFile $tmpFile -match=para,Plane=$plane
        } else {
            file copy -force $weightFile $tmpFile
        }
        exec sddsxref $tmpFile $differenceFile -pipe=out -leave=TagName -match=TagName -nowarning \
            | sddsconvert -pipe -del=col,TagName \
            | sddsprocess -pipe "-redef=col,%s,%s Weight *,select=*,exclude=Weight" \
            | sddsconvert -pipe=in $tmpFile.1 -del=col,Weight
        exec sddsxref $tmpFile.1 $differenceFile $tmpFile.2 -take=TagName
        file copy -force $tmpFile.2 $differenceFile
        lappend deleteFiles $tmpFile $tmpFile.1 $tmpFile.2
    }
    
    #------ Adjusting for corrector weights...
    if [file exists $weightCorrectorsFile] {
        OutputStatusMessage "Applying corrector weight..."
        set correctorList [join [exec sdds2stream $weightCorrectorsFile -col=Corrector] ]
        set weightList [join [exec sdds2stream $weightCorrectorsFile -col=Weight] ]
        set columnList [join [exec sddsquery -col $differenceFile] ]
        set redefLine ""
        foreach corrector $correctorList weight $weightList {
            if {[lsearch -exact $columnList $corrector] != -1} {
                append redefLine " \"-redef=col,$corrector,$corrector $weight *\" "
            }
        }
        if [string length $redefLine] {
            eval exec sddsprocess -nowarning $differenceFile $redefLine
            lappend deleteFiles ${differenceFile}~
        }
    }

    #------ Applying sigmaWeightFile:
    if $useSigmaWeightFile {
	switch -exact $plane {
	    X {set sigmaWeightFile $sigmaWeightFileX}
	    Y {set sigmaWeightFile $sigmaWeightFileY}
	    default {return -code error "Wrong plane variable: $plane"}
	}
	if [file exists $sigmaWeightFile] {
	    set tmpFile $tmpDir/[APSTmpString].sigmaweight
	    #------ sddsmatrixop deletes all string columns.
	    exec sddsmatrixop $differenceFile -pipe=out -push=$sigmaWeightFile -divide=hadamard \
		| sddsxref -pipe=in $differenceFile $tmpFile -take=TagName
	    file copy -force $tmpFile $differenceFile
	    lappend deleteFiles $tmpFile
	}
    }
    catch {eval file delete $deleteFiles}
}

#-------------------------------------------------------------------------------------------------------------------------

proc CalculateResidualRms { args } {
    
    global dispColumn thetaKickX thetaKickY dispFitWeight offDiagWeight options
    APSParseArguments {vectorX vectorY measuredFileX measuredFileY vectorQuads quadConstraintWeight coupledMatrix}
    
    set vectorFileList [list $vectorX $vectorY]
    set thetaKickList [list $thetaKickX $thetaKickY]
    set measuredFileList [list $measuredFileX $measuredFileY]
    set planeList [list X Y]
    foreach vectorFile $vectorFileList thetaKick $thetaKickList measuredFile $measuredFileList plane $planeList {
        if [exec sddscombine $vectorFile -pipe=out -merge | sdds2stream -pipe=in -rows=bare] {
            #------ Undoing offDiagWeight before calculating rms...
            if {[exec sdds2stream $vectorFile -npages=bare] == 2} {
                if {$plane == "X"} { set pageNum 2 } else { set pageNum 1 }
                exec sddsprocess $vectorFile $vectorFile.tmp \
                    "-redef=col,%s,i_page $pageNum == pop pop ? %s $offDiagWeight / : %s \$,select=*,exclude=Rootname"
            } else {
                file copy -force $vectorFile $vectorFile.tmp
            }
            #------ This is for dispersion-dependent processing...
            if {[lsearch -exact [exec sddsquery -col $measuredFile] $dispColumn] != -1 } {
                #------ Orbit rms in mm (transformation from m/rad to mm).
                set stdDevZ [exec sddsprocess $vectorFile.tmp -pipe=out "-match=col,Rootname=*${dispColumn}*,!" \
                                 | sddsprocess -pipe \
                                 "-redef=col,MatrixDifference,MatrixDifference $thetaKick 1000 * *" \
                                 | tee $vectorFile.Zmm \
                                 | sddsprocess -pipe -process=MatrixDifference,standardDeviation,stdDev \
                                 | sdds2stream -pipe=in -para=stdDev]
                #------ Dispersion rms in m.
                set stdDevD [exec sddsprocess $vectorFile.tmp -pipe=out "-match=col,Rootname=*${dispColumn}*" \
                                 | sddsprocess -pipe \
                                 "-redef=col,MatrixDifference,MatrixDifference $dispFitWeight /" \
                                 | tee $vectorFile.Dm \
                                 | sddsprocess -pipe -process=MatrixDifference,standardDeviation,stdDev \
                                 | sdds2stream -pipe=in -para=stdDev]
                set stdDev$plane [list $stdDevZ $stdDevD]
                exec sddscombine $vectorFile.Zmm $vectorFile.Dm $vectorFile.mm -merge
                file delete $vectorFile.Zmm $vectorFile.Dm
            } else {
                exec sddsprocess $vectorFile.tmp $vectorFile.mm \
                    "-redef=col,MatrixDifference,MatrixDifference $thetaKick 1000 * *"
                if [catch {exec sddsprocess $vectorFile.mm -pipe=out \
                               -process=MatrixDifference,standardDeviation,stdDev \
                               | sdds2stream -pipe=in -para=stdDev} stdDevZ] {
                    return -code error "Error processing file $vectorFile.mm: $stdDevZ"
                }
                set stdDev$plane [list $stdDevZ]
            }
            file delete $vectorFile.tmp
        } else {
            set stdDev$plane Nan
            file copy -force $vectorFile $vectorFile.mm
        }
    }
    set combineFiles "$vectorX.mm $vectorY.mm"

    #------ Calculating quad constraint std...
    if [file exists $vectorQuads] {
	#------ Find stdev without quads first:
	if [catch {eval exec sddscombine $combineFiles -pipe=out -merge \
		       | sddsprocess -pipe -process=MatrixDifference,standardDeviation,stdDev \
		       | sdds2stream -pipe=in -para=stdDev} stdDev1] {
	    return -code error "Error combining files (1) ($combineFiles): $stdDev1"
	}
	#------ Divide by weight (actually need to involve quadConstraintFile here):
	if [catch {exec sddsprocess $vectorQuads -pipe=out \
		       "-redef=col,MatrixDifference,MatrixDifference $quadConstraintWeight /" \
		       -process=MatrixDifference,standardDeviation,StdDev \
		       | tee $vectorQuads.mm \
		       | sdds2stream -pipe=in -para=StdDev} stdDevQ] {
	    return -code error "Error calculating quad vector: $stdDevQ"
	}
	lappend combineFiles $vectorQuads.mm
    } else {
	set stdDevQ 0
    }
    
    if [catch {eval exec sddscombine $combineFiles -pipe=out -merge \
		   | sddsprocess -pipe -process=MatrixDifference,standardDeviation,stdDev \
		   | sdds2stream -pipe=in -para=stdDev} stdDev] {
	return -code error "Error combining files (2) ($combineFiles): $stdDev"
    }
    eval file delete $combineFiles
    if ![info exists stdDev1] {set stdDev1 $stdDev}
    
    set returnList ""
    if {$coupledMatrix == 2} {
	set returnList [concat $returnList "rmsYX [lindex [lindex $stdDevY 0] 0] rmsYY [lindex [lindex $stdDevY 0] 1]"]
	set returnList [concat $returnList "rmsXX [lindex [lindex $stdDevX 0] 0] rmsXY [lindex [lindex $stdDevX 0] 1]"]
	if $options(fitDispersion) {
	    set returnList [concat $returnList "rmsDX [lindex [lindex $stdDevX 1] 0] rmsDY [lindex [lindex $stdDevX 1] 1]"]
	}
    } else {
	set returnList [concat $returnList "rmsYY $stdDevY"]
	if $options(fitDispersion) {
	    set returnList [concat $returnList "rmsXX [lindex $stdDevX 0] rmsDX [lindex $stdDevX 1]"]
	} else {
	    set returnList [concat $returnList "rmsXX [lindex $stdDevX 0]"]
	}
    }
    set returnList [concat $returnList "rmsTotal $stdDev rmsTotalWQ $stdDev1 rmsQ $stdDevQ"]
    return $returnList
}

#-------------------------------------------------------------------------------------------------------------------------

proc MakeTunesColumn {args} {
    APSParseArguments {filename}
    global tmpDir nuxMeasured nuyMeasured twiFile tuneFitWeight
    set calculatedTunes [exec sdds2stream -para=nux,nuy $twiFile]
    if {[lsearch -glob $calculatedTunes *NaN*] != -1} {
        return -code error "Error: Tunes are unstable!"
    }
    exec sddsmakedataset -pipe=out \
	-col=Rootname,type=string -data=Snux,Snuy \
	-col=CalcTunes,type=double -data=[join $calculatedTunes ,] \
	-col=MeasTunes,type=double -data=$nuxMeasured,$nuyMeasured \
	| sddsprocess -pipe "-redef=col,MatrixDifference,MeasTunes CalcTunes - $tuneFitWeight *" \
	| sddsconvert -pipe=in $filename -del=col,CalcTunes,MeasTunes
    OutputStatusMessage [format "%s %10.4f %10.4f" "Calculated Tunes are " [lindex $calculatedTunes 0] [lindex $calculatedTunes 1]]
}

#-------------------------------------------------------------------------------------------------------------------------
#-------------------------------------------------------------------------------------------------------------------------
#------ Response Matrix Derivative related Procedures --------------------------------------------------------------------
#-------------------------------------------------------------------------------------------------------------------------
#-------------------------------------------------------------------------------------------------------------------------

proc BuildResponseMatrixDeriv { args } {
    
    global elementListVar elementListRM rmDerivElementFile allElementsFile
    global latticeIterFile latticeRMFile dispersionFile twiFile variableFile rmDifferenceFile
    global options kickRMD fixedOrbitLength directRMCalc dispColumn dispFitWeight bRho 
    global waitTimeLong waitIntervalLong numberOfQsubTasksLong useDoubleOrbits
    global workDir tmpDir locoBinDir qsubCommandLong submissionPauseLong
    global hYes vYes closedOrbit qsubRespProcCommandLong queueSystemNameLong
    global rmdQuadDelta rmdTiltDelta remotePath verbose useDirectRMDCalc
    global twissRefElement twissRefFile malignElement specialElementsInputFile

    set verboseLocal 1
    APSParseArguments {continuePrevious abortFile rmDerivFile useQsubLong hrmFile vrmFile lteFile \
			   coupledMatrix acceleratorCode beamlineName solutionParamFileList \
			   useKickFiles kickFileX kickFileY hrmFile2plane vrmFile2plane verboseLocal}
    
    if ![file exists $workDir/dqs] {exec mkdir $workDir/dqs}
    
    OutputStatusMessage "Doing the orbit response matrix derivative..."

    set uniqueMatrixCode [APSTmpString]
    if [catch {MakeRMDerivElementListFile -varList $elementListVar -rmList $elementListRM -outputFile $rmDerivElementFile \
		   -uniqueMatrixCode $uniqueMatrixCode -coupledMatrix $coupledMatrix -rmDerivFile $rmDerivFile} result] {
	return -code error "MakeElementListFile: $result"
    }

    if $verboseLocal {set verboseString ""} else {set verboseString " -logFileOnly 1 "}

    set deleteFiles ""
    set matrixDoneList ""
    set matrixErrorList ""
    set matrixArrayList ""
    set hrmFile00 $workDir/computed.hrm.000
    set vrmFile00 $workDir/computed.vrm.000
    exec sddscombine $hrmFile $hrmFile00 -merge -overWrite
    exec sddscombine $vrmFile $vrmFile00 -merge -overWrite
    lappend deleteFiles $hrmFile00 $vrmFile00
    set tmpRoot $tmpDir/[APSTmpString]-RMD
    set matrixAbortFile $tmpRoot.abort
    switch -exact -- $coupledMatrix {
        0 { 
	    set indexList [list useQuads useCorrs useBPMs useEnergy] 
	    set arrayIndexList [list quads hCorr vCorr xBpm yBpm energy]
	}
        1 { 
	    set indexList [list useSkew useCorrsTilt useBPMsTilt] 
	    set arrayIndexList [list skew hCorrTilt vCorrTilt xBpmTilt yBpmTilt]
	}
        2 { 
	    set indexList [list useQuads useSkew useCorrs useBPMs useEnergy useCorrsTilt useBPMsTilt useBPMsCoef useSpecElem]
	    set arrayIndexList [list quads hCorr vCorr xBpm yBpm energy skew hCorrTilt vCorrTilt bpmTilt bpmCoef special]
	}
    }
    if ![file exists $workDir/dqs] {exec mkdir $workDir/dqs}
    
    foreach {directCalcList indirectCalcList} [list "" ""] {}
    foreach index $indexList {
	if $options($index) {
	    switch -exact -- $index {
		useQuads      {lappend directCalcList quads}
		useSkew       {lappend directCalcList skew}
		useCorrs      {if {$useDirectRMDCalc || $useKickFiles} {lappend directCalcList hCorr vCorr} else {lappend indirectCalcList hCorr vCorr}}
		useBPMs       {if $useDirectRMDCalc {lappend directCalcList xBpm yBpm} else {lappend indirectCalcList xBpm yBpm}}
		useCorrsTilt  {if {$useDirectRMDCalc || $useKickFiles} {lappend directCalcList hCorrTilt vCorrTilt} else {lappend indirectCalcList hCorrTilt vCorrTilt}}
		useBPMsTilt   {if $useDirectRMDCalc {lappend directCalcList bpmTilt} else {lappend indirectCalcList bpmTilt}}
		useBPMsCoef   {if $useDirectRMDCalc {lappend directCalcList bpmCoef} else {lappend indirectCalcList bpmCoef}}
		fitPathLength {lappend directCalcList fitPathLength}
		useSpecElem   {lappend directCalcList special}
	    }
	}
    }
    
    #------ Running direct RMD calculations:
    set splitTasksLocal $numberOfQsubTasksLong 
    foreach varType $directCalcList {
	set ext $varType
	set colPrefix ""
	set colSuffix ""
	set parameterDelta 1e-4
	set localVarFile $variableFile.$varType
	set statusMessage "Doing direct ${varType}..."
	switch -exact -- $varType {
	    special {
		set varList [join [exec sdds2stream $specialElementsInputFile -col=TagName] ]
		set localVarFile $variableFile.special
		set splitTasksLocal 1
	    }
	    quads {
		set varList [lindex $elementListVar 0]
		set localVarFile $variableFile.Quad
		set parameterDelta $rmdQuadDelta
	    }
	    skew {
		set varList [lindex $elementListVar 6]
		set localVarFile $variableFile.Skew
		set parameterDelta $rmdTiltDelta
		set colPrefix a
	    }
	    hCorr {set varList [lindex $elementListVar 1]}
	    vCorr {set varList [lindex $elementListVar 2]}
	    hCorrTilt {
		set varList [lindex $elementListVar 7]
		set colSuffix "#tilt"
	    }
	    vCorrTilt {
		set varList [lindex $elementListVar 8]
		set colSuffix "#tilt"
	    }
	    xBpm {set varList [lindex $elementListVar 3]}
	    yBpm {set varList [lindex $elementListVar 4]}
	    bpmTilt {
		set varList [lindex $elementListVar 9]
		set colSuffix "#tilt"
	    }
	    bpmCoef {
		set varList [lindex $elementListVar 10]
		set colSuffix "#coef"
	    }
	    fitPathLength {
		set ext dz
		set varList "MALIN#1"
		set localVarFile $variableFile.$ext
		set parameterDelta 1e-6
	    }
	    default {return -code error "Error in index: varType= $varType"}
	}
	eval OutputStatusMessage \"$statusMessage\" $verboseString
	if [string length $remotePath] {set remoteQueue 1} else {set remoteQueue 0}
	set queueDirName $workDir/dqs/$ext
	set localTmpDir $tmpDir/$ext
	set remotePath1 $remotePath/$ext
	set matrixFile $tmpRoot.$ext
	set rootTaskName $ext
	#------ Need "dummy" here because in string scriptParameters zero length variable means just empty space.
	set possibleEmptyVars [list lteFile beamlineName solutionParamFileList twissRefFile twissRefElement malignElement kickFileX kickFileY]
	foreach var $possibleEmptyVars {if ![string length [set $var]] {set $var dummy}}

	set scriptName $locoBinDir/locoRunQuads
	set scriptParameters    "-parameterDelta $parameterDelta -acceleratorCode $acceleratorCode "
	append scriptParameters "-tuneInclude $options(useTune) -kick $kickRMD -bRho $bRho "
	append scriptParameters "-fixedOrbitLength $fixedOrbitLength -fitDispersion $options(fitDispersion) "
	append scriptParameters "-dispColumn $dispColumn -dispFitWeight $dispFitWeight "
	append scriptParameters "-coupledMatrix $coupledMatrix "
	append scriptParameters "-directRMCalc $directRMCalc "
	append scriptParameters "-makeBaselineCalc 1 "
	append scriptParameters "-lteFile $lteFile -beamlineName $beamlineName "
	append scriptParameters "-computedRMX $hrmFile00 -computedRMY $vrmFile00 -closedOrbit $closedOrbit "
	append scriptParameters "-useDoubleOrbits $useDoubleOrbits -latticeParamFileList \"$solutionParamFileList\" "
	append scriptParameters "-twissRefFile $twissRefFile -twissRefElement $twissRefElement -malignElement $malignElement "
	append scriptParameters "-useKickFiles $useKickFiles -kickFileX $kickFileX -kickFileY $kickFileY "
	append scriptParameters "-localTmpDir $localTmpDir -variableFile $localVarFile "
	if {[string first "bpm" $varType] != -1 || [string first "Bpm" $varType] != -1} {append scriptParameters "-bpmMode 1 "}
	if [string match special $varType] {append scriptParameters "-specialElementsInputFile $specialElementsInputFile "}
	if ![llength $varList] {
	    OutputStatusMessage "BuildResponseMatrixDeriv Warning: No variables specified for $index option (elementListVar: $elementListVar)"
	} else {
	    set doneFile $matrixFile.done
	    set errorFile $matrixFile.error
	    set optionsFile $matrixFile.options
	    if [catch {CalculateOrbitRMD -varList $varList \
			   -splitTasks $splitTasksLocal -colPrefix $colPrefix -colSuffix $colSuffix \
			   -scriptName $scriptName -scriptParameters $scriptParameters \
			   -matrixFile $matrixFile -queueDirName $queueDirName -rootTaskName $rootTaskName \
			   -useQsub $useQsubLong -continuePrevious $continuePrevious -usePopupWindow 0 \
			   -waitTime $waitTimeLong -qsubCommand $qsubCommandLong \
			   -waitInterval $waitIntervalLong -submissionPause $submissionPauseLong \
			   -waitForCompletion 0 -doneFile $doneFile -errorFile $errorFile \
			   -matrixAbortFile $matrixAbortFile -optionsFile $optionsFile -remotePath $remotePath1 \
			   -lteFile $lteFile -beamlineName $beamlineName \
			   -useKickFiles $useKickFiles -kickFileX $kickFileX -kickFileY $kickFileY \
			   -solutionParamFileList $solutionParamFileList -remoteQueue $remoteQueue \
			   -qsubRespProcCommand $qsubRespProcCommandLong -queueSystemName $queueSystemNameLong \
			   -abortFile $abortFile -verboseLocal $verboseLocal} result] {
		return -code error "CalculateOrbitRMD: $result"
	    }
	}
	lappend matrixDoneList $doneFile; lappend deleteFiles $doneFile
	lappend matrixErrorList $errorFile
	lappend matrixArrayList $ext $matrixFile
    }

    #------ Running indirect RMD calculations:
    foreach varType $indirectCalcList {
	set colPrefix ""
	set colSuffix ""
	set ext $varType
	set parameterDelta 1e-4
	switch -regexp -- $varType {
	    Corr$ {
		set statusMessage "Doing indirect $varType correctors..."
		set rmFileString " -hrmFile $hrmFile -vrmFile $vrmFile "
	    }
	    CorrTilt$ {
		set statusMessage "Doing indirect $varType corrector tilts..."
		set colSuffix "#tilt"
		set tiltedHRMFile $tmpRoot.tiltedHRM
		set tiltedVRMFile $tmpRoot.tiltedVRM
		set rmFileString " -hrmFile $tiltedHRMFile -vrmFile $tiltedVRMFile "
		lappend deleteFiles $tiltedHRMFile $tiltedVRMFile
		#------ Preparing files for QuickRMDeriv -- need response matrix where all correctors are tilted by delta
		if [catch {CalculateTiltedCorrResponse -lteFile $lteFile -beamlineName $beamlineName -delta $parameterDelta \
			       -latticeParamFileList $solutionParamFileList \
			       -coupledMatrix $coupledMatrix -tiltedHRMFile $tiltedHRMFile \
			       -tiltedVRMFile $tiltedVRMFile} result] {
		    return -code error "CalculateTiltedCorrResponse: $result"
		}
	    }
	    Bpm$ {
		set statusMessage "Doing indirect $varType calibrations..."
		set rmFileString " -hrmFile $hrmFile -vrmFile $vrmFile -hrmFile00 $hrmFile "
	    }
	    bpmTilt|bpmCoef {
		set statusMessage "Doing indirect $varType..."
		set rmFileString " -hrmFile $hrmFile2plane -vrmFile $vrmFile2plane -hrmFile00 $hrmFile "
		set colSuffix "[string range $varType 3 6]"
	    }
	    default {return -code error "Error in index: varType= $varType"}
  	}
	eval OutputStatusMessage \"$statusMessage\" $verboseString
	set matrixFile $tmpRoot.$ext
	lappend matrixArrayList $ext $matrixFile
	if [catch {eval QuickRMDeriv $rmFileString -mode $ext -outputFile $matrixFile} result] {
	    return -code error "QuickRMDeriv: $result"
	}
	OutputStatusMessage "QuickRMDeriv is done." -logFileOnly 1
    }
if 0 {
    if $options(useSpecElem) {
	if [exec sdds2stream $specialElementsInputFile -rows=bare] {
	    OutputStatusMessage "Doing special elements..."
	    set matrixFile $tmpRoot.special
	    lappend matrixArrayList special $matrixFile
#	    set specialElementList [join [exec sdds2stream $specialElementsInputFile -col=ElementName] ]
#	    set specialElementScripts [join [exec sdds2stream $specialElementsInputFile -col=ScriptName] ]
#	    set specialElementDeltas [join [exec sdds2stream $specialElementsInputFile -col=Delta] ]
#	    set specialElementSuffixes [join [exec sdds2stream $specialElementsInputFile -col=NameSuffix] ]
	    
	    set localTmpDir $tmpDir
	    if ![file exists $localTmpDir] {exec mkdir $localTmpDir}
	    if ![file exists $workDir/dqs/special] {exec mkdir $workDir/dqs/special}
	    set scriptName $locoBinDir/locoRunQuads
	    set totalPath $workDir/dqs/special
	    if ![file exists $totalPath] {exec mkdir $totalPath}
	    set localMatrixFile $totalPath/responseMatrix.sdds
	    set scriptParameters    "-acceleratorCode $acceleratorCode "
	    append scriptParameters "-tuneInclude $options(useTune) -kick $kickRMD -bRho $bRho "
	    append scriptParameters "-fixedOrbitLength $fixedOrbitLength -fitDispersion $options(fitDispersion) "
	    append scriptParameters "-dispColumn $dispColumn -dispFitWeight $dispFitWeight "
	    append scriptParameters "-coupledMatrix $coupledMatrix "
	    append scriptParameters "-directRMCalc $directRMCalc "
	    append scriptParameters "-makeBaselineCalc 1 "
	    append scriptParameters "-lteFile $lteFile -beamlineName $beamlineName "
	    append scriptParameters "-computedRMX $hrmFile00 -computedRMY $vrmFile00 -closedOrbit $closedOrbit "
	    append scriptParameters "-useDoubleOrbits $useDoubleOrbits -latticeParamFileList \"$solutionParamFileList\" "
	    append scriptParameters "-twissRefFile $twissRefFile -twissRefElement $twissRefElement -malignElement $malignElement "
	    if $useKickFiles {append scriptParameters "-useKickFiles $useKickFiles -kickFileX $kickFileX -kickFileY $kickFileY "} else {append scriptParameters "-useKickFiles $useKickFiles"}
	    append scriptParameters "-localTmpDir $localTmpDir -specialElementsInputFile $specialElementsInputFile "

	    set scriptParametersFile $totalPath/scriptParameters
	    if [catch {open $scriptParametersFile w} fid] { 
		return -code error "$fid" 
	    }
	    puts $fid "-path $totalPath $scriptParameters"
	    close $fid
	    if [catch {exec $mainScript -argMode file -scriptParametersFile $scriptParametersFile -verbose $verbose} result] {
		return -code error "$mainScript: $result"
	    }
	    if {[llength $specialElementMatrixList] == 1} {
		file copy -force $specialElementMatrixList $tmpRoot.special
	    } else {
		eval exec sddsxref $specialElementMatrixList $tmpRoot.special -take=* -leave=Rootname -nowarning
	    }
	    lappend matrixArrayList special $matrixFile
	    lappend arrayIndexList special
	}
    }
}
    #------ Waiting for the matrix done files to appear:
    if [llength $matrixDoneList] {
	eval OutputStatusMessage \"Waiting for background jobs to complete...\" $verboseString
	set done 0
	set waitInterval [expr $waitIntervalLong * 1000]
	while {$done == 0} {
	    set done 1
	    update
	    after $waitInterval
	    #------ Checking if the done files exist...
	    foreach filename $matrixDoneList {
		if ![file exists $filename] { set done 0 }
	    }
	    #------ Checking if error files have appeared...
	    foreach filename $matrixErrorList {
		set existingFiles ""
		if [file exists $filename] { lappend existingFiles $filename }
		if [llength $existingFiles] {
		    return -code error "Error running background jobs: see files $existingFiles"
		}
	    }
	    #------ Aborting if requested...
	    if [file exists $abortFile] {
		Fit_WriteToFile -filename $matrixAbortFile -accessMode w -line "Error: Interrupted by user."
		return -code error "Interrupted by user."
	    }
	}
    }
    
    eval OutputStatusMessage \"Final combining...\" $verboseString
    file delete $rmDerivFile
    if 1 {
	#------ Make properly sorted matrixFile list:
	array set matrixArray $matrixArrayList
	set matrixFileList ""
	set existingIndexList [array names matrixArray]
	foreach arrayIndex $arrayIndexList {
	    if {[lsearch -exact $existingIndexList $arrayIndex] != -1} {
		lappend matrixFileList $matrixArray($arrayIndex)
	    }
	}
	set constraintLines [exec sdds2stream $rmDifferenceFile -para=ConstraintLines]
	OutputStatusMessage "List of combined files: $matrixFileList" -logFileOnly 1
	eval exec sddsprocess $rmDifferenceFile -pipe=out -clip=0,$constraintLines \
	    | sddsconvert -pipe -retain=col,Rootname \
	    | sddsxref -pipe $matrixFileList -take=* -leave=Rootname \
	    | sddsprocess -pipe=in $rmDerivFile -print=para,MatrixCode,$uniqueMatrixCode \
	    -print=para,ElementFile,[file tail $rmDerivElementFile] -def=para,CoupledMatrixAnalysis,$coupledMatrix
    } else {
	set tmpDerivFile $tmpDir/[file tail $rmDerivFile]
	file delete $rmDerivFile $tmpDerivFile
	if {[llength $matrixArrayList] == 2} {
	    file copy -force $matrixFileList $tmpDerivFile
	} else {
	    array set matrixArray $matrixArrayList
	    set matrixFileList ""
	    set existingIndexList [array names matrixArray]
	    foreach arrayIndex $arrayIndexList {
		if {[lsearch -exact $existingIndexList $arrayIndex] != -1} {
		    lappend matrixFileList $matrixArray($arrayIndex)
		}
	    }
	    OutputStatusMessage "List of combined files: $matrixFileList" -logFileOnly 1
	    if [catch {eval exec nice -19 sddsxref $matrixFileList $tmpDerivFile \
			   "-take=* -leave=Rootname" -nowarning} result] {
		return -code error "Error during final combining: $result"
	    }
	}
	exec sddsprocess $tmpDerivFile $rmDerivFile -print=para,MatrixCode,$uniqueMatrixCode \
	    -print=para,ElementFile,[file tail $rmDerivElementFile] -def=para,CoupledMatrixAnalysis,$coupledMatrix
    }
    set deleteFiles [concat $deleteFiles $matrixFileList]
    eval file delete $deleteFiles
    eval OutputStatusMessage \"The response matrix derivative is done.\" $verboseString
    return 0
}

#-------------------------------------------------------------------------------------------------------------------------

proc MakeRMDerivElementListFile {args} {
    global tmpDir options variableTypeList rmRootList fullOptionsList
    APSParseArguments {outputFile varList rmList uniqueMatrixCode coupledMatrix rmDerivFile}
    set fullOptionsList [list useQuads useCorrs useCorrs useBPMs useBPMs useEnergy useSkew useCorrsTilt useCorrsTilt \
			     useBPMsTilt useBPMsTilt]
    set varRootList [list Quad HCorr VCorr HBPM VBPM Energy Skew HCorrTilt VCorrTilt HBPMTilt VBPMTilt]
    set rmRootList  [list HCorr VCorr HBPM VBPM]
    set tmpRoot $tmpDir/[APSTmpString]-rmderivElements
    set fileList ""
    foreach rootname $variableTypeList option $fullOptionsList elementList $varList varRoot $varRootList {
	if $options($option) {
	    exec sddsmakedataset $tmpRoot.VAR.$rootname -para=VarType,type=string -data=Var \
		-para=VarName,type=string -data=$rootname \
		-para=VarRootname,type=string -data=$varRoot \
		-col=ElementName,type=string -data=[join $elementList ,]
	    lappend fileList $tmpRoot.VAR.$rootname
	} else {
	    exec sddsmakedataset $tmpRoot.VAR.$rootname -para=VarType,type=string -data=Var \
		-para=VarName,type=string -data=$rootname \
		-para=VarRootname,type=string -data=$varRoot \
		-col=ElementName,type=string -data
	    lappend fileList $tmpRoot.VAR.$rootname
	}
    }
    foreach rootname $rmRootList elementList $rmList {
	exec sddsmakedataset $tmpRoot.RM.$rootname -para=VarType,type=string -data=RM \
	    -para=VarName,type=string -data=$rootname \
	    -para=VarRootname,type=string -data=$rootname \
	    -col=ElementName,type=string -data=[join $elementList ,]
	lappend fileList $tmpRoot.RM.$rootname
    }
    eval exec sddscombine $fileList -pipe=out \
	| sddsprocess -pipe=in $outputFile -print=para,MatrixCode,$uniqueMatrixCode -nowarning \
	"\"-print=para,options,[array get options]\"" -def=para,CoupledMatrixAnalysis,$coupledMatrix \
	-print=para,RMDFile,[file tail $rmDerivFile]
    eval file delete $fileList
}

#-------------------------------------------------------------------------------------------------------------------------

proc ChangeColumnNames {args} {
    APSParseArguments {matrixFile colPrefix colSuffix}
    exec sddsconvert $matrixFile -pipe=out "-editName=col,*,a i%$colPrefix% e i%$colSuffix%" \
	| sddsconvert -pipe=in $matrixFile.edit -rename=col,${colPrefix}Rootname${colSuffix}=Rootname
    file copy -force $matrixFile.edit $matrixFile
    file delete $matrixFile.edit
}

#-------------------------------------------------------------------------------------------------------------------------

proc RecombineResponseMatrixDeriv {args} {
    global tmpDir env locoBinDir
    APSParseArguments {input output constraintWeight constraintFile variableFile coupledMatrix useQuadConstraints}
    OutputStatusMessage "[exec date +%D--%T]: Recombining RMD with quad constraints using source file $input..."
    set tmpRoot $tmpDir/[APSTmpString]-recomb
    set deleteFiles ""
    set varNumberList [join [exec sdds2stream $variableFile.sdds -rows=bare] ]
    set quadNumber [lindex $varNumberList 0]
    set skewNumber [lindex $varNumberList 6]
    if {$coupledMatrix == 2 && $skewNumber && $useQuadConstraints == 2} {set doSkew 1} else {set doSkew 0}
    if $doSkew {
	exec sddsprocess $variableFile.Skew $tmpRoot.skew "-reedit=col,TagName,i/a/"
	exec sddscombine $variableFile.Quad $tmpRoot.skew $tmpRoot.QuadSkew -merge -overWrite
	lappend deleteFiles $tmpRoot.skew
	set quadVariableFile $tmpRoot.QuadSkew
	lappend deleteFiles $tmpRoot.QuadSkew
	#------ Calculate number of columns before skew columns:
	set addColumnsBefore [expr [lindex $varNumberList 0] + [lindex $varNumberList 1] + [lindex $varNumberList 2] \
				  + [lindex $varNumberList 3] + [lindex $varNumberList 4] + [lindex $varNumberList 5]] 
    } else {
	set quadVariableFile $variableFile.Quad
    }

    set doConstrVector 1
    if [file exists $constraintFile] {
	exec sddsconvert $quadVariableFile -pipe=out -retain=col,TagName \
	    | sddsxref -pipe $constraintFile -take=CW -match=TagName -fillIn -nowarning \
	    | sddsprocess -pipe=in $tmpRoot.constrVector \
	    "-redef=col,CW,CW 0 == ? $constraintWeight : $constraintWeight CW * \$" 
	set doConstrVector 0
    }
    if $doConstrVector {
	exec sddsconvert $quadVariableFile -pipe=out -retain=col,TagName \
	    | sddsprocess -pipe=in $tmpRoot.constrVector -def=col,CW,$constraintWeight
    }
    lappend deleteFiles $tmpRoot.constrVector
    if $doSkew {
	if [catch {exec sddsrespmatrixderivative $input $tmpRoot.rmderiv1 -mode=quad -addRowsBefore=0 -addRowsAfter=$quadNumber} result] {
	    return -code error "Error calculating constraint rows (1): $result"
	}
	if [catch {exec sddsrespmatrixderivative $input $tmpRoot.rmderiv2 -mode=quad -addRowsBefore=$addColumnsBefore -addRowsAfter=$skewNumber} result] {
	    return -code error "Error calculating constraint rows (2): $result"
	}
	if [catch {exec sddscombine $tmpRoot.rmderiv1 $tmpRoot.rmderiv2 -pipe=out -merge \
		       | sddsxref -pipe $tmpRoot.constrVector -take=CW,TagName -rename=col,TagName=Rootname1 \
		       | sddsprocess -pipe "-redef=col,%s,%s CW *,select=*,exclude=Rootname*" \
		       | sddsconvert -pipe -del=col,CW,Rootname \
		       | sddsconvert -pipe=in $tmpRoot.constraints -rename=col,Rootname1=Rootname} result] {
	    return -code error "Error calculating constraint rows: $result"
	}
	lappend deleteFiles $tmpRoot.rmderiv1 $tmpRoot.rmderiv2
    } else {
	if [catch {exec sddsrespmatrixderivative $input -pipe=out -mode=quad -addRowsBefore=0 -addRowsAfter=$quadNumber \
		       | sddsxref -pipe $tmpRoot.constrVector -take=CW,TagName -rename=col,TagName=Rootname1 \
		       | sddsprocess -pipe "-redef=col,%s,%s CW *,select=*,exclude=Rootname*" \
		       | sddsconvert -pipe -del=col,CW,Rootname \
		       | sddsconvert -pipe=in $tmpRoot.constraints -rename=col,Rootname1=Rootname} result] {
	    return -code error "Error calculating constraint rows: $result"
	}
    }
    lappend deleteFiles $tmpRoot.constraints
    exec sddscombine $input $tmpRoot.constraints -pipe=out -merge \
	| sddsprocess -pipe=in $output -print=para,SourceFile,$input
    eval file delete $deleteFiles
}

#-------------------------------------------------------------------------------------------------------------------------

proc CalculateOrbitRMD { args } {

    global locoBinDir workDir definitionFile closedOrbit env verbose
    set waitForCompletion 0

    set verboseLocal 1
    APSParseArguments {varList splitTasks scriptName scriptParameters matrixFile queueDirName useQsub qsubCommand \
				 continuePrevious rootTaskName usePopupWindow waitTime waitInterval colPrefix colSuffix \
				 submissionPause waitForCompletion remotePath doneFile matrixAbortFile optionsFile \
				 errorFile lteFile beamlineName solutionParamFileList remoteQueue \
				 abortFile qsubRespProcCommand queueSystemName useKickFiles kickFileX kickFileY verboseLocal}

    if {!$remoteQueue || !$useQsub} {
	#------ Running jobs on local computer (same computer as where the GUI runs)
	if ![file exists $queueDirName] {exec mkdir $queueDirName}
	set varListFile $queueDirName/$rootTaskName.varList.sdds
	if $useQsub {
	    if $verboseLocal {OutputStatusMessage "Starting background job..."}
	} else {
	    #------ The qsubCommand in the inititalDefinitions file could have "&" at the end. We need to redefine:
	    set qsubCommand "exec \$scriptFile"
	}
	set command "exec $locoBinDir/locoRemoteCalculateDerivative -mode localRun -workDir $workDir \
            -varList \"$varList\" -splitTasks $splitTasks -scriptName $scriptName \
             -scriptParameters \"[ProtectQuoteSigns $scriptParameters]\" \
            -matrixFile $matrixFile -queueDirName $queueDirName -rootTaskName $rootTaskName \
            -useQsub $useQsub -continuePrevious $continuePrevious -waitTime $waitTime -queueSystemName $queueSystemName \
            -qsubCommand \"[ProtectDollarSigns $qsubCommand]\" -qsubRespProcCommand \"[ProtectDollarSigns $qsubRespProcCommand]\" \
            -waitInterval $waitInterval -submissionPause $submissionPause \
            -colPrefix \"$colPrefix\" -colSuffix \"$colSuffix\" -doneFile $doneFile -errorFile $errorFile -verbose $verbose"
	OutputStatusMessage $command -logFileOnly 1
	if $useQsub {eval $command &} else {eval $command}
	#------ If you submit 2 jobs (quad and skew), they might get errors without pause
	after 1000
    } else {
	#------ Running jobs on remote host (different from computer running GUI)
	set remoteWorkDir [string range $remotePath [expr [string first : $remotePath] + 1] [string length $remotePath]]
	set remoteHostname [string range $remotePath 0 [expr [string first : $remotePath] - 1]]
	set copyFiles ""

	set parameterFileList ""
	if [llength $solutionParamFileList] { 
	    foreach paramFile $solutionParamFileList {
		lappend parameterFileList $remoteWorkDir/[file tail $paramFile]
		lappend copyFiles $paramFile
	    }
	}

	#------ Reading variables in scriptOptions and change it:
	array set optionsArray $scriptParameters
	set optionsArray(-latticeParamFileList) $parameterFileList
	set changeNameList [list lteFile twissRefFile computedRMX computedRMY variableFile]
	if $useKickFiles {lappend changeNameList kickFileX kickFileY}
	foreach name $changeNameList {
	    #------ Order of next three lines is important
	    set $name $optionsArray(-$name)
	    lappend copyFiles [set $name]
	    set optionsArray(-$name) $remoteWorkDir/[file tail $optionsArray(-$name)]
	}
	switch -exact $remoteHostname {
	    weed {set optionsArray(-localTmpDir) /tmp/$env(USER)/[file tail $optionsArray(-localTmpDir)]}
	    orthros { }
	    default {return -code error "Remote host $remoteHostname is not supported."}
	}

	catch {file delete $doneFile $matrixAbortFile $optionsFile}
	
	set optionsList [list splitTasks useQsub continuePrevious rootTaskName usePopupWindow waitTime waitInterval submissionPause \
			    queueSystemName]
	if [catch {open $optionsFile w} fid] { 
	    return -code error "$fid" 
	}
	foreach option $optionsList {
	    puts $fid "set $option [set $option]"
	}
	puts $fid "set queueDirName $remoteWorkDir"
	puts $fid "set scriptName $remoteWorkDir/[file tail $scriptName]"
	puts $fid "set matrixFile $remoteWorkDir/[file tail $matrixFile]"
	puts $fid "set computedRMX $remoteWorkDir/[file tail $computedRMX]"
	puts $fid "set computedRMY $remoteWorkDir/[file tail $computedRMY]"
	puts $fid "set varList \" $varList \""
	puts $fid "set scriptParameters \" [array get optionsArray] \""
	#--- Insert "-v LOCO_BINDIR_OUTSIDE=$remoteWorkDir" into qsubCommand:
	set qsubCommand1 [exec editstring -stream $qsubCommand "-edit=%@-V@-V -v LOCO_BINDIR_OUTSIDE=$remoteWorkDir@"]
	puts $fid "set qsubCommand \" [ProtectDollarSigns $qsubCommand1] \""
	puts $fid "set qsubRespProcCommand \" [ProtectDollarSigns $qsubRespProcCommand] \""
	#------ Have to do prefix-suffix inside remote program because the waiting in the main program is after the 
	#------ entire loop is done.
	if [string length $colPrefix] { puts $fid "set colPrefix $colPrefix" } else { puts $fid "set colPrefix \"\"" }
	if [string length $colSuffix] { puts $fid "set colPrefix $colSuffix" } else { puts $fid "set colSuffix \"\"" }
	close $fid
	
	lappend copyFiles $locoBinDir/locoRemoteCalculateDerivative $locoBinDir/locoCalculateResponseMatrix \
	    $scriptName $definitionFile $optionsFile $locoBinDir/locoCommonProcedures

	OutputStatusMessage "Starting remote job..."
	set command "exec $locoBinDir/locoRemoteCalculateDerivative -workDir $workDir -remotePath $remotePath -optionsFile $optionsFile \
	    -mode copy -matrixFile $matrixFile -copyFiles \"$copyFiles\" -doneFile $doneFile -abortFile $matrixAbortFile \
	    -errorFile $errorFile -verbose $verbose"
	OutputStatusMessage $command -logFileOnly 1
	eval $command &
	#------ If you submit 2 jobs (quad and skew), they might get errors without pause
	after 1000
    }
    #------ Waiting for the doneFile if required...
    if $waitForCompletion {
	while {[file exists $doneFile] != 1} {
	    after 1000
	    update
	    if [file exists $abortFile] {
		exec echo ABORT > $matrixAbortFile
		return -code error "Interrupted by user."
	    }
	}
    }
}

#-------------------------------------------------------------------------------------------------------------------------
proc ProtectDollarSigns {string} {
    set newString ""
    set stringLength [string length $string]
    for {set i 0} {$i < $stringLength} {incr i} {
	if {[string compare [string range $string $i $i] "$"] == 0} {
	    append newString "\\\$"
	} else {
	    append newString [string range $string $i $i]
	}
    }
    return $newString
}
#-------------------------------------------------------------------------------------------------------------------------
proc ProtectQuoteSigns {string} {
    set newString ""
    set stringLength [string length $string]
    for {set i 0} {$i < $stringLength} {incr i} {
	if {[string compare [string range $string $i $i] "\""] == 0} {
	    append newString "\\\""
	} else {
	    append newString [string range $string $i $i]
	}
    }
    return $newString
}
#-------------------------------------------------------------------------------------------------------------------------
proc AddTwoTuneLines {args} {
    global tmpDir
    APSParseArguments {matrixFile}
    set tmpRoot $tmpDir/[APSTmpString]-add2lines
    exec sddsmakedataset $tmpRoot.rootname -col=Rootname,type=string -data=nux,nuy
    exec sddsprocess $matrixFile -pipe=out -clip=2,0,inv "-redef=col,%s,0,select=*,exclude=Rootname" \
	| sddsxref -pipe=in $tmpRoot.rootname $tmpRoot.2 -replace=col,Rootname -leave=*
    exec sddscombine $matrixFile $tmpRoot.2 $tmpRoot -merge -overwrite
    file copy -force $tmpRoot $matrixFile
    file delete $tmpRoot.rootname $tmpRoot.2 $tmpRoot
}
#-------------------------------------------------------------------------------------------------------------------------
# Corrector derivative for corrector NAME has only one nonzero column that is equal to the NAME column of the response matrix
# BPM gain derivative for BPM NAME has only one nonzero row that is equal to the NAME row of the response matrix (part X or Y)
# BPM tilt/coef derivative for BPM NAME has only 2 nonzero rows (one in X, one in Y) that are equal to the NAME row of the 
# response matrix but where X and Y parts of the matrix are switched.
# First, hrmFile is transformed into a column, then vrmFile is transformed into a column, then H column is put on top of the V column
proc QuickRMDeriv {args} {
    global tmpDir allElementsFile options
    APSParseArguments {hrmFile vrmFile mode outputFile hrmFile00}
    set tmpRoot $tmpDir/[APSTmpString]-QRMD
    set localRMFile $tmpRoot.localrm
    if {[string first "Corr" $mode] == -1} {set elementType bpm} else {set elementType corr}

    set suffix ""
    switch -exact $mode {
	hCorrTilt {set suffix "#tilt"}
	vCorrTilt {set suffix "#tilt"}
	bpmTilt   {set suffix "#tilt"}
	bpmCoef   {set suffix "#coef"}
    }
    set elementList [join [exec sdds2stream $allElementsFile.$mode -col=TagName] ]

    exec sddscombine $allElementsFile.xBpm $allElementsFile.yBpm $tmpRoot.bpms -overwrite
    foreach rmFile [list $hrmFile $vrmFile] outFile [list $tmpRoot.out.h $tmpRoot.out.v] {
	set fileList ""
	file copy -force $rmFile $localRMFile
	set rowNamecolumn TagName
	#--- If mode=tilt or coef, flip the response matrix pages
	if {[string compare $mode bpmTilt] == 0 || [string compare $mode bpmCoef] == 0} {
	    exec sddsprocess $rmFile -pipe=out -redef=para,PageNumber,i_page "-reedit=col,TagName,%@#X@@ %@#Y@@" \
		| sddssort -pipe -para=PageNumber,decreasing \
		| sddsxref -pipe $rmFile -take=TagName -rename=col,TagName=TagNamePlane -nowarning \
		| sddsselect -pipe=in $hrmFile00 $localRMFile -match=TagNamePlane=TagName
	    set rowNamecolumn TagNamePlane
	    if {[string compare $mode bpmTilt] == 0} {
		exec sddsprocess $localRMFile "-redef=col,%s,PageNumber 1 == ? %s -1 * : %s \$,select=*,exclude=*Name*" -nowarning
		file delete ${localRMFile}~
	    }
	}

	#--- Go element by element
	foreach element $elementList {
	    if {[string compare $elementType "bpm"] == 0} {
		#--- Matrix with non-zero row at the element position for BPMs
		exec sddsprocess $localRMFile -pipe=out "-redef=col,%s,TagName \"$element\" streq ? %s : 0 \$,select=*,exclude=*Name*" \
		    | sddscombine -pipe -merge \
		    | sddsmatrix2column -pipe=in $tmpRoot.$element$suffix -dataColumn=$element$suffix -major=row -rowNameColumn=$rowNamecolumn
	    } else {
		#--- Matrix with non-zero column at the element position for correctors
		exec sddsconvert $localRMFile -pipe=out -del=col,*Name \
		    | sddsprocess -pipe "-redef=col,%s,0,select=*,exclude=$element" \
		    | sddscombine -pipe -merge \
		    | sddsmatrix2column -pipe=in $tmpRoot.$element$suffix -dataColumn=$element$suffix -major=row
	    }
	    lappend fileList $tmpRoot.$element$suffix
	}
	eval exec sddsxref $fileList $outFile -take=* -leave=Rootname
	eval file delete $fileList $tmpRoot.rootname $tmpRoot.localrm
    }

    eval exec sddscombine $tmpRoot.out.h $tmpRoot.out.v $outputFile -merge -over
    if $options(useTune) {if [catch {AddTwoTuneLines -matrixFile $outputFile} result] {return -code error "AddTwoTuneLines $result"}}
    file delete $tmpRoot.out.h $tmpRoot.out.v $tmpRoot.bpms 
}
#-------------------------------------------------------------------------------------------------------------------------
#------ If not all variables are used in the fit (but used in the RM), filter them out from the matrix file.

proc RemoveNotUsedColumns { args } {
    APSParseArguments { varList rmList matrixFile }
    if {[llength $varList] != [llength $rmList]} {
        set removeList ""
        foreach rmElement $rmList {
            if {[lsearch -exact $varList $rmElement] == -1} {
                lappend removeList $rmElement
            }
        }
        exec sddsconvert $matrixFile -nowarning -del=col,[join $removeList ,]
        file delete $matrixFile~
    }
}

#-------------------------------------------------------------------------------------------------------------------------

proc CalculateTiltedCorrResponse { args } {

    global tmpDir workDir 
    global variableFile acceleratorCode hrmIterFile vrmIterFile solutionParamFileList
    APSParseArguments {lteFile beamlineName latticeParamFileList coupledMatrix tiltedHRMFile tiltedVRMFile delta}
    
    #------ tmpRoot is local here so that nodes knew the file:
    set tmpRoot $workDir/[APSTmpString]-crossRM
    
    #------ Preparing files for calculating RM derivative with respect to corrector tilts:
    #------ We calculate response matrix where all correctors are tilted, take difference from
    #------ initial RM and then divide it by delta. This file is then used by QuickRMDeriv
    #------ to build RM derivative.

    if [catch {join [exec sdds2stream $variableFile.hCorrTilt -col=TagName] } hCorrList] {
        return -code error "Error reading $variableFile.hCorrTilt: $hCorrList"
    }
    if [catch {join [exec sdds2stream $variableFile.vCorrTilt -col=TagName] } vCorrList] {
        return -code error "Error reading $variableFile.vCorrTilt: $vCorrList"
    }
    #--- Make parameter file where all correctors are tilted:
    set tiltParamFile $tmpRoot.tiltParam
    lappend deleteFiles $tiltParamFile
    exec sddsmakedataset -pipe=out -col=TagName,type=string -data=[join [concat $hCorrList $vCorrList] ,] \
	| sddsprocess -pipe=in $tiltParamFile -print=col,ElementParameter,TILT -def=col,ParameterValue,$delta \
	-print=col,ParameterMode,absolute
    if [catch {Fit_MakeElementNameFromTagName -input $tiltParamFile} result] {
	return -code error "Fit_MakeElementNameFromTagName: $result"
    }
    set localParamFileList [concat $solutionParamFileList $tiltParamFile]
    set latticeOptions "-lteFile $lteFile -beamlineName $beamlineName -latticeParamFileList \"$localParamFileList\""

    if [catch {eval CalculateResponseMatrix $latticeOptions \
		   {-acceleratorCode $acceleratorCode \
			-hrmFile $tmpRoot.tiltCorrXRM -vrmFile $tmpRoot.tiltCorrYRM -hCorrList $hCorrList -vCorrList $vCorrList \
			-workDir $workDir -twiFile $tmpRoot.twi -coupledMatrix $coupledMatrix}} result] {
        return -code error "CalculateResponseMatrix: $result"
    }
    lappend deleteFiles $tmpRoot.tiltCorrXRM $tmpRoot.tiltCorrYRM $tmpRoot.twi

    #------ Filtering out not used bpms (or other elements) and calculating (M-Mo)/delta and combining in one page for rmDeriv.c...
    foreach tiltCorrFile [list $tiltedHRMFile $tiltedVRMFile] iterFile [list $hrmIterFile $vrmIterFile] \
	ext [list tiltCorrXRM tiltCorrYRM] {
	    if [catch {exec sddsselect $tmpRoot.$ext $iterFile -pipe=out -match=TagName \
			   | sddschanges -pipe -base=$iterFile -changes=exclude=*Name,* -copy=TagName,BPMName -parallelPages \
			   | sddsprocess -pipe "-redef=col,%s,%s $delta /,select=ChangeIn*" \
			   | sddsconvert -pipe=in $tiltCorrFile "-edit=col,ChangeIn*,8d"} result] {
		return -code error "Error processing corrector cross response (files: $tmpRoot.$ext and $iterFile): $result"
	    }
	}
    eval file delete $deleteFiles

}

#-------------------------------------------------------------------------------------------------------------------------

proc FilterResponseMatrixDeriv { args } {
    global tmpDir dispColumn options
    set tmpRoot $tmpDir/[APSTmpString]-filtMatr
    APSParseArguments {inputFile outputFile elementListVar elementListRM optionsList coupledMatrix fitDispersion \
			   useTune useKickFiles kickFileX kickFileY}

    OutputStatusMessage "[exec date +%D--%T]: Filtering response matrix derivative (using source $inputFile)..."
    if $useKickFiles {OutputStatusMessage "Warning: for useKickFiles=1, filtering is done for columns only."}

    #------ Checking that we are not overwriting inputFile (that inputFile and outputFile are different files). 
    #------ Cannot just compare filenames because of disk mounting at APS.
    if ![file exists $inputFile] {
	return -code error "Error: File $inputFile does not exist."
    }
    if [file exists $outputFile] {
	file rename -force $outputFile $outputFile.tmpcopy
	if ![file exists $inputFile] {
	    file rename -force $outputFile.tmpcopy $outputFile
	    return -code error "Error: Input and output files are the same: $inputFile --- $outputFile"
	}
	file rename -force $outputFile.tmpcopy $outputFile
    }

    #------ Inserting $dispColumn at the end of H corr list (index 1 in the big list):
    if {$fitDispersion && $options(useCorrs)} {
	set elementListVar [lreplace $elementListVar 1 1 \
			     [linsert [lindex $elementListVar 1] [llength [lindex $elementListVar 1]] $dispColumn]]
    }

    set filterColumns 0
    set filterRows 0
    set columnListOld [exec sddsquery -col $inputFile]
    if [catch {BuildColumnList -elementList $elementListVar -coupledMatrix $coupledMatrix} columnListNew] {
	return -code error "BuildColumnList: $columnListNew"
    }
    #------ Checking for elements that are not in original derivative matrix:
    set missingElements ""
    set error 0
    foreach column $columnListNew {
	if {[lsearch -exact $columnListOld $column] == -1} {
	    set error 1
	    lappend missingElements $column
	}
    }
    if $error {return -code error "Error: Element that is not in original matrix is used for filtering: $missingElements."}
    if {[llength $columnListNew] < [llength $columnListOld]} {
	set filterColumns 1
	set removeColumnList ""
	foreach column $columnListOld {
	    if {[lsearch -exact $columnListNew $column] == -1} {
		lappend removeColumnList $column
	    }
	}
	OutputStatusMessage "Filtering: number of columns to remove: [llength $removeColumnList]" -logFileOnly 1
    }

    #------ Building lists of new and old rows:
    set rowListOld [exec sdds2stream $inputFile -col=Rootname]
    if [catch {BuildRowList -elementList $elementListRM -coupledMatrix $coupledMatrix \
		   -fitDispersion $fitDispersion -useTune $useTune} rowListNew] {
	return -code error "BuildRowList: $rowListNew"
    }

    #------ Building file with rows to be removed in the Rootname column:
    set absentRows ""
    if {[llength $rowListNew] < [llength $rowListOld]} {
	set filterRows 1
	foreach element $rowListOld {
	    if {[lsearch -exact $rowListNew $element] == -1} { lappend absentRows $element }
	}
	set rowLimit 6000
	if {[llength $absentRows] > $rowLimit} {
	    set splits [expr [llength $absentRows] / $rowLimit + 1]
	    set newList [Fit_SplitList -List $absentRows -splitTasks $splits]
	    set fileList ""; set counter 0
	    foreach smallList $newList {
		if [catch {exec sddsmakedataset $tmpRoot.$counter -col=Rootname,type=string -data=[join $smallList ,]} result] {
		    return -code error "Error running sddsmakedataset (1) (list length [llength $smallList]): $result"
		}
		lappend fileList $tmpRoot.$counter
		incr counter
	    }
	    eval exec sddscombine $fileList $tmpRoot -merge
	    eval file delete $fileList
	} else {
	    if [catch {exec sddsmakedataset $tmpRoot -col=Rootname,type=string -data=[join $absentRows ,]} result] {
		return -code error "Error running sddsmakedataset (2) (list length [llength $absentRows]): $result"
	    }
	}
    }
    OutputStatusMessage "Filtering: Input file is $inputFile; output file is $outputFile. Old rows: [llength $rowListOld]; new rows:[llength $rowListNew]; absent rows: [llength $absentRows]." -logFileOnly 1

    set tmpOutputFile $tmpDir/[file tail $outputFile]
    file delete $tmpOutputFile
    #------ Do not filter rows in useKickFiles=1 (don't have time to do it).
    if $useKickFiles {set filterRows 0}
    #------ Applying filtering commands:
    if {$filterColumns && !$filterRows} {
	exec sddsconvert $inputFile $tmpOutputFile -del=col,[join $removeColumnList ,]
    }
    if {$filterColumns && $filterRows} {
	exec sddsconvert $inputFile -pipe=out -del=col,[join $removeColumnList ,] \
	    | sddsselect -pipe=in $tmpRoot $tmpOutputFile -match=Rootname -invert
    }
    if {!$filterColumns && $filterRows} {
	exec sddsselect $inputFile $tmpRoot $tmpOutputFile -match=Rootname -invert
    }
    file delete $tmpRoot
    if {$filterColumns || $filterRows} {
	OutputStatusMessage "[exec date +%D--%T]: Filtering is done."
    } else {
	OutputStatusMessage "[exec date +%D--%T]: $inputFile has no columns to filter. Copying..."
	file copy -force $inputFile $tmpOutputFile
    }
    file copy -force $tmpOutputFile $outputFile
    file delete $tmpOutputFile
}

#-------------------------------------------------------------------------------------------------------------------------

proc BuildColumnList { args } {
    APSParseArguments { elementList coupledMatrix}
    #------ Adding X or Y to bpm names (only when coupledMatrix is 0 or 1):
    switch -exact -- $coupledMatrix {
	0 {
	    foreach index [list 3 4] suffix [list X Y] {
		set elementList [lreplace $elementList $index $index \
				     [EditListElements -elementList [lindex $elementList $index] -suffix $suffix]]
	    }
	}
	1 {
	    foreach bpmIndex [list 3 4] plane [list Y X] {
		set elementList [lreplace $elementList $index $index \
				     [EditListElements -elementList [lindex $elementList $index] -suffix $suffix]]
	    }
	}
	2 {
	    #------ Adding prefix "a" to corrector and bpm:
	    foreach index [list 7 8 9 10] suffix [list #tilt #tilt #tilt #coef] {
		set elementList [lreplace $elementList $index $index \
				     [EditListElements -elementList [lindex $elementList $index] -suffix $suffix]]
	    }
	}
    }
    #------ Adding E to corrector names for energy fitting:
    set index 5
    set suffix E
    set elementList [lreplace $elementList $index $index \
			 [EditListElements -elementList [lindex $elementList $index] -suffix $suffix]]
    #------ Adding a before skew quad names:
    set index 6
    set prefix a
    set elementList [lreplace $elementList $index $index \
			 [EditListElements -elementList [lindex $elementList $index] -prefix $prefix]]
    return [concat Rootname [join $elementList ]]
}

#-------------------------------------------------------------------------------------------------------------------------

proc EditListElements { args } {
    set prefix ""
    set suffix ""
    APSParseArguments {elementList prefix suffix}
    set localList ""
    foreach element $elementList {
	lappend localList ${prefix}${element}${suffix}
    }
    return $localList
}

#-------------------------------------------------------------------------------------------------------------------------

proc BuildRowList { args } {
    global dispColumn
    APSParseArguments {elementList coupledMatrix fitDispersion useTune}
    #------ Inserting $dispColumn at the end of the H corr list (position 0 in the big list):
    if $fitDispersion {
	set elementList [lreplace $elementList 0 0 \
			     [linsert [lindex $elementList 0] [llength [lindex $elementList 0]] $dispColumn]]
    }
    switch -exact -- $coupledMatrix {
	0 {
	    set rowList ""
	    foreach bpmIndex [list 2 3] corrIndex [list 0 1] {
		foreach bpm [lindex $elementList $bpmIndex] {
		    foreach corr [lindex $elementList $corrIndex] {
			lappend rowList $bpm$corr
		    }
		}
	    }
	}
	1 {
	    set rowList ""
	    foreach bpmIndex [list 3 2] corrIndex [list 0 1] {
		foreach bpm [lindex $elementList $bpmIndex] {
		    foreach corr [lindex $elementList $corrIndex] {
			lappend rowList $bpm$corr
		    }
		}
	    }
	}
	2 {
	    set rowList ""
	    foreach corrIndex [list 0 1] {
		foreach bpm [concat [lindex $elementList 2] [lindex $elementList 3]] {
		    foreach corr [lindex $elementList $corrIndex] {
			lappend rowList $bpm$corr
		    }
		}
	    }
	}
    }
    if $useTune {
	set rowList [concat $rowList Snux Snuy]
    }
    return $rowList
}

#-------------------------------------------------------------------------------------------------------------------------

proc CalculateRMDerivInverse { args } {

    global tmpDir

    set verbose 1
    APSParseArguments { SVnumber rmDerivFile rmDerivInverseFile rmDifferenceFile programName verbose}

    if $verbose {set verboseString ""} else {set verboseString " -logFileOnly 1 "}
    eval OutputStatusMessage \"[exec date +%D--%T]: Building the inverse matrix using $programName...\" $verboseString
    catch {eval file delete [glob ${rmDerivInverseFile}*]}

    global svdProgramFullName
    if [info exists svdProgramFullName] {
	set executable $svdProgramFullName
    } else {
	set executable sddspseudoinverse
    }
    if [catch {GetNumberOfCores} cores] {OutputStatusMessage "Warning: Error reading number of cores"; set cores 1}
    OutputStatusMessage "Using the following program for inverse: $executable, with $cores threads." -logFileOnly 1
    eval OutputStatusMessage \"Using $cores threads.\" $verboseString
    set tikhonov 1
    if $tikhonov {
	eval OutputStatusMessage \"With Tikhonov regularization.\" $verboseString
	if [catch {exec $executable $rmDerivFile -pipe=out -tikhonov=svn=$SVnumber \
		       -sFile=$rmDerivFile.SV -economy -thread=$cores \
		       | sddsconvert -pipe=in $rmDerivInverseFile -del=col,OldColumnNames -del=para,*} result] {
	    return -code error "Error doing sddspseudoinverse: $result"
	}
    } else {
	if [catch {exec $executable $rmDerivFile -pipe=out -largestSingularValues=$SVnumber \
		       -sFile=$rmDerivFile.SV -economy -thread=$cores \
		       | sddsconvert -pipe=in $rmDerivInverseFile -del=col,OldColumnNames -del=para,*} result] {
	    return -code error "Error doing sddspseudoinverse: $result"
	}
    }
    exec sddsconvert $rmDerivFile.SV -nowarning -rename=col,SingularValues=SV
    file delete $rmDerivFile.SV~
    file stat $rmDerivInverseFile statArray
    eval OutputStatusMessage \"RM Inverse calculated. The size is $statArray(size).\" $verboseString

}

#-------------------------------------------------------------------------------------------------------------------------
proc GetNumberOfCores {args} {
    set threadsPerCore [lindex [exec lscpu | fgrep Thread] end]
    set threads [exec nproc --all]
    set cores [expr $threads / $threadsPerCore]
    return $cores
}
#-------------------------------------------------------------------------------------------------------------------------

#-------------------------------------------------------------------------------------------------------------------------
#-------------------------------------------------------------------------------------------------------------------------
#------ Miscelaneous procedures used during iterations like output, model update, etc. --------------------------------------
#-------------------------------------------------------------------------------------------------------------------------
#-------------------------------------------------------------------------------------------------------------------------

proc IterationResultsOutput {args} {

    global useQuadConstraints
    set vectorStdDev ""
    APSParseArguments {rmsList fitDispersion vectorStdDev coupledMatrix iteration sddsFile bigIter}

    set outputLine ""
    set lineVectorStdDev ""
    array set rmsArray $rmsList
    array set sddsArray [list iteration $iteration stDev 0.0 rmsTotal 0.0 rmsTotalWQ 0.0 rmsDX 0.0 rmsDY 0.0 \
			     rmsXX 0.0 rmsXY 0.0 rmsYX 0.0 rmsYY 0.0 rmsQ 0.0]

    foreach index [array names rmsArray] {set sddsArray($index) $rmsArray($index)}

    if [string length $vectorStdDev] {
        set lineVectorStdDev [format "%s %9.3e" "Var. stdev:" $vectorStdDev]
	if $useQuadConstraints {append lineVectorStdDev [format "%s %9.3e " ". Total res. error: " $rmsArray(rmsTotal)]}
        append lineVectorStdDev "\n"
	set sddsArray(stDev) $vectorStdDev
    }

    if {$coupledMatrix == 2} {
        set lineX [format "%9.3e,%9.3e" $rmsArray(rmsXX) $rmsArray(rmsXY)]
        set lineY [format "%9.3e,%9.3e" $rmsArray(rmsYX) $rmsArray(rmsYY)]
    } else {
        set lineX [format "%9.3e" $rmsArray(rmsXX)]
        set lineY [format "%9.3e" $rmsArray(rmsYY)]
    }
    if $fitDispersion {
        if {$coupledMatrix == 2} {
            set lineDisp "D(m): [format "%9.3e,%9.3e" $rmsArray(rmsDX) $rmsArray(rmsDY)])"
        } else {
            set lineDisp "D(m): [format "%9.3e" $rmsArray(rmsDX)])"
        }
    } else {
        set lineDisp ")"
    }
    if $useQuadConstraints {
	set lineQuadConstraints " Quads: [format "%9.3e" $rmsArray(rmsQ)]"
	set totalRmsOutput $rmsArray(rmsTotalWQ)
    } else {
	set lineQuadConstraints ""
	set totalRmsOutput $rmsArray(rmsTotal)
    }
    set lineRes [format "%s %9.3e %s %s %s %s " "Res. error: " $totalRmsOutput " (X(mm):" $lineX "Y(mm):" $lineY]
    append outputLine $lineVectorStdDev $lineRes $lineDisp $lineQuadConstraints
    
    OutputStatusMessage $outputLine

    if [catch {OutputSddsFile -sddsList [array get sddsArray] -outputFile $sddsFile -bigIter $bigIter} result] {
	return -code error "OutputSddsFile: $result"
    }
}

#-------------------------------------------------------------------------------------------------------------------------

proc OutputSddsFile {args} {
    APSParseArguments {sddsList outputFile bigIter}
    array set sddsArray $sddsList
    set deleteFiles ""

    #--- Read overall iteration number:
    if ![file exists $outputFile] {
	set outputIter 0
    } else {
	if {$bigIter == 0} {
	    set outputIter $sddsArray(iteration)
	} else {
	    if [catch {exec sddscombine $outputFile -pipe=out -merge \
			   | sddsprocess -pipe -filter=col,BigIteration,0,[expr $bigIter - 1] -clip=0,1,inv \
			   | sdds2stream -pipe=out -col=Iteration} iterations0] {
		return -code error "Error reading iterations0: $iterations0"
	    }
	    set outputIter [expr $sddsArray(iteration) + $iterations0]
	}
    }

    #--- Make one-line log file:
    exec sddsmakedataset $outputFile.tmp \
	-col=Iteration,type=long -data=$outputIter \
	-col=BigIteration,type=long -data=$bigIter \
	-col=StDev,type=double -data=$sddsArray(stDev) \
	-col=TotalRms,type=double -data=$sddsArray(rmsTotal) \
	-col=TotalRmsWQ,type=double -data=$sddsArray(rmsTotalWQ) \
	-col=RmsXX,type=double -data=$sddsArray(rmsXX) \
	-col=RmsXY,type=double -data=$sddsArray(rmsXY) \
	-col=RmsYX,type=double -data=$sddsArray(rmsYX) \
	-col=RmsYY,type=double -data=$sddsArray(rmsYY) \
	-col=RmsDx,type=double -data=$sddsArray(rmsDX) \
	-col=RmsDy,type=double -data=$sddsArray(rmsDY) \
	-col=RmsQ,type=double -data=$sddsArray(rmsQ)
    lappend deleteFiles $outputFile.tmp

    #--- Make one- or multi-page file:
    if ![file exists $outputFile] {
	file copy -force $outputFile.tmp $outputFile
    } else {
	exec sddscombine $outputFile $outputFile.tmp -pipe=out -merge \
	    | sddsbreak -pipe=in $outputFile.tmp1 -change=BigIteration
	file copy -force $outputFile.tmp1 $outputFile
	lappend deleteFiles $outputFile.tmp1
    }

    eval file delete $deleteFiles
}

#-------------------------------------------------------------------------------------------------------------------------

proc UpdateModel {args} {
    
    APSParseArguments { elementsUpdateFile vectorOutFile corFraction variableFile allElementsFile }
    
    global workDir tmpDir options dispColumn acceleratorCode keepTunesOnTarget tuneCorrectionFile
    global indexList variableTypeList adjustAverageGains specialElementsInputFile specialParamFile

    set tmpRoot $tmpDir/[APSTmpString]-updateModel
    set deleteFiles ""
    #------ Preparing update file...
    if [catch { exec sddsxref $elementsUpdateFile $vectorOutFile -pipe=out -take=MatrixDifference \
                    | sddsconvert -pipe -del=para,PageID \
                    | sddsbreak -pipe -changeof=VariableType \
                    | sddsprocess -pipe -process=VariableType,first,PageID \
                    | sddsconvert -pipe=in $tmpRoot.Update -del=col,ParameterValue \
                    -rename=col,MatrixDifference=VarAdjustment} result] {
        return -code error "Error creating update file ($elementsUpdateFile; $vectorOutFile): $result"
    }
    set existingVarList [join [exec sdds2stream $tmpRoot.Update -para=PageID] ]
    #------ Adjusting average gains...
    if {[info exists adjustAverageGains] && $options(useBPMs)} {
        if [catch {AdjustAverageGains -updateFile $tmpRoot.Update -adjustAverageGains $adjustAverageGains \
                       -dispColumn $dispColumn -dataColumn VarAdjustment} result] {
            return -code error "AdjustAverageGains: $result"
        }
    }
    #------ Update variable and allElements files...
    if {![string compare $acceleratorCode elegant] && $options(useSkew)} {
	#----- Remove K1 rows from elements.Skew before updating...
	exec sddsprocess $allElementsFile.Skew -match=col,ElementParameter=*TILT -nowarning
	lappend deleteFiles $allElementsFile.Skew~
    }
    set variableFileList ""
    set allElementsFileList ""
    foreach index [linsert $indexList 100 special] variableType [linsert $variableTypeList 100 special] {
        if {[lsearch -exact $existingVarList $variableType] != -1} {
            exec sddsprocess $tmpRoot.Update $tmpRoot.$variableType -nowarning \
                -match=para,PageID=$variableType
            if ![exec sdds2stream $tmpRoot.$variableType -rows=bare] {
                OutputStatusMessage "Warning: $index option is on, but $elementsUpdateFile file has zero rows of this \
                                                    variable type."
            }
            foreach rootFile [list $variableFile $allElementsFile] {
                if [catch {exec sddsxref $rootFile.$variableType $tmpRoot.$variableType -pipe=out -nowarning \
                               -take=VarAdjustment -match=TagName -fillIn \
                               | sddsprocess -pipe \
                               "-redefine=col,ParameterValue,ParameterValue VarAdjustment $corFraction * + " \
                               | sddsconvert -pipe=in $tmpRoot.$variableType.new -del=col,VarAdjustment } result] {
                    return -code error "Error updating allElements file for $variableType type: $result"
                }
                file copy -force $tmpRoot.$variableType.new $rootFile.$variableType
                lappend deleteFiles $tmpRoot.$variableType.new
            }
        }
        lappend deleteFiles $tmpRoot.$variableType 
        lappend variableFileList $variableFile.$variableType
        lappend allElementsFileList $allElementsFile.$variableType
    }

    #------ Adjust tunes, change *.Quad files:
    if $keepTunesOnTarget {
	set applyCorr [exec sdds2stream $tuneCorrectionFile -para=ApplyCorrection]
	if $applyCorr {
	    foreach rootFile [list $variableFile $allElementsFile] {
		if [catch {exec sddsxref $rootFile.Quad $tuneCorrectionFile -pipe=out -nowarning -take=ParameterValue -match=TagName \
			       -rename=col,ParameterValue=TuneCorrection -fillIn \
			       | sddsprocess -pipe "-redefine=col,ParameterValue,ParameterValue TuneCorrection + " \
			       | sddsconvert -pipe=in $tmpRoot.tuneCorr -del=col,TuneCorrection} result] {
		    return -code error "Error updating allElements file for $variableType type: $result"
		}
                file copy -force $tmpRoot.tuneCorr $rootFile.Quad
		lappend deleteFiles $tmpRoot.tuneCorr
	    }
	}
    }

    #------ Final combining of updated files
    eval exec sddscombine $variableFileList $variableFile.sdds -overWrite
    eval exec sddscombine $allElementsFileList -pipe=out | sddsconvert -pipe=in $allElementsFile.sdds -del=col,ElementOccurence
    exec sddscombine $variableFile.sdds $elementsUpdateFile -merge -overWrite
    lappend deleteFiles ${elementsUpdateFile}~ $tmpRoot.Update
    
    if $options(useSpecElem) {
	if [catch {MakeRealParamFileFromSpecialVariables -specialElementsInputFile $specialElementsInputFile -specialElementsFile $allElementsFile.special \
		       -specialParamFile $specialParamFile} result] {
	    return -code error "MakeRealParamFileFromSpecialVariables: $result"
	}
    }
    eval file delete $deleteFiles
}

#-------------------------------------------------------------------------------------------------------------------------
#--- specialElementsInputFile -- special element input file with elements and their description, etc
#--- specialElementsFile -- elements.special for LOCO
#--- specialParamFile -- real elegant parameter file
proc MakeRealParamFileFromSpecialVariables {args} {
    APSParseArguments {specialElementsInputFile specialElementsFile specialParamFile}
    set tmpRoot [APSTmpString]-specVar
    if [exec sdds2stream $specialElementsFile -rows=bare] {
	set specialFile $tmpRoot.special
	if [catch {join [exec sddsbreak $specialElementsFile -pipe=out -changeof=ElementParameter \
			     | sddsprocess -pipe -process=ElementParameter,first,SpecialType \
			     | sdds2stream -pipe=in -para=SpecialType] } specialVarTypeList] {
	    return -code error "Error reading special var file $specialElementsFile: $specialVarTypeList"
	}

	set specialFileList ""
	foreach specialVarType $specialVarTypeList {
	    if [catch {HandleSpecialVariableType -specialElementsInputFile $specialElementsInputFile -specialElementsFile $specialElementsFile \
			   -specialVarType $specialVarType -outputFile $specialFile.$specialVarType} result] {
		return -code error "HandleSpecialVariableType: $result"
	    }
	    lappend specialFileList $specialFile.$specialVarType
	}
	if {[llength $specialFileList] == 1} {
	    file copy -force $specialFileList $specialParamFile
	} else {
	    eval exec sddscombine $specialFileList $specialParamFile -overWrite
	}
	eval file delete $specialFileList
    }
}
#-------------------------------------------------------------------------------------------------------------------------

proc AdjustAverageGains { args } {
    APSParseArguments { updateFile adjustAverageGains dispColumn dataColumn}
    set pageIDList [join [exec sdds2stream $updateFile -para=PageID] ]
    set nSigma 3
    set averageList [join [exec sddsoutlier $updateFile -pipe=out -col=$dataColumn -stDev=$nSigma \
			       | sddsprocess -pipe -process=$dataColumn,average,Aver \
			       | sdds2stream -pipe=in -para=Aver] ]

    set adjNameList [list hCorr xBpm vCorr yBpm]
    set presentVars [exec sdds2stream $updateFile -para=PageID]
    foreach adjName $adjNameList {
	if {[lsearch -exact $presentVars $adjName] == -1} {
	    puts stdout "Warning: will not do average gain adjustment"
	    return
	}
    }
    switch -exact -- $adjustAverageGains {
        dispersion {
	    if [catch {exec sddsprocess $updateFile -pipe=out -match=para,PageID=hCorr \
			   | sddsprocess -pipe -match=col,TagName=$dispColumn \
			   | sdds2stream -pipe=in -col=$dataColumn} gainAdjustmentX] {
		return -code error "Error getting dispersion gain calibration: $gainAdjustmentX"
	    }
	    if ![string length $gainAdjustmentX] {set gainAdjustmentX 0}
	    set gainAdjustmentY [lindex $averageList [lsearch -exact $pageIDList vCorr]]
	    set gainList [list [expr $gainAdjustmentX * -1] $gainAdjustmentX [expr $gainAdjustmentY * -1] $gainAdjustmentY]
	}
	averageBpm {
	    set gainAdjustmentX [lindex $averageList [lsearch -exact $pageIDList xBpm]]
	    set gainAdjustmentY [lindex $averageList [lsearch -exact $pageIDList yBpm]]
	    set gainList [list $gainAdjustmentX  [expr $gainAdjustmentX * -1] $gainAdjustmentY [expr $gainAdjustmentY * -1]]
	}
	averageCor {
	    set gainAdjustmentX [lindex $averageList [lsearch -exact $pageIDList hCorr]]
	    set gainAdjustmentY [lindex $averageList [lsearch -exact $pageIDList vCorr]]
	    set gainList [list [expr $gainAdjustmentX * -1] $gainAdjustmentX [expr $gainAdjustmentY * -1] $gainAdjustmentY]
	}
	default {
            set gainList [list 0 0 0 0]
	}
    }
    foreach pageID $pageIDList {
	set index [lsearch -exact $adjNameList $pageID]
	if {$index != -1} { lappend dataList [lindex $gainList $index] } else { lappend dataList 0 }
    }
    exec sddsmakedataset -pipe=out -col=GainAdj,type=double -data=[join $dataList ,] \
	| sddsbreak -pipe -rowLimit=1 \
	| sddsprocess -pipe=in $updateFile.gainAdj -process=GainAdj,first,GainAdjustment
    exec sddsxref $updateFile $updateFile.gainAdj -pipe=out -leave=* -transfer=para,GainAdjustment \
	| sddsprocess -pipe=in $updateFile.tmp "-redef=col,$dataColumn,$dataColumn GainAdjustment +"
    file copy -force $updateFile.tmp $updateFile
    file delete $updateFile.gainAdj $updateFile.tmp
}

#-------------------------------------------------------------------------------------------------------------------------

proc StopExecution {args} {
    
    global workDir badPointsLevel variableArchiveFile twiFile
    global elementOrderFile rmDerivInverseFile
    APSParseArguments {resultsFile allElementsFile hrmIterFile vrmIterFile rmsList badPointsIterCounter}

    array set rmsArray $rmsList
    set residualError $rmsArray(rmsTotal)
    if {[lsearch -exact [array name rmsArray] rmsTotalWQ] == -1} {
	set residualErrorWQ $residualError
    } else {
	set residualErrorWQ $rmsArray(rmsTotalWQ)
    }
    unset badPointsLevel

    #------ If iterations=0, the file does not exist:
    if [file exists $variableArchiveFile] {
	#--- Finding best iteration (Residual error is without quads in the file):
	if [catch {exec sddscollapse $variableArchiveFile -pipe=out \
                     | sddsbreak -pipe -change=Iteration \
                     | sddsprocess -pipe -clip=1,0,inv \
                     | sddscombine -pipe -merge \
                     | sddsprocess -pipe -process=ResidualError,min,BestIteration,position,functionof=Iteration \
                     | sdds2stream -pipe=in -para=BestIteration} bestIteration] {
	    return -code error "Error reading best iteration: $bestIteration"
	}
	#--- Filtering the best iteration results from the *.iterations files:
	exec sddsprocess $variableArchiveFile $resultsFile -nowarning \
          -filter=para,Iteration,$bestIteration,$bestIteration \
          -redef=para,ResidualError,$residualError -redef=para,ResidualErrorWQ,$residualErrorWQ \
          -redef=para,BadPointsIterCounter,$badPointsIterCounter
	exec sddsprocess $hrmIterFile.iterations $hrmIterFile.after -filter=para,Iteration,$bestIteration,$bestIteration
	exec sddsprocess $vrmIterFile.iterations $vrmIterFile.after -filter=para,Iteration,$bestIteration,$bestIteration
	exec sddsprocess $twiFile.iterations $twiFile.after -filter=para,Iteration,$bestIteration,$bestIteration
    } else {
	file copy -force $hrmIterFile.before $hrmIterFile.after
	file copy -force $vrmIterFile.before $vrmIterFile.after
	file copy -force $allElementsFile.sdds $resultsFile
    }

    file delete $workDir/vector_out $workDir/computed.Xmatrix $workDir/computed.Ymatrix
    file delete $elementOrderFile
    catch {eval file delete [glob -nocomplain -- $workDir/*~]}
}

#-------------------------------------------------------------------------------------------------------------------------
#-------------------------------------------------------------------------------------------------------------------------
#-------------------------------------------------------------------------------------------------------------------------
