use of org.vcell.util.ISize in project vcell by virtualcell.
the class RunRefSimulationFastOp method runRefSimulation.
private RowColumnResultSet runRefSimulation(ROI cellROI, ROI[] imageDataROIs, UShortImage psf, FloatImage initRefConc, double experimentalRecoveryTime, LocalWorkspace localWorkspace, ClientTaskStatusSupport progressListener) throws Exception {
User owner = LocalWorkspace.getDefaultOwner();
KeyValue simKey = LocalWorkspace.createNewKeyValue();
//
// save first image from normalized time series as the initial concentration field data
//
ExternalDataInfo initialConcentrationExtData = createNewExternalDataInfo(localWorkspace, INITCONC_EXTDATA_NAME);
Extent extent = initRefConc.getExtent();
Origin origin = initRefConc.getOrigin();
ISize isize = new ISize(initRefConc.getNumX(), initRefConc.getNumY(), initRefConc.getNumZ());
saveExternalData(initRefConc, INITCONC_EXTDATA_VARNAME, initialConcentrationExtData.getExternalDataIdentifier(), localWorkspace);
FieldFunctionArguments initConditionFFA = new FieldFunctionArguments(INITCONC_EXTDATA_NAME, INITCONC_EXTDATA_VARNAME, new Expression(0.0), VariableType.VOLUME);
//
// save ROIs as a multivariate field data
//
ExternalDataInfo roiExtData = createNewExternalDataInfo(localWorkspace, ROI_EXTDATA_NAME);
saveROIsAsExternalData(imageDataROIs, localWorkspace, roiExtData.getExternalDataIdentifier());
ArrayList<FieldFunctionArguments> roiFFAs = new ArrayList<FieldFunctionArguments>();
for (ROI roi : imageDataROIs) {
roiFFAs.add(new FieldFunctionArguments(ROI_EXTDATA_NAME, ROI_MASK_NAME_PREFIX + roi.getROIName(), new Expression(0.0), VariableType.VOLUME));
}
//
// save PSF as a field data
//
ExternalDataInfo psfExtData = createNewExternalDataInfo(localWorkspace, PSF_EXTDATA_NAME);
savePsfAsExternalData(psf, PSF_EXTDATA_VARNAME, psfExtData.getExternalDataIdentifier(), localWorkspace);
FieldFunctionArguments psfFFA = new FieldFunctionArguments(PSF_EXTDATA_NAME, PSF_EXTDATA_VARNAME, new Expression(0.0), VariableType.VOLUME);
TimeBounds timeBounds = getEstimatedRefTimeBound(experimentalRecoveryTime);
double timeStepVal = REFERENCE_DIFF_DELTAT;
Expression chirpedDiffusionRate = new Expression(REFERENCE_DIFF_RATE_COEFFICIENT + "*(t+" + REFERENCE_DIFF_DELTAT + ")");
BioModel bioModel = createRefSimBioModel(simKey, owner, origin, extent, cellROI, timeStepVal, timeBounds, VAR_NAME, new Expression(initConditionFFA.infix()), psfFFA, chirpedDiffusionRate);
if (progressListener != null) {
progressListener.setMessage("Running Reference Simulation...");
}
// run simulation
Simulation simulation = bioModel.getSimulation(0);
ROIDataGenerator roiDataGenerator = getROIDataGenerator(localWorkspace, imageDataROIs);
simulation.getMathDescription().getPostProcessingBlock().addDataGenerator(roiDataGenerator);
runFVSolverStandalone(new File(localWorkspace.getDefaultSimDataDirectory()), simulation, initialConcentrationExtData.getExternalDataIdentifier(), roiExtData.getExternalDataIdentifier(), psfExtData.getExternalDataIdentifier(), progressListener, true);
KeyValue referenceSimKeyValue = simulation.getVersion().getVersionKey();
VCSimulationIdentifier vcSimID = new VCSimulationIdentifier(referenceSimKeyValue, LocalWorkspace.getDefaultOwner());
VCSimulationDataIdentifier vcSimDataID = new VCSimulationDataIdentifier(vcSimID, 0);
File hdf5File = new File(localWorkspace.getDefaultSimDataDirectory(), vcSimDataID.getID() + SimDataConstants.DATA_PROCESSING_OUTPUT_EXTENSION_HDF5);
// get post processing info (time points, variable sizes)
DataOperation.DataProcessingOutputInfoOP dataOperationInfo = new DataOperation.DataProcessingOutputInfoOP(null, /*no vcDataIdentifier OK*/
false, null);
DataOperationResults.DataProcessingOutputInfo dataProcessingOutputInfo = (DataOperationResults.DataProcessingOutputInfo) DataSetControllerImpl.getDataProcessingOutput(dataOperationInfo, hdf5File);
// get post processing data
DataOperation.DataProcessingOutputDataValuesOP dataOperationDataValues = new DataOperation.DataProcessingOutputDataValuesOP(null, /*no vcDataIdentifier OK*/
ROI_EXTDATA_NAME, TimePointHelper.createAllTimeTimePointHelper(), DataIndexHelper.createSliceDataIndexHelper(0), null, null);
DataOperationResults.DataProcessingOutputDataValues dataProcessingOutputDataValues = (DataOperationResults.DataProcessingOutputDataValues) DataSetControllerImpl.getDataProcessingOutput(dataOperationDataValues, hdf5File);
//
// delete the simulation files
//
//
// remove reference simulation files and field data files
//
File userDir = new File(localWorkspace.getDefaultSimDataDirectory());
File[] oldSimFilesToDelete = getSimulationFileNames(userDir, referenceSimKeyValue);
for (int i = 0; oldSimFilesToDelete != null && i < oldSimFilesToDelete.length; i++) {
oldSimFilesToDelete[i].delete();
}
deleteCanonicalExternalData(localWorkspace, initialConcentrationExtData.getExternalDataIdentifier());
deleteCanonicalExternalData(localWorkspace, roiExtData.getExternalDataIdentifier());
deleteCanonicalExternalData(localWorkspace, roiExtData.getExternalDataIdentifier());
// get ref sim time points (distorted by time dilation acceleration)
double[] rawRefDataTimePoints = dataProcessingOutputInfo.getVariableTimePoints();
// get shifted time points
double[] correctedRefDataTimePoints = shiftTimeForBaseDiffRate(rawRefDataTimePoints);
double[][] refData = dataProcessingOutputDataValues.getDataValues();
//
// for rowColumnResultSet with { "t", "roi1", .... , "roiN" } for reference data
//
int numROIs = imageDataROIs.length;
String[] columnNames = new String[numROIs + 1];
columnNames[0] = "t";
for (int i = 0; i < numROIs; i++) {
columnNames[i + 1] = imageDataROIs[i].getROIName();
}
RowColumnResultSet reducedData = new RowColumnResultSet(columnNames);
for (int i = 0; i < correctedRefDataTimePoints.length; i++) {
double[] row = new double[numROIs + 1];
row[0] = correctedRefDataTimePoints[i];
double[] roiData = refData[i];
for (int j = 0; j < numROIs; j++) {
// roiData[0] is the average over the cell .. postbleach this shouldn't change for pure diffusion
row[j + 1] = roiData[j + 1];
}
reducedData.addRow(row);
}
return reducedData;
}
use of org.vcell.util.ISize in project vcell by virtualcell.
the class RunRefSimulationFastOp method savePsfAsExternalData.
private void savePsfAsExternalData(UShortImage psf, String varName, ExternalDataIdentifier newPsfExtDataID, LocalWorkspace localWorkspace) throws IOException, ImageException, ObjectNotFoundException {
Origin origin = new Origin(0, 0, 0);
Extent ext = new Extent(1, 1, 1);
ISize isize = new ISize(3, 3, 1);
VCImageUncompressed vcImage = new VCImageUncompressed(null, new byte[isize.getXYZ()], ext, isize.getX(), isize.getY(), isize.getZ());
RegionImage regionImage = new RegionImage(vcImage, 0, null, null, RegionImage.NO_SMOOTHING);
CartesianMesh cartesianMesh = CartesianMesh.createSimpleCartesianMesh(origin, ext, isize, regionImage);
int NumTimePoints = 1;
// 8 rois integrated into 1 image
int NumChannels = 1;
short[][][] pixData = new short[NumTimePoints][NumChannels][1];
pixData[0][0] = psf.getPixels();
FieldDataFileOperationSpec fdos = new FieldDataFileOperationSpec();
fdos.opType = FieldDataFileOperationSpec.FDOS_ADD;
fdos.cartesianMesh = cartesianMesh;
fdos.shortSpecData = pixData;
fdos.specEDI = newPsfExtDataID;
fdos.varNames = new String[] { varName };
fdos.owner = LocalWorkspace.getDefaultOwner();
fdos.times = new double[] { 0.0 };
fdos.variableTypes = new VariableType[] { VariableType.VOLUME };
fdos.origin = origin;
fdos.extent = ext;
fdos.isize = isize;
localWorkspace.getDataSetControllerImpl().fieldDataFileOperation(fdos);
}
use of org.vcell.util.ISize in project vcell by virtualcell.
the class FitTimeSeries method fitToGaussian.
static GaussianFitResults fitToGaussian(double init_center_i, double init_center_j, double init_radius2, FloatImage image) {
//
// do some optimization on the image (fitting to a Gaussian)
// set initial guesses from ROI operation.
//
ISize imageSize = image.getISize();
final int num_i = imageSize.getX();
final int num_j = imageSize.getY();
final float[] floatPixels = image.getFloatPixels();
//
// initial guess based on previous fit of ROI
// do gaussian fit in index space for center and standard deviation (later to translate it back to world coordinates)
//
final int window_size = (int) Math.sqrt(init_radius2) * 4;
// final int window_min_i = 0; // (int) Math.max(0, Math.floor(init_center_i - window_size/2));
// final int window_max_i = num_i-1; // (int) Math.min(num_i-1, Math.ceil(init_center_i + window_size/2));
// final int window_min_j = 0; // (int) Math.max(0, Math.floor(init_center_j - window_size/2));
// final int window_max_j = num_j-1; // (int) Math.min(num_j-1, Math.ceil(init_center_j + window_size/2));
final int window_min_i = (int) Math.max(0, Math.floor(init_center_i - window_size / 2));
final int window_max_i = (int) Math.min(num_i - 1, Math.ceil(init_center_i + window_size / 2));
final int window_min_j = (int) Math.max(0, Math.floor(init_center_j - window_size / 2));
final int window_max_j = (int) Math.min(num_j - 1, Math.ceil(init_center_j + window_size / 2));
final int PARAM_INDEX_CENTER_I = 0;
final int PARAM_INDEX_CENTER_J = 1;
final int PARAM_INDEX_K = 2;
final int PARAM_INDEX_HIGH = 3;
final int PARAM_INDEX_RADIUS_SQUARED = 4;
final int NUM_PARAMETERS = 5;
double[] initParameters = new double[NUM_PARAMETERS];
initParameters[PARAM_INDEX_CENTER_I] = init_center_i;
initParameters[PARAM_INDEX_CENTER_J] = init_center_j;
initParameters[PARAM_INDEX_HIGH] = 1.0;
initParameters[PARAM_INDEX_K] = 10;
initParameters[PARAM_INDEX_RADIUS_SQUARED] = init_radius2;
PowellOptimizer optimizer = new PowellOptimizer(1e-4, 1e-1);
MultivariateFunction func = new MultivariateFunction() {
@Override
public double value(double[] point) {
double center_i = point[PARAM_INDEX_CENTER_I];
double center_j = point[PARAM_INDEX_CENTER_J];
double high = point[PARAM_INDEX_HIGH];
double K = point[PARAM_INDEX_K];
double radius2 = point[PARAM_INDEX_RADIUS_SQUARED];
double error2 = 0;
for (int j = window_min_j; j <= window_max_j; j++) {
// double y = j - center_j;
double y = j;
for (int i = window_min_i; i <= window_max_i; i++) {
// double x = i - center_i;
double x = i;
double modelValue = high - FastMath.exp(-K * FastMath.exp(-2 * (x * x + y * y) / radius2));
double imageValue = floatPixels[j * num_i + i];
double error = modelValue - imageValue;
error2 += error * error;
}
}
System.out.println(new GaussianFitResults(center_i, center_j, radius2, K, high, error2));
return error2;
}
};
PointValuePair pvp = optimizer.optimize(new ObjectiveFunction(func), new InitialGuess(initParameters), new MaxEval(100000), GoalType.MINIMIZE);
double[] fittedParamValues = pvp.getPoint();
double fitted_center_i = fittedParamValues[PARAM_INDEX_CENTER_I];
double fitted_center_j = fittedParamValues[PARAM_INDEX_CENTER_J];
double fitted_radius2 = fittedParamValues[PARAM_INDEX_RADIUS_SQUARED];
double fitted_K = fittedParamValues[PARAM_INDEX_K];
double fitted_high = fittedParamValues[PARAM_INDEX_HIGH];
double objectiveFunctionValue = pvp.getValue();
return new GaussianFitResults(fitted_center_i, fitted_center_j, fitted_radius2, fitted_K, fitted_high, objectiveFunctionValue);
}
use of org.vcell.util.ISize in project vcell by virtualcell.
the class FitTimeSeries method fit.
public FitBleachSpotOpResults fit(ROI bleachROI, FloatImage normImage) {
//
// find initial guess by centroid from bleach ROI and total area of bleach ROI (assuming a circle of radius R)
//
int nonzeroPixelsCount = bleachROI.getNonzeroPixelsCount();
ISize size = bleachROI.getISize();
if (size.getZ() > 1) {
throw new RuntimeException("expecting 2D bleach region ROI");
}
short[] pixels = bleachROI.getBinaryPixelsXYZ(1);
long numPixelsInROI = 0;
Extent extent = bleachROI.getRoiImages()[0].getExtent();
double totalX_um = 0.0;
double totalY_um = 0.0;
double total_I = 0.0;
double total_J = 0.0;
int pixelIndex = 0;
for (int i = 0; i < size.getX(); i++) {
double x = extent.getX() * (i + 0.5) / (size.getX() - 1);
for (int j = 0; j < size.getY(); j++) {
if (pixels[pixelIndex] != 0) {
double y = extent.getY() * (j + 0.5) / (size.getY() - 1);
totalX_um += x;
totalY_um += y;
total_I += i;
total_J += j;
numPixelsInROI++;
}
pixelIndex++;
}
}
Origin origin = bleachROI.getRoiImages()[0].getOrigin();
double roiCenterX_um = origin.getX() + totalX_um / numPixelsInROI;
double roiCenterY_um = origin.getY() + totalY_um / numPixelsInROI;
double roiCenterI_pixelscale = total_I / numPixelsInROI;
double roiCenterJ_pixelscale = total_J / numPixelsInROI;
// Area = PI * R^2
// R = sqrt(Area/PI)
double roiBleachSpotArea_um2 = (extent.getX() * extent.getY() * numPixelsInROI) / (size.getX() * size.getY());
double roiBleachRadius_um = Math.sqrt(roiBleachSpotArea_um2 / Math.PI);
double roiBleachRadius_pixelscale = Math.sqrt(numPixelsInROI / Math.PI);
FitBleachSpotOpResults results = new FitBleachSpotOpResults();
results.bleachRadius_ROI = roiBleachRadius_um;
results.centerX_ROI = roiCenterX_um;
results.centerY_ROI = roiCenterY_um;
GaussianFitResults gfresults = fitToGaussian(roiCenterI_pixelscale, roiCenterJ_pixelscale, roiBleachRadius_pixelscale * roiBleachRadius_pixelscale, normImage);
results.bleachFactorK_GaussianFit = gfresults.K;
results.bleachRadius_GaussianFit = Math.sqrt(gfresults.radius2);
results.centerX_GaussianFit = origin.getX() + extent.getX() * (gfresults.centerI + 0.5) / size.getX();
results.centerY_GaussianFit = origin.getY() + extent.getY() * (gfresults.centerJ + 0.5) / size.getY();
return results;
}
use of org.vcell.util.ISize in project vcell by virtualcell.
the class ImageJHelper method vcellSendImage.
public static void vcellSendImage(final Component requester, final PDEDataContext pdeDataContext, SubVolume subvolume, Hashtable<SampledCurve, int[]>[] membraneTables, String description, double[] timePoints, String[] channelDescriptions, int[] colormap) throws Exception {
// xyz, 1 time, 1 var
final ImageJConnection[] imageJConnectionArr = new ImageJConnection[1];
AsynchClientTask sendImageTask = new AsynchClientTask("Send image to ImageJ...", AsynchClientTask.TASKTYPE_NONSWING_BLOCKING) {
@Override
public void run(Hashtable<String, Object> hashTable) throws Exception {
try {
ImageJConnection imageJConnection = new ImageJConnection(ExternalCommunicator.IMAGEJ);
imageJConnectionArr[0] = imageJConnection;
imageJConnection.openConnection(VCellImageJCommands.vcellSendImage, description);
// send size of the standard 5 dimensions in this order (width, height, nChannels, nSlices, nFrames)
ISize xyzSize = pdeDataContext.getCartesianMesh().getISize();
Extent extent = pdeDataContext.getCartesianMesh().getExtent();
BasicStackDimensions basicStackDimensions = new BasicStackDimensions(xyzSize.getX(), xyzSize.getY(), xyzSize.getZ(), 1, 1);
sendData0(imageJConnection, new HyperStackHelper(basicStackDimensions, extent, true, Float.class.getSimpleName(), null, new int[] { subvolume.getHandle() }, timePoints, channelDescriptions, colormap), pdeDataContext.getDataValues(), "'" + pdeDataContext.getVariableName() + "'" + pdeDataContext.getTimePoint());
sendVolumeDomain0(imageJConnection, pdeDataContext.getCartesianMesh(), null, description);
sendMembraneOutline(imageJConnection, membraneTables);
} catch (Exception e) {
if (e instanceof UserCancelException) {
throw e;
}
e.printStackTrace();
hashTable.put("imagejerror", e);
} finally {
try {
if (imageJConnectionArr[0] != null) {
imageJConnectionArr[0].closeConnection();
}
} catch (Exception e) {
e.printStackTrace();
}
}
}
};
ClientTaskDispatcher.dispatch(requester, new Hashtable<>(), new AsynchClientTask[] { sendImageTask }, false, true, new ProgressDialogListener() {
@Override
public void cancelButton_actionPerformed(EventObject newEvent) {
if (imageJConnectionArr[0] != null) {
imageJConnectionArr[0].closeConnection();
}
}
});
}
Aggregations